Kvantový aproximační optimalizační algoritmus
Odhadovaná doba použití: 22 minut na procesoru Heron r3 (POZNÁMKA: Jedná se pouze o odhad. Skutečná doba běhu se může lišit.)
Pozadí
Tento tutoriál ukazuje, jak implementovat Kvantový aproximační optimalizační algoritmus (QAOA) – hybridní (kvantově-klasickou) iterativní metodu – v kontextu Qiskit vzorů. Nejprve vyřešíš problém maximálního řezu (neboli Max-Cut) pro malý graf a poté se naučíš, jak ho spustit v užitkové škále. Všechna hardwarová spuštění v tutoriálu by měla proběhnout v časovém limitu volně dostupného Open Plánu.
Problém Max-Cut je optimalizační problém, který je obtížné vyřešit (konkrétněji jde o NP-těžký problém) s řadou různých aplikací v oblasti clusterování, síťové vědy a statistické fyziky. Tento tutoriál uvažuje graf uzlů propojených hranami a snaží se rozdělit uzly do dvou množin tak, aby byl maximalizován počet hran, které tento řez překračuje.

Požadavky
Před zahájením tohoto tutoriálu se ujisti, že máš nainstalováno následující:
- Qiskit SDK v1.0 nebo novější s podporou vizualizace
- Qiskit Runtime v0.22 nebo novější (
pip install qiskit-ibm-runtime)
Kromě toho budeš potřebovat přístup k instanci na IBM Quantum Platform. Vezmi na vědomí, že tento tutoriál nelze spustit na Open Plánu, protože používá workloady prostřednictvím sessions, které jsou dostupné pouze s přístupem Premium Plánu.
Nastavení
# Added by doQumentation — required packages for this notebook
!pip install -q matplotlib numpy qiskit qiskit-ibm-runtime rustworkx scipy
import matplotlib
import matplotlib.pyplot as plt
import rustworkx as rx
from rustworkx.visualization import mpl_draw as draw_graph
import numpy as np
from scipy.optimize import minimize
from collections import defaultdict
from typing import Sequence
from qiskit.quantum_info import SparsePauliOp
from qiskit.circuit.library import QAOAAnsatz
from qiskit.transpiler.preset_passmanagers import generate_preset_pass_manager
from qiskit_ibm_runtime import QiskitRuntimeService
from qiskit_ibm_runtime import Session, EstimatorV2 as Estimator
from qiskit_ibm_runtime import SamplerV2 as Sampler
Část I. QAOA v malém měřítku
První část tohoto tutoriálu používá problém Max-Cut v malém měřítku k ilustraci kroků pro řešení optimalizačního problému pomocí kvantového počítače.
Aby sis lépe porozuměl kontextu před mapováním tohoto problému na kvantový algoritmus, lze lépe pochopit, jak se problém Max-Cut stává klasickým kombinatorickým optimalizačním problémem, když nejprve uvažujeme minimalizaci funkce
kde vstup je vektor, jehož složky odpovídají jednotlivým uzlům grafu. Pak omez každou z těchto složek na hodnotu nebo (které představují zahrnutí nebo nezahrnutí do řezu). Tento příklad v malém měřítku používá graf s uzly.
Mohl bys napsat funkci dvojice uzlů , která indikuje, zda odpovídající hrana leží v řezu. Například funkce je 1, pouze pokud jedno z nebo je 1 (což znamená, že hrana leží v řezu), a jinak je nulová. Problém maximalizace hran v řezu lze formulovat jako
což lze přepsat jako minimalizaci ve tvaru
Minimum v tomto případě nastane, když je počet hran překračovaných řezem maximální. Jak vidíš, dosud tu není nic, co by se týkalo kvantového počítání. Je potřeba přeformulovat tento problém do podoby, které kvantový počítač rozumí. Inicializuj svůj problém vytvořením grafu s uzly.
n = 5
graph = rx.PyGraph()
graph.add_nodes_from(np.arange(0, n, 1))
edge_list = [
(0, 1, 1.0),
(0, 2, 1.0),
(0, 4, 1.0),
(1, 2, 1.0),
(2, 3, 1.0),
(3, 4, 1.0),
]
graph.add_edges_from(edge_list)
draw_graph(graph, node_size=600, with_labels=True)
Krok 1: Mapování klasických vstupů na kvantový problém
Prvním krokem vzoru je mapování klasického problému (grafu) na kvantové Circuit a operátory. K tomu jsou potřeba tři hlavní kroky:
- Využít řadu matematických přeformulování k reprezentaci tohoto problému pomocí notace problémů Kvadratické neomezené binární optimalizace (QUBO).
- Přepsat optimalizační problém jako Hamiltonián, jehož základní stav odpovídá řešení minimalizujícímu nákladovou funkci.
- Vytvořit kvantový Circuit, který připraví základní stav tohoto Hamiltoniánu procesem podobným kvantovému žíhání.
Poznámka: V metodologii QAOA chceš nakonec mít operátor (Hamiltonián), který reprezentuje nákladovou funkci našeho hybridního algoritmu, a také parametrizovaný Circuit (Ansatz), který reprezentuje kvantové stavy s kandidátními řešeními problému. Z těchto kandidátních stavů lze vzorkovat a poté je vyhodnocovat pomocí nákladové funkce.
Graf → optimalizační problém
Prvním krokem mapování je změna notace. Následující výraz vyjadřuje problém v QUBO notaci:
kde je matice reálných čísel, odpovídá počtu uzlů v grafu, je vektor binárních proměnných zavedených výše a označuje transpozici vektoru .
Maximize
-2*x_0*x_1 - 2*x_0*x_2 - 2*x_0*x_4 - 2*x_1*x_2 - 2*x_2*x_3 - 2*x_3*x_4 + 3*x_0
+ 2*x_1 + 3*x_2 + 2*x_3 + 2*x_4
Subject to
No constraints
Binary variables (5)
x_0 x_1 x_2 x_3 x_4
Optimalizační problém → Hamiltonián
Poté lze QUBO problém přeformulovat jako Hamiltonián (zde matice reprezentující energii systému):
Kroky přeformulování z QAOA problému na Hamiltonián
Abychom ukázali, jak lze QAOA problém přepsat tímto způsobem, nejprve nahraď binární proměnné novou sadou proměnných pomocí
Zde vidíš, že pokud je rovno , pak musí být . Když se nahradí v optimalizačním problému (), lze získat ekvivalentní formulaci.
Pokud nyní definujeme , odstraníme prefaktor a konstantní člen , dostaneme dvě ekvivalentní formulace stejného optimalizačního problému.
Zde závisí na . Všimni si, že pro získání jsme vypustili faktor 1/4 a konstantní posun , které nehrají roli v optimalizaci.
Nyní, pro získání kvantové formulace problému, povyš proměnné na Pauliho matici , jako je matice ve tvaru
Když tyto matice dosadíš do výše uvedeného optimalizačního problému, získáš následující Hamiltonián
Také si vzpomeň, že matice jsou vloženy do výpočetního prostoru kvantového počítače, tedy do Hilbertova prostoru velikosti . Proto bys měl chápat výrazy jako jako tenzorový součin vložený do Hilbertova prostoru . Například v problému s pěti rozhodovacími proměnnými je výraz chápán jako , kde je jednotková matice .
Tento Hamiltonián se nazývá Hamiltonián nákladové funkce. Má tu vlastnost, že jeho základní stav odpovídá řešení, které minimalizuje nákladovou funkci . Aby sis tedy vyřešil optimalizační problém, potřebuješ nyní připravit základní stav (nebo stav s vysokým překryvem s ním) na kvantovém počítači. Vzorkování z tohoto stavu pak s vysokou pravděpodobností přinese řešení . Nyní uvažme Hamiltonián pro problém Max-Cut. Přiřaď každý vrchol grafu ke Qubitu ve stavu nebo , kde hodnota označuje množinu, do které vrchol patří. Cílem problému je maximalizovat počet hran , pro které platí a , nebo naopak. Pokud přiřadíme operátor každému Qubitu, kde
pak hrana patří do řezu, pokud vlastní hodnota ; jinými slovy, Qubity přiřazené k a jsou různé. Podobně nepatří do řezu, pokud vlastní hodnota . Všimni si, že nás nezajímá přesný stav Qubitu přiřazený každému vrcholu, ale pouze to, zda jsou na hraně stejné nebo různé. Problém Max-Cut vyžaduje najít přiřazení Qubitů na vrcholech tak, aby byla minimalizována vlastní hodnota následujícího Hamiltoniánu
Jinými slovy, pro všechna v problému Max-Cut. Hodnota označuje váhu hrany. V tomto tutoriálu uvažujeme neváhovaný graf, tedy pro všechna .
def build_max_cut_paulis(
graph: rx.PyGraph,
) -> list[tuple[str, list[int], float]]:
"""Convert the graph to Pauli list.
This function does the inverse of `build_max_cut_graph`
"""
pauli_list = []
for edge in list(graph.edge_list()):
weight = graph.get_edge_data(edge[0], edge[1])
pauli_list.append(("ZZ", [edge[0], edge[1]], weight))
return pauli_list
max_cut_paulis = build_max_cut_paulis(graph)
cost_hamiltonian = SparsePauliOp.from_sparse_list(max_cut_paulis, n)
print("Cost Function Hamiltonian:", cost_hamiltonian)
Cost Function Hamiltonian: SparsePauliOp(['IIIZZ', 'IIZIZ', 'ZIIIZ', 'IIZZI', 'IZZII', 'ZZIII'],
coeffs=[1.+0.j, 1.+0.j, 1.+0.j, 1.+0.j, 1.+0.j, 1.+0.j])
Hamiltonián → kvantový Circuit
Hamiltonián obsahuje kvantovou definici tvého problému. Nyní můžeš vytvořit kvantový Circuit, který pomůže vzorkovat dobrá řešení z kvantového počítače. QAOA je inspirován kvantovým žíháním a aplikuje střídavé vrstvy operátorů v kvantovém Circuit.
Obecná myšlenka je začít v základním stavu známého systému, výše, a poté řídit systém do základního stavu nákladového operátoru, o který máš zájem. To se provádí aplikací operátorů a s úhly a .
Kvantový Circuit, který vygeneruješ, je parametrizován hodnotami a , takže můžeš zkoušet různé hodnoty a a vzorkovat z výsledného stavu.
V tomto případě vyzkoušíš příklad s jednou vrstvou QAOA obsahující dva parametry: a .
circuit = QAOAAnsatz(cost_operator=cost_hamiltonian, reps=2)
circuit.measure_all()
circuit.draw("mpl")
circuit.parameters
ParameterView([ParameterVectorElement(β[0]), ParameterVectorElement(β[1]), ParameterVectorElement(γ[0]), ParameterVectorElement(γ[1])])
Krok 2: Optimalizace problému pro spuštění na kvantovém hardwaru
Circuit výše obsahuje řadu abstrakcí užitečných pro přemýšlení o kvantových algoritmech, ale na hardwaru je nelze spustit. Aby bylo možné spustit na QPU, musí Circuit projít řadou operací, které tvoří krok transpilace nebo optimalizace Circuit vzoru.
Knihovna Qiskit nabízí řadu transpilačních průchodů, které pokrývají širokou škálu transformací Circuit. Je potřeba zajistit, aby byl tvůj Circuit optimalizován pro daný účel.
Transpilace může zahrnovat několik kroků, jako jsou:
- Počáteční mapování Qubitů v Circuit (jako jsou rozhodovací proměnné) na fyzické Qubity na zařízení.
- Rozvinutí instrukcí v kvantovém Circuit do nativních instrukcí hardwaru, kterým Backend rozumí.
- Směrování Qubitů v Circuit, které spolu interagují, na fyzické Qubity, které jsou vzájemně sousední.
- Potlačení chyb přidáním jednoQubitových Gate pro potlačení šumu pomocí dynamického oddělování.
Více informací o transpilaci najdeš v naší dokumentaci.
Následující kód transformuje a optimalizuje abstraktní Circuit do formátu připraveného ke spuštění na jednom ze zařízení přístupných přes cloud pomocí Qiskit IBM Runtime služby.
service = QiskitRuntimeService()
backend = service.least_busy(
operational=True, simulator=False, min_num_qubits=127
)
print(backend)
# Create pass manager for transpilation
pm = generate_preset_pass_manager(optimization_level=3, backend=backend)
candidate_circuit = pm.run(circuit)
candidate_circuit.draw("mpl", fold=False, idle_wires=False)
<IBMBackend('test_heron_pok_1')>

Krok 3: Spuštění pomocí Qiskit primitiv
V pracovním postupu QAOA jsou optimální parametry QAOA nalezeny v iterační optimalizační smyčce, která spouští řadu vyhodnocení Circuit a používá klasický optimalizátor k nalezení optimálních parametrů a . Tato spouštěcí smyčka se provádí pomocí následujících kroků:
- Definuj počáteční parametry
- Vytvoř novou
Sessionobsahující optimalizační smyčku a primitivum použité pro vzorkování Circuit - Jakmile je nalezena optimální sada parametrů, spusť Circuit naposledy pro získání finálního rozdělení, které bude použito v kroku po-zpracování.
Definuj Circuit s počátečními parametry
Začínáme s libovolně zvolenými parametry.
initial_gamma = np.pi
initial_beta = np.pi / 2
init_params = [initial_beta, initial_beta, initial_gamma, initial_gamma]
Definuj Backend a spouštěcí primitivum
Použij Qiskit Runtime primitiva pro interakci s IBM® Backendy. Dvě primitiva jsou Sampler a Estimator a volba primitiva závisí na tom, jaký typ měření chceš spustit na kvantovém počítači. Pro minimalizaci použij Estimator, protože měření nákladové funkce je jednoduše střední hodnota .
Spustit
Primitiva nabízejí různé režimy spuštění pro plánování workloadů na kvantových zařízeních a pracovní postup QAOA běží iterativně v Session.

Funkci nákladů založenou na Sampleru můžeš zapojit do minimalizační rutiny SciPy pro nalezení optimálních parametrů.
def cost_func_estimator(params, ansatz, hamiltonian, estimator):
# transform the observable defined on virtual qubits to
# an observable defined on all physical qubits
isa_hamiltonian = hamiltonian.apply_layout(ansatz.layout)
pub = (ansatz, isa_hamiltonian, params)
job = estimator.run([pub])
results = job.result()[0]
cost = results.data.evs
objective_func_vals.append(cost)
return cost
objective_func_vals = [] # Global variable
with Session(backend=backend) as session:
# If using qiskit-ibm-runtime<0.24.0, change `mode=` to `session=`
estimator = Estimator(mode=session)
estimator.options.default_shots = 1000
# Set simple error suppression/mitigation options
estimator.options.dynamical_decoupling.enable = True
estimator.options.dynamical_decoupling.sequence_type = "XY4"
estimator.options.twirling.enable_gates = True
estimator.options.twirling.num_randomizations = "auto"
result = minimize(
cost_func_estimator,
init_params,
args=(candidate_circuit, cost_hamiltonian, estimator),
method="COBYLA",
tol=1e-2,
)
print(result)
message: Return from COBYLA because the trust region radius reaches its lower bound.
success: True
status: 0
fun: -1.6295230263157894
x: [ 1.530e+00 1.439e+00 4.071e+00 4.434e+00]
nfev: 26
maxcv: 0.0
Optimalizátoru se podařilo snížit náklady a najít lepší parametry pro Circuit.
plt.figure(figsize=(12, 6))
plt.plot(objective_func_vals)
plt.xlabel("Iteration")
plt.ylabel("Cost")
plt.show()
Jakmile jsi nalezl optimální parametry pro Circuit, můžeš je přiřadit a vzorkovat finální rozdělení získané s optimalizovanými parametry. Zde by mělo být použito primitivum Sampler, protože jde o rozdělení pravděpodobnosti bitřetězcových měření, která odpovídají optimálnímu řezu grafu.
Poznámka: To znamená připravit kvantový stav v počítači a poté ho změřit. Měření zkolabuje stav do jediného výpočetního základního stavu – například 010101110000... – který odpovídá kandidátnímu řešení našeho původního optimalizačního problému ( nebo v závislosti na úloze).
optimized_circuit = candidate_circuit.assign_parameters(result.x)
optimized_circuit.draw("mpl", fold=False, idle_wires=False)

# If using qiskit-ibm-runtime<0.24.0, change `mode=` to `backend=`
sampler = Sampler(mode=backend)
sampler.options.default_shots = 10000
# Set simple error suppression/mitigation options
sampler.options.dynamical_decoupling.enable = True
sampler.options.dynamical_decoupling.sequence_type = "XY4"
sampler.options.twirling.enable_gates = True
sampler.options.twirling.num_randomizations = "auto"
pub = (optimized_circuit,)
job = sampler.run([pub], shots=int(1e4))
counts_int = job.result()[0].data.meas.get_int_counts()
counts_bin = job.result()[0].data.meas.get_counts()
shots = sum(counts_int.values())
final_distribution_int = {key: val / shots for key, val in counts_int.items()}
final_distribution_bin = {key: val / shots for key, val in counts_bin.items()}
print(final_distribution_int)
{28: 0.0328, 11: 0.0343, 2: 0.0296, 25: 0.0308, 16: 0.0303, 27: 0.0302, 13: 0.0323, 7: 0.0312, 4: 0.0296, 9: 0.0295, 26: 0.0321, 30: 0.031, 23: 0.0324, 31: 0.0303, 21: 0.0335, 15: 0.0317, 12: 0.0309, 29: 0.0297, 3: 0.0313, 5: 0.0312, 6: 0.0274, 10: 0.0329, 22: 0.0353, 0: 0.0315, 20: 0.0326, 8: 0.0322, 14: 0.0306, 17: 0.0295, 18: 0.0279, 1: 0.0325, 24: 0.0334, 19: 0.0295}
Krok 4: Po-zpracování a vrácení výsledku v požadovaném klasickém formátu
Krok po-zpracování interpretuje výstup vzorkování a vrací řešení původního problému. V tomto případě tě zajímá bitřetězec s nejvyšší pravděpodobností, protože ten určuje optimální řez. Symetrie v problému umožňují čtyři možná řešení a proces vzorkování vrátí jedno z nich s mírně vyšší pravděpodobností, ale v níže zobrazeném rozdělení vidíš, že čtyři z bitřetězců jsou zřetelně pravděpodobnější než ostatní.
# auxiliary functions to sample most likely bitstring
def to_bitstring(integer, num_bits):
result = np.binary_repr(integer, width=num_bits)
return [int(digit) for digit in result]
keys = list(final_distribution_int.keys())
values = list(final_distribution_int.values())
most_likely = keys[np.argmax(np.abs(values))]
most_likely_bitstring = to_bitstring(most_likely, len(graph))
most_likely_bitstring.reverse()
print("Result bitstring:", most_likely_bitstring)
Result bitstring: [0, 1, 1, 0, 1]
matplotlib.rcParams.update({"font.size": 10})
final_bits = final_distribution_bin
values = np.abs(list(final_bits.values()))
top_4_values = sorted(values, reverse=True)[:4]
positions = []
for value in top_4_values:
positions.append(np.where(values == value)[0])
fig = plt.figure(figsize=(11, 6))
ax = fig.add_subplot(1, 1, 1)
plt.xticks(rotation=45)
plt.title("Result Distribution")
plt.xlabel("Bitstrings (reversed)")
plt.ylabel("Probability")
ax.bar(list(final_bits.keys()), list(final_bits.values()), color="tab:grey")
for p in positions:
ax.get_children()[int(p[0])].set_color("tab:purple")
plt.show()
Vizualizace nejlepšího řezu
Z optimálního bitřetězce pak můžeš vizualizovat tento řez na původním grafu.
# auxiliary function to plot graphs
def plot_result(G, x):
colors = ["tab:grey" if i == 0 else "tab:purple" for i in x]
pos, _default_axes = rx.spring_layout(G), plt.axes(frameon=True)
rx.visualization.mpl_draw(
G, node_color=colors, node_size=100, alpha=0.8, pos=pos
)
plot_result(graph, most_likely_bitstring)
A vypočítej hodnotu řezu:
def evaluate_sample(x: Sequence[int], graph: rx.PyGraph) -> float:
assert len(x) == len(
list(graph.nodes())
), "The length of x must coincide with the number of nodes in the graph."
return sum(
x[u] * (1 - x[v]) + x[v] * (1 - x[u])
for u, v in list(graph.edge_list())
)
cut_value = evaluate_sample(most_likely_bitstring, graph)
print("The value of the cut is:", cut_value)
The value of the cut is: 5
Část II. Škálování na větší problémy!
Máš přístup k mnoha zařízením s více než 100 qubity na platformě IBM Quantum®. Vyber jedno, na kterém vyřešíš Max-Cut na váženém grafu se 100 uzly. Jedná se o problém v „utility-scale" měřítku. Kroky k sestavení pracovního postupu jsou stejné jako výše, ale s mnohem větším grafem.
n = 100 # Number of nodes in graph
graph_100 = rx.PyGraph()
graph_100.add_nodes_from(np.arange(0, n, 1))
elist = []
for edge in backend.coupling_map:
if edge[0] < n and edge[1] < n:
elist.append((edge[0], edge[1], 1.0))
graph_100.add_edges_from(elist)
draw_graph(graph_100, node_size=200, with_labels=True, width=1)
