Přeskočit na hlavní obsah

Kvantový aproximační optimalizační algoritmus

Odhad využití: 22 minut na procesoru Heron r3 (POZNÁMKA: Jedná se pouze o odhad. Skutečná doba běhu se může lišit.)

Výsledky učení

Po dokončení tohoto tutoriálu můžeš očekávat, že porozumíš následujícím informacím:

  • Jak mapovat klasický kombinatorický optimalizační problém (max-cut) na kvantový Hamiltonián
  • Jak implementovat a spustit Kvantový aproximační optimalizační algoritmus (QAOA) pomocí Qiskit Runtime sessions
  • Jak škálovat pracovní postup QAOA z malého příkladu na simulátoru na spuštění na hardwaru v užitkové škále

Předpoklady

Doporučuje se, aby ses seznámil s těmito tématy:

Pozadí

Kvantový aproximační optimalizační algoritmus (QAOA) je hybridní kvantově-klasická iterativní metoda pro řešení kombinatorických optimalizačních problémů. V tomto tutoriálu použiješ QAOA k vyřešení problému maximálního řezu (max-cut) — NP-těžkého optimalizačního problému s aplikacemi v clusterování, síťové vědě a statistické fyzice. Pro daný graf uzlů propojených hranami je cílem rozdělit uzly do dvou množin tak, aby byl maximalizován počet hran, které řez překračují.

Ilustrace problému max-cut

Od klasické optimalizace ke kvantovým Circuit

Max-cut lze vyjádřit jako klasický binární optimalizační problém. Každému uzlu je přiřazena binární proměnná xi{0,1}x_i \in \{0, 1\} označující, do které množiny patří. Cílem je maximalizovat počet hran, jejichž koncové body jsou v různých množinách:

maxx{0,1}n(i,j)xi+xj2xixj.\max_{x \in \{0,1\}^n} \sum_{(i,j)} x_i + x_j - 2x_ix_j.

To je ekvivalentně problém kvadratické neomezené binární optimalizace (QUBO) ve tvaru minxxTQx\min_x\, x^T Q x. Pomocí standardní substituce proměnných (xi(1Zi)/2x_i \to (1 - Z_i)/2) lze QUBO přepsat jako nákladový Hamiltonián, jehož základní stav kóduje optimální řešení. Obecně má tento Hamiltonián jak kvadratické, tak lineární členy:

HC=ijQijZiZj+ibiZi.H_C = \sum_{ij} Q_{ij} \, Z_i Z_j + \sum_i b_i \, Z_i.

Pro neváhovaný problém max-cut, který zde uvažujeme, lineární koeficienty vymizí (bi=0b_i = 0) a Qij=1Q_{ij} = 1 pro každou hranu, čímž zůstává jednodušší tvar HC=(i,j)EZiZjH_C = \sum_{(i,j) \in E} Z_i Z_j, který sestavíš v kódu níže. Obecnější tvar výše je to, co bys potřeboval k přizpůsobení tohoto pracovního postupu váhovaným grafům nebo jiným problémům vyjádřitelným pomocí QUBO.

Jak QAOA funguje

QAOA připravuje kandidátní řešení aplikací střídavých vrstev dvou operátorů na počáteční superpoziční stav Hn0H^{\otimes n}|0\rangle: nákladového operátoru eiγkHCe^{-i\gamma_k H_C} a mixujícího operátoru eiβkHme^{-i\beta_k H_m}. Úhly γk\gamma_k a βk\beta_k jsou optimalizovány v klasické zpětnovazební smyčce; kvantový počítač vyhodnocuje nákladovou funkci a klasický optimalizátor aktualizuje parametry až do konvergence. Tato iterační smyčka běží uvnitř Qiskit Runtime session, která udržuje kvantové zařízení rezervované napříč iteracemi pro nižší latenci.

Schéma Circuit s vrstvami QAOA

Pro hlubší zpracování teorie QAOA, včetně úplného odvození QUBO-na-Hamiltonián, viz modul kurzu QAOA.

V tomto tutoriálu nejprve vyřešíš max-cut na malém grafu s pěti uzly a poté škáluješ stejný pracovní postup na problém se 100 uzly v užitkové škále na reálném hardwaru. Poznámka k přístupu dle plánu: Tento tutoriál používá Qiskit Runtime sessions, které jsou dostupné pouze v rámci Premium Plánu. Pokud jsi na Open Plánu, nemůžeš tento tutoriál spustit tak, jak je napsán; místo toho budeš muset vyměnit Session za režim úloh (job mode) (tedy odeslat každou iteraci jako nezávislou úlohu místo zabalení optimalizační smyčky do with Session(...)). Pracovní postup stále poběží, ale každá iterace utrpí plnou latenci fronty místo opětovného použití rezervovaného zařízení. Více informací viz Přehled dostupných plánů.

Požadavky

Před zahájením tohoto tutoriálu se ujisti, že máš nainstalováno následující:

  • Qiskit SDK v2.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.

Nastavení

# Added by doQumentation — required packages for this notebook
!pip install -q matplotlib numpy qiskit qiskit-ibm-runtime rustworkx scipy
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

Příklad v malém měřítku

Tato sekce provádí každým krokem pracovního postupu QAOA na malé instanci max-cut s pěti uzly. Navzdory označení "malé měřítko" se tento příklad stále spouští na reálném hardwaru IBM Quantum — kód vybere Backend se 127 nebo více Qubity a spustí tam Circuit. Inicializuj svůj problém vytvořením grafu s n=5n=5 uzly.

n_small = 5

graph = rx.PyGraph()
graph.add_nodes_from(np.arange(0, n_small, 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)

Výstup předchozí buňky kódu

Krok 1: Mapování klasických vstupů na kvantový problém

Mapuj klasický graf na kvantové Circuit a operátory. Jak je popsáno v sekci Pozadí, pro neváhovaný max-cut se nákladový Hamiltonián redukuje na HC=(i,j)EZiZjH_C = \sum_{(i,j) \in E} Z_i Z_j a QAOA používá parametrizovaný ansatz Circuit k přípravě kandidátních základních stavů HCH_C.

Sestavení nákladového Hamiltoniánu

Převeď hrany grafu na Pauliho členy ZiZjZ_iZ_j k sestavení HCH_C (odvození viz Pozadí).

def build_max_cut_paulis(
graph: rx.PyGraph,
) -> list[tuple[str, list[int], float]]:
"""Convert graph edges to a list of ZZ Pauli terms.

The returned list is in the sparse format expected by
``SparsePauliOp.from_sparse_list``: each element is
``(pauli_string, qubit_indices, coefficient)``.
"""
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_small)
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])

Sestavení QAOA ansatz Circuit

Použij QAOAAnsatz k sestavení parametrizovaného QAOA Circuit z nákladového Hamiltoniánu. Zde používáme reps=2 (dvě vrstvy QAOA, čtyři parametry: β0,β1,γ0,γ1\beta_0, \beta_1, \gamma_0, \gamma_1).

circuit = QAOAAnsatz(cost_operator=cost_hamiltonian, reps=2)
circuit.measure_all()

circuit.draw("mpl")

Výstup předchozí buňky kódu

circuit.parameters
ParameterView([ParameterVectorElement(β[0]), ParameterVectorElement(β[1]), ParameterVectorElement(γ[0]), ParameterVectorElement(γ[1])])

Krok 2: Optimalizace problému pro spuštění na kvantovém hardwaru

Transpiluj abstraktní Circuit na instrukce nativní pro hardware. Tento krok řeší mapování Qubitů, dekompozici Gate, směrování a potlačení chyb. Více informací viz dokumentace k transpilaci.

service = QiskitRuntimeService()
backend = service.least_busy(
operational=True, simulator=False, min_num_qubits=127
)
print(backend)

# Create pass manager for transpilation. Level 3 is the most aggressive
# preset: slower to transpile, but produces shorter circuits that are
# more robust to hardware noise.
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('ibm_pittsburgh')>

Výstup předchozí buňky kódu

Krok 3: Spuštění pomocí Qiskit primitiv

Optimalizační smyčka QAOA běží uvnitř Qiskit Runtime session, aby zařízení zůstalo rezervované napříč iteracemi. Estimator vyhodnocuje HC\langle H_C \rangle v každém kroku a klasický optimalizátor (COBYLA) aktualizuje parametry až do konvergence.

Ilustrace zobrazující chování režimů Single job, Batch a Session runtime. Definuj počáteční parametry a spusť optimalizační smyčku:

# QAOA doesn't prescribe principled default angles — any bounded choice
# works as a warm start for problems this small. beta and gamma are
# periodic (beta in [0, pi] and gamma in [0, 2*pi] modulo the underlying
# Pauli-rotation periods), and pi/2 and pi are just midpoints of those
# ranges. For harder problems you would typically warm start from known
# good angles or transfer parameters from smaller instances.
initial_gamma = np.pi
initial_beta = np.pi / 2
init_params = [initial_beta, initial_beta, initial_gamma, initial_gamma]
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"
estimator.options.environment.job_tags = ["TUT_QAOA"]

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: -2.0402211719947774
x: [ 3.041e+00 1.212e+00 2.081e+00 4.471e+00]
nfev: 36
maxcv: 0.0

Optimalizátoru se podařilo snížit náklady a najít lepší parametry pro Circuit.

Hladce klesající křivka, která se ustálí na plató, je charakteristickým znakem konvergence. Zašuměná, nemonotónní křivka obvykle naznačuje, že něco výše v procesu vyžaduje pozornost; běžnými příčinami jsou příliš málo shotů na vyhodnocení (vysoký rozptyl estimatoru), špatné počáteční parametry nebo Circuit, jehož hloubka je ovládána hardwarovým šumem. COBYLA nepoužívá derivace a je poměrně robustní vůči mírnému šumu, ale když šum zaplaví skutečná zlepšení nákladů na krok, jeho lineárně-aproximační model už nedokáže rozlišit skutečný sestup od náhodného chvění a optimalizátor začne bloudit.

plt.figure(figsize=(12, 6))
plt.plot(objective_func_vals)
plt.xlabel("Iteration")
plt.ylabel("Cost")
plt.show()

Výstup předchozí buňky kódu

Přiřaď optimalizované parametry a navzorkuj finální rozdělení pomocí primitiva Sampler.

optimized_circuit = candidate_circuit.assign_parameters(result.x)
optimized_circuit.draw("mpl", fold=False, idle_wires=False)

Výstup předchozí buňky kódu

# 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"

sampler.options.environment.job_tags = ["TUT_QAOA"]

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)
{18: 0.039, 5: 0.0665, 20: 0.0973, 29: 0.0063, 9: 0.0899, 13: 0.0379, 2: 0.0047, 1: 0.0153, 11: 0.0932, 14: 0.0327, 12: 0.0314, 25: 0.0193, 21: 0.0398, 6: 0.0224, 4: 0.0197, 10: 0.0387, 3: 0.0181, 26: 0.07, 17: 0.0327, 19: 0.0332, 22: 0.0914, 24: 0.007, 0: 0.0033, 8: 0.0066, 30: 0.0158, 28: 0.0169, 27: 0.0222, 16: 0.0073, 7: 0.0057, 23: 0.0062, 15: 0.0054, 31: 0.0041}

Krok 4: Po-zpracování a vrácení výsledku v požadovaném klasickém formátu

Extrahuj nejpravděpodobnější bitřetězec z navzorkovaného rozdělení. Ten reprezentuje nejlepší řez, který QAOA našel.

# 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, 0, 1, 0, 1]
plt.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()

Výstup předchozí buňky kódu

Vizualizace nejlepšího řezu

Z optimálního bitřetězce pak můžeš tento řez vizualizovat 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)

Výstup předchozí buňky kódu

Nyní 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

Pro takto malý graf je skutečné optimum snadné získat hrubou silou, takže si můžeš výsledky ověřit porovnáním výsledku QAOA s přesnou odpovědí.

# Classical baseline: enumerate all 2**n_small bitstrings and take the best cut.
def brute_force_max_cut(graph: rx.PyGraph) -> tuple[int, list[int]]:
n = len(list(graph.nodes()))
best_cut = -1
best_x: list[int] = []
for i in range(2**n):
x = [(i >> k) & 1 for k in range(n)]
cut = evaluate_sample(x, graph)
if cut > best_cut:
best_cut = int(cut)
best_x = x
return best_cut, best_x

classical_best, classical_x = brute_force_max_cut(graph)
print(f"Classical optimum (brute force): {classical_best}")
print(f"QAOA cut value: {cut_value}")
Classical optimum (brute force): 5
QAOA cut value: 5

Příklad na hardwaru ve velkém měřítku

Na IBM Quantum Platform máš přístup k mnoha zařízením s více než 100 Qubity. Vyber si jedno, na kterém vyřešíš max-cut na váhovaném grafu se 100 uzly. Jedná se o problém v "užitkové škále". Pracovní postup sleduje stejné kroky jako výše, aplikované na mnohem větší graf.

Kompletní pracovní postup v užitkové škále

Všechny čtyři kroky jsou ukázány níže, aplikované na graf se 100 uzly. Struktura je stejná jako u průvodce v malém měřítku: mapuj, transpiluj, spusť, po-zpracuj — ale s větším problémem a rozděleno do čtyř buněk níže pro přehlednost.

# Precomputed parity lookup table: _PARITY[b] = +1 if popcount(b) is even, else -1.
# We use this to vectorize expectation-value evaluation across all Pauli terms.
_PARITY = np.array(
[-1 if bin(i).count("1") % 2 else 1 for i in range(256)],
dtype=np.complex128,
)

def evaluate_sparse_pauli(state: int, observable: SparsePauliOp) -> complex:
"""Expectation value of a SparsePauliOp on a single computational-basis state.

For a Z-only observable (which QAOA cost Hamiltonians are, after the
QUBO-to-Hamiltonian mapping), the eigenvalue of each Pauli term on a
computational-basis state is simply (-1)**popcount(z_mask AND state),
i.e., the parity of the bitwise-AND of the term's Z-support and the
measured bitstring.

This routine packs the Z-support of every Pauli term into bytes, ANDs
them against the measured state in a single vectorized op, and looks up
the parity in _PARITY. For a 100-qubit / ~hundreds-of-terms Hamiltonian
over 10_000 samples, this is dramatically faster than calling
SparsePauliOp.expectation_value per sample.
"""
packed_uint8 = np.packbits(observable.paulis.z, axis=1, bitorder="little")
state_bytes = np.frombuffer(
state.to_bytes(packed_uint8.shape[1], "little"), dtype=np.uint8
)
reduced = np.bitwise_xor.reduce(packed_uint8 & state_bytes, axis=1)
return np.sum(observable.coeffs * _PARITY[reduced])

def best_solution(samples, hamiltonian):
"""Return the sampled bitstring (as int) with the lowest Hamiltonian cost."""
min_cost = float("inf")
min_sol = None
for bit_str in samples.keys():
candidate_sol = int(bit_str)
fval = evaluate_sparse_pauli(candidate_sol, hamiltonian).real
if fval <= min_cost:
min_cost = fval
min_sol = candidate_sol
return min_sol

def _plot_cdf(objective_values: dict, ax, color):
x_vals = sorted(objective_values.keys(), reverse=True)
y_vals = np.cumsum([objective_values[x] for x in x_vals])
ax.plot(x_vals, y_vals, color=color)

def plot_cdf(dist, ax, title):
_plot_cdf(dist, ax, "C1")
ax.vlines(min(list(dist.keys())), 0, 1, "C1", linestyle="--")
ax.set_title(title)
ax.set_xlabel("Objective function value")
ax.set_ylabel("Cumulative distribution function")
ax.grid(alpha=0.3)

def samples_to_objective_values(samples, hamiltonian):
"""Convert the samples to values of the objective function."""
objective_values = defaultdict(float)
for bit_str, prob in samples.items():
candidate_sol = int(bit_str)
fval = evaluate_sparse_pauli(candidate_sol, hamiltonian).real
objective_values[fval] += prob
return objective_values

Krok 1: Sestav graf, nákladový Hamiltonián a ansatz.

# Step 1: build the 100-node graph, cost Hamiltonian, and QAOA ansatz.
n_large = 100
graph_100 = rx.PyGraph()
graph_100.add_nodes_from(np.arange(0, n_large, 1))
elist = []
for edge in backend.coupling_map:
if edge[0] < n_large and edge[1] < n_large:
elist.append((edge[0], edge[1], 1.0))
graph_100.add_edges_from(elist)

max_cut_paulis_100 = build_max_cut_paulis(graph_100)
cost_hamiltonian_100 = SparsePauliOp.from_sparse_list(
max_cut_paulis_100, n_large
)

circuit_100 = QAOAAnsatz(cost_operator=cost_hamiltonian_100, reps=1)
circuit_100.measure_all()

Krok 2: Transpiluj pro vybraný hardwarový Backend.

# Step 2: transpile for hardware.
pm = generate_preset_pass_manager(optimization_level=3, backend=backend)
candidate_circuit_100 = pm.run(circuit_100)

Krok 3: Spusť optimalizační smyčku QAOA uvnitř session a poté navzorkuj.

# Step 3: run the QAOA optimization loop on the device, then sample the
# final distribution with the optimized parameters.
initial_gamma = np.pi
initial_beta = np.pi / 2
init_params = [initial_beta, initial_gamma]

objective_func_vals = [] # Global variable
with Session(backend=backend) as 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"
estimator.options.environment.job_tags = ["TUT_QAOA"]

result = minimize(
cost_func_estimator,
init_params,
args=(candidate_circuit_100, cost_hamiltonian_100, estimator),
method="COBYLA",
)
print(result)

# Assign optimal parameters and sample the final distribution.
optimized_circuit_100 = candidate_circuit_100.assign_parameters(result.x)

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"

# Add a unique tag to the job execution
sampler.options.environment.job_tags = ["TUT_QAOA"]

pub = (optimized_circuit_100,)
job = sampler.run([pub], shots=int(1e4))

counts_int = job.result()[0].data.meas.get_int_counts()
shots = sum(counts_int.values())
final_distribution_100_int = {
key: val / shots for key, val in counts_int.items()
}
message: Return from COBYLA because the trust region radius reaches its lower bound.
success: True
status: 0
fun: -17.172689238986344
x: [ 2.574e+00 4.166e+00]
nfev: 28
maxcv: 0.0

Krok 4: Po-zpracuj navzorkované rozdělení k extrakci nejlepšího řezu.

# Step 4: find the best-cost sample and evaluate its cut value.
best_sol_100 = best_solution(final_distribution_100_int, cost_hamiltonian_100)
best_sol_bitstring_100 = to_bitstring(int(best_sol_100), len(graph_100))
best_sol_bitstring_100.reverse()

print("Result bitstring:", best_sol_bitstring_100)

cut_value_100 = evaluate_sample(best_sol_bitstring_100, graph_100)
print("The value of the cut is:", cut_value_100)
Result bitstring: [1, 1, 0, 1, 1, 0, 1, 1, 1, 0, 1, 0, 1, 1, 1, 0, 1, 1, 1, 1, 0, 0, 1, 0, 0, 1, 1, 0, 1, 0, 1, 0, 1, 1, 0, 1, 1, 0, 0, 0, 0, 0, 1, 1, 0, 1, 0, 0, 0, 1, 0, 0, 1, 0, 1, 0, 0, 1, 1, 1, 1, 1, 0, 1, 0, 1, 0, 0, 0, 1, 0, 1, 0, 1, 0, 0, 0, 0, 1, 0, 0, 1, 0, 1, 1, 1, 1, 0, 1, 0, 1, 0, 1, 1, 1, 0, 1, 1, 1, 0]
The value of the cut is: 156

Zkontroluj, že náklady minimalizované v optimalizační smyčce zkonvergovaly, a vizualizuj výsledky.

# Plot convergence
plt.figure(figsize=(12, 6))
plt.plot(objective_func_vals)
plt.xlabel("Iteration")
plt.ylabel("Cost")
plt.show()

# Visualize the cut
plot_result(graph_100, best_sol_bitstring_100)

# Plot cumulative distribution function
result_dist = samples_to_objective_values(
final_distribution_100_int, cost_hamiltonian_100
)
fig, ax = plt.subplots(1, 1, figsize=(8, 6))
plot_cdf(result_dist, ax, backend.name)

Výstup předchozí buňky kódu

Výstup předchozí buňky kódu

Výstup předchozí buňky kódu

Další kroky

Pokud tě tato práce zaujala, mohly by tě zajímat následující materiály:

Doporučení