Demo 2：QAOA MaxCut 4-node

第 8 章 量子机器学习 · 本章演示 · 代码演示总索引

可视化：双图（热图+top bitstring） —— 14×14 网格+爬山轨迹

对应教程：第 8 章 量子机器学习 中相关小节 · 依赖安装：pip install "unified-quantum[visualization,simulation]"

实现要点

核心：K4 是 4 节点完全图，3-正则，含 6 条边；3-正则图上 QAOA p=1 的 最坏情形近似比下界为 r_1 = 0.6924（Farhi et al., 2014）。 这里在 (γ, β) 网格上扫描 p=1 QAOA 期望成本 ⟨H_C⟩，画出等高线 + 经典爬山优化轨迹；K4 这个特定实例上 QAOA p=1 能达到 ~0.92 的 实际近似比（远好于最坏情形保证）。

完整源码

demos/ch08-qml/02_qaoa_maxcut_4node.py

"""
ch08-qml demo 2: QAOA on the 4-node MaxCut (K4)

对应教程：chenzhaoyun.com/quantum-tutorial/ch08-qml/qaoa-tutorial.html
核心：K4 是 4 节点完全图，3-正则，含 6 条边；3-正则图上 QAOA p=1 的
     最坏情形近似比下界为 r_1 = 0.6924（Farhi et al., 2014）。
     这里在 (γ, β) 网格上扫描 p=1 QAOA 期望成本 ⟨H_C⟩，画出等高线 +
     经典爬山优化轨迹；K4 这个特定实例上 QAOA p=1 能达到 ~0.92 的
     实际近似比（远好于最坏情形保证）。

依赖：pip install unified-quantum[visualization]
运行：python demos/ch08-qml/02_qaoa_maxcut_4node.py
"""
from __future__ import annotations

from pathlib import Path

import matplotlib
matplotlib.use("Agg")
import matplotlib.pyplot as plt
import numpy as np

from uniqc import qaoa_ansatz, submit_task, wait_for_result

OUT = Path(__file__).parent

# === 1. K4 graph and problem Hamiltonian ===
# Edges of K4: every pair of the 4 nodes
EDGES = [(0, 1), (0, 2), (0, 3), (1, 2), (1, 3), (2, 3)]
N_QUBITS = 4
N_EDGES = len(EDGES)
OPT_CUT = 4  # K4 MaxCut optimal: partition {0,1}|{2,3} → 4 cut edges

# Cost Hamiltonian H_C = sum_{(i,j)∈E} (1 - Z_i Z_j)/2
# (1 - Z_i Z_j)/2 contributes +1 iff cut. The constant 1/2 per edge shifts
# energy by |E|/2; we drop the constant and encode with negative Z_i Z_j
# coefficients so that QAOA maximises <H_C>.
H_C = [(f"Z{i}Z{j}", -0.5) for i, j in EDGES]

SHOTS = 1500


def evaluate_qaoa(gamma: float, beta: float, shots: int = SHOTS) -> tuple[float, dict]:
    """Build the p=1 QAOA circuit, measure, return (mean_cut_count, raw_counts)."""
    c = qaoa_ansatz(H_C, p=1, gammas=np.array([gamma]), betas=np.array([beta]))
    for q in range(N_QUBITS):
        c.measure(q)
    r = wait_for_result(submit_task(circuit=c, backend="dummy:local:simulator", shots=shots))
    counts = r.counts
    total = sum(counts.values())
    if total <= 0:
        return 0.0, counts
    cut_count = 0.0
    for i, j in EDGES:
        p_same = 0.0
        for bs, cnt in counts.items():
            bs_str = bs.zfill(N_QUBITS)
            qi = int(bs_str[N_QUBITS - 1 - i])
            qj = int(bs_str[N_QUBITS - 1 - j])
            if qi == qj:
                p_same += cnt
        p_same /= total
        cut_count += 1.0 - p_same
    return cut_count, counts


# === 2. γ × β grid scan (one quantum submission per cell) ===
print("ch08-qml-demo-2: QAOA p=1 on K4 MaxCut")
print(f"  nodes={N_QUBITS}, edges={N_EDGES}, optimal cut={OPT_CUT}")
print(f"  shots per QAOA run: {SHOTS}")

G = np.linspace(0.05, 1.4, 14)  # γ
B = np.linspace(0.05, 0.95, 14)  # β
cut_grid = np.full((len(G), len(B)), np.nan)
total_evals = 0
for gi, g in enumerate(G):
    for bi, b in enumerate(B):
        c, _ = evaluate_qaoa(g, b)
        cut_grid[gi, bi] = c
        total_evals += 1

# === 3. classical hill-climb starting from the grid optimum ===
gi_best, bi_best = np.unravel_index(int(np.nanargmax(cut_grid)), cut_grid.shape)
g_best, b_best = float(G[gi_best]), float(B[bi_best])
c_best, best_counts = evaluate_qaoa(g_best, b_best, shots=8000)

# A short COBYLA-like refinement: gradient-free random walk
rng = np.random.default_rng(0)
trajectory = [(g_best, b_best, c_best)]
for step in range(20):
    sigma = max(0.03, 0.12 * (0.6 ** step))
    g_try = g_best + rng.normal(0, sigma)
    b_try = b_best + rng.normal(0, sigma)
    g_try = max(0.01, min(np.pi, g_try))
    b_try = max(0.01, min(np.pi, b_try))
    c_try, _ = evaluate_qaoa(g_try, b_try, shots=4000)
    if c_try > c_best:
        g_best, b_best, c_best = g_try, b_try, c_try
        best_counts = _
    trajectory.append((g_best, b_best, c_best))

# === 4. Build bitstring histogram for the best run ===
hist = {bs.zfill(N_QUBITS): cnt for bs, cnt in best_counts.items()}
total_best = sum(hist.values())
hist_sorted = sorted(hist.items(), key=lambda kv: -kv[1])
top_states = hist_sorted[:8]

# best observed single sample's cut
def cut_of(bs: str) -> int:
    bs = bs.zfill(N_QUBITS)
    return sum(1 for i, j in EDGES if int(bs[N_QUBITS - 1 - i]) != int(bs[N_QUBITS - 1 - j]))

best_sample_cut = max(cut_of(bs) for bs, _ in top_states)

# === 5. matplotlib 可视化 ===
fig = plt.figure(figsize=(12.5, 5.4))
gs = fig.add_gridspec(1, 2, width_ratios=[1.05, 1.0])

# Left: heatmap of ⟨cuts⟩ with optimal marker
ax1 = fig.add_subplot(gs[0, 0])
im = ax1.imshow(
    cut_grid, origin="lower", aspect="auto", cmap="viridis",
    extent=[B[0] - (B[1] - B[0]) / 2, B[-1] + (B[1] - B[0]) / 2,
            G[0] - (G[1] - G[0]) / 2, G[-1] + (G[1] - G[0]) / 2],
)
ax1.plot([b_best], [g_best], "*", color="#ffd700", ms=18, mec="black", mew=0.9,
         label=fr"optimum $\gamma^*={g_best:.3f}$, $\beta^*={b_best:.3f}$")
cs = ax1.contour(B, G, cut_grid, levels=[1.0, 2.0, 3.0, 3.5], colors="white",
                 linewidths=0.8, alpha=0.55)
ax1.clabel(cs, inline=True, fontsize=8, fmt="%.1f")
# trajectory line
traj_g = [t[0] for t in trajectory]
traj_b = [t[1] for t in trajectory]
ax1.plot(traj_b, traj_g, "-o", color="#ff6f61", ms=3, lw=1.0, alpha=0.75,
         label="classical hill-climb trajectory")
ax1.set_xlabel(r"mixer angle $\beta$")
ax1.set_ylabel(r"cost angle $\gamma$")
ax1.set_title(r"p=1 QAOA $\langle H_C\rangle$ landscape on K4")
cbar = plt.colorbar(im, ax=ax1)
cbar.set_label(r"mean # cuts (out of 4)")
ax1.legend(loc="lower left", fontsize=8)

# Right: histogram of top bitstrings at (γ*, β*)
ax2 = fig.add_subplot(gs[0, 1])
labels = [bs for bs, _ in top_states]
probs = [cnt / total_best for _, cnt in top_states]
cuts = [cut_of(bs) for bs, _ in top_states]
colors = ["#2ca02c" if c == OPT_CUT else "#1f77b4" if c == OPT_CUT - 1 else "#d62728" for c in cuts]
bars = ax2.bar(range(len(labels)), probs, color=colors, ec="black", lw=0.5)
for k, (lab, p, c) in enumerate(zip(labels, probs, cuts)):
    ax2.text(k, p + 0.012, f"{lab}\ncut={c}", ha="center", va="bottom",
             fontsize=7.5)
ax2.set_xticks(range(len(labels)))
ax2.set_xticklabels(labels, rotation=0, fontsize=8)
ax2.set_ylabel("sampling probability")
ax2.set_title(fr"top bitstrings at $(\gamma^*,\beta^*)$ — best cut in samples = {best_sample_cut}/{OPT_CUT}")
ax2.set_ylim(0, max(probs) * 1.45)
ax2.grid(alpha=0.3, axis="y")

plt.suptitle("ch08-qml demo 2: QAOA p=1 on K4 MaxCut", y=1.02)
plt.tight_layout()
out_path = OUT / "02_qaoa_maxcut_4node.png"
plt.savefig(out_path, dpi=120, bbox_inches="tight")
plt.close()

# === 6. stdout 关键统计 ===
print(f"  total QAOA circuit evaluations: {total_evals + 1 + (len(trajectory) - 1)}")
print(f"  grid scan size: {len(G)}x{len(B)} = {len(G)*len(B)} cells")
print(f"  grid best      : cuts={float(cut_grid.max()):.3f}/{OPT_CUT} "
      f"at γ={G[gi_best]:.3f}, β={B[bi_best]:.3f}")
print(f"  after hill-climb: γ*={g_best:.4f}, β*={b_best:.4f}")
print(f"  mean cut count  : {c_best:.4f} / {OPT_CUT}")
print(f"  approximation r : {c_best/OPT_CUT:.4f}  (3-regular worst-case guarantee = 0.6924)")
print(f"  best single sample cut: {best_sample_cut} / {OPT_CUT}")
print(f"  top 3 bitstrings:")
for bs, cnt in hist_sorted[:3]:
    print(f"    |{bs}>  → cut={cut_of(bs)}, prob={cnt/total_best:.4f}")
print(f"  PNG -> {out_path}")

运行结果

运行方式

# 仓库根目录执行；输出 PNG 到 demos/ch08-qml/02_qaoa_maxcut_4node.png，关键统计到 stdout
.venv/bin/python3 demos/ch08-qml/02_qaoa_maxcut_4node.py

---

← 本章演示 · ↑ 代码演示总索引 · ⌂ 教程门户