Quantenalgorithmen: Variationale Quantenalgorithmen

1. Einführung

Dieses Tutorial gibt einen Überblick über einen hybriden quantenklassischen Algorithmus, wobei der Schwerpunkt auf dem Variations-Quanten-Eigensolver (VQE) und dem Quanten-Approximationsalgorithmus (QAOA) liegt. Das Hauptziel dieser Algorithmen ist die Lösung von Optimierungsproblemen durch den Einsatz von Quantenschaltungen mit parametrisierten Quantengattern.

Trotz der Fortschritte in der Quanteninformatik ist es aufgrund des Rauschens in aktuellen Quantengeräten schwierig, aussagekräftige Ergebnisse aus tiefen Quantenschaltungen zu gewinnen. Um diese Herausforderung zu meistern, verfolgen VQE und QAOA einen hybriden quantenklassischen Ansatz, bei dem relativ kurze Quantenschaltungen mit Hilfe von Quantenberechnungen iterativ ausgeführt und die Parameter der parametrisierten Zielquantenschaltungen mit Hilfe klassischer Berechnungen optimiert werden.

QAOA hat das Potenzial, optimale Lösungen für die Zielprobleme auf der Ebene der Versorgungsunternehmen zu liefern, dank der Anwendung verschiedener Techniken zur Fehlerbegrenzung und -unterdrückung. VQE hat viele Anwendungen (wie die Quantenchemie), bei denen es weniger skalierbar ist. Es gibt jedoch eine Reihe von eigenwertbezogenen Ansätzen, die die VQE ergänzen und erweitern, darunter die Krylov-Unterraumdiagonalisierung und die Sampling-basierte Quantendiagonalisierung (SQD). Das Verständnis der VQE ist ein wichtiger erster Schritt, um die breite Palette der entstandenen klassisch-quantischen Hybridalgorithmen zu verstehen.

Dieses Modul beschreibt die grundlegenden Konzepte und die Implementierung von VQE und QAOA. In weiteren Tutorien werden fortgeschrittene Themen und Techniken zur Skalierung dieser Algorithmen behandelt.

Sie benötigen die folgende Bibliothek in Ihrer Umgebung, um dieses Notizbuch auszuführen. Wenn Sie es noch nicht installiert haben, können Sie es installieren, indem Sie den Kommentar entfernen und die folgende Zelle ausführen.

# % pip install 'qiskit[visualization]' qiskit-ibm-runtime

2. Berechnung des minimalen Eigenwerts eines einfachen Hamilton-Operators

Wir beginnen mit der Anwendung von VQE auf einen sehr einfachen Fall, um zu sehen, wie es funktioniert. Wir werden den minimalen Eigenwert der Pauli-Matrix $Z$ mit VQE berechnen. Wir beginnen mit dem Import einiger allgemeiner Pakete.

import numpy as np
from qiskit.circuit import ParameterVector, QuantumCircuit
from qiskit.primitives import StatevectorEstimator, StatevectorSampler
from qiskit.quantum_info import SparsePauliOp
from scipy.optimize import minimize

Wir definieren nun den interessierenden Operator und betrachten ihn in Matrixform.

op = SparsePauliOp("Z")
op.to_matrix()

Output:

array([[ 1.+0.j,  0.+0.j],
       [ 0.+0.j, -1.+0.j]])

Es ist einfach, die Eigenwerte klassisch zu ermitteln, so dass wir unsere Arbeit überprüfen können. Dies könnte sich als schwierig erweisen, wenn wir uns in Richtung Nutzwert entwickeln. Hier verwenden wir numpy.

# compute eigenvalues with numpy
result = np.linalg.eigh(op.to_matrix())
print("Eigenvalues:", result.eigenvalues)

Output:

Eigenvalues: [-1.  1.]

Um Eigenwerte mit Hilfe eines Quantenvariationsalgorithmus zu erhalten, konstruieren wir einen Schaltkreis mit Gattern, die Variationsparameter annehmen:

# define a variational form
param = ParameterVector("a", 3)
qc = QuantumCircuit(1, 1)
qc.u(param[0], param[1], param[2], 0)
qc_estimator = qc.copy()
qc.measure(0, 0)
qc.draw("mpl")

Output:

Wenn wir den Erwartungswert eines Operators (wie $Z$ ) schätzen wollen, sollten wir Estimator verwenden. Wenn wir uns die Zustände des Systems ansehen wollen, verwenden wir Sampler.

sampler = StatevectorSampler()
estimator = StatevectorEstimator()

Wir können die Anzahl der Bitstrings 0 und 1 mit zufälligen Parameterwerten [1, 2, 3] mit Sampler berechnen.

# compute counts of bitstrings with random parameter values by Sampler
result = sampler.run([(qc, [1, 2, 3])]).result()
counts = result[0].data.c.get_counts()
counts

Output:

{'0': 783, '1': 241}

Wir wissen, dass wir den Erwartungswert von Z durch $\langle Z \rangle = p_0 - p_1$ mit Wahrscheinlichkeiten $\{0: p_0, 1: p_1\}$ berechnen können.

# compute the expectation value of Z based on the counts
(counts.get("0", 0) - counts.get("1", 0)) / sum(counts.values())

Output:

0.529296875

Diese Schaltung funktionierte, aber die gewählten Parameterwerte entsprachen nicht einem Zustand mit sehr niedriger Energie (oder niedrigem Eigenwert). Der erhaltene Eigenwert ist um einiges höher als das Minimum. Das Ergebnis ist ähnlich, wenn der Schätzer verwendet wird.

Beachten Sie, dass Estimator Quantenschaltungen ohne Messungen verwendet.

result = estimator.run([(qc_estimator, op, [1, 2, 3])]).result()
result[0].data.evs

Output:

array(0.54030231)

Wir müssen die Parameter durchsuchen und diejenigen finden, die den niedrigsten Eigenwert ergeben. Wir erstellen eine Funktion, die die Parameterwerte der Variationsform aufnimmt und den Erwartungswert $\langle Z \rangle$ zurückgibt.

# define a cost function to look for the minimum eigenvalue of Z
def cost(x):
    result = sampler.run([(qc, x)]).result()
    counts = result[0].data.c.get_counts()
    expval = (counts.get("0", 0) - counts.get("1", 0)) / sum(counts.values())
    # the following line shows the trajectory of the optimization
    print(expval, counts)
    return expval

Wenden wir die Funktion SciPy's minimize an, um den minimalen Eigenwert von Z zu finden.

# minimize the cost function with scipy's minimize
min_result = minimize(cost, [0, 0, 0], method="COBYLA", tol=1e-8)
min_result

Output:

1.0 {'0': 1024}
0.494140625 {'0': 765, '1': 259}
0.466796875 {'0': 751, '1': 273}
0.564453125 {'0': 801, '1': 223}
-0.4296875 {'1': 732, '0': 292}
-0.984375 {'1': 1016, '0': 8}
-0.8984375 {'1': 972, '0': 52}
-0.990234375 {'1': 1019, '0': 5}
-0.892578125 {'1': 969, '0': 55}
-0.986328125 {'1': 1017, '0': 7}
-0.861328125 {'1': 953, '0': 71}
-1.0 {'1': 1024}
-0.982421875 {'1': 1015, '0': 9}
-0.99609375 {'1': 1022, '0': 2}
-0.986328125 {'1': 1017, '0': 7}
-1.0 {'1': 1024}
-0.990234375 {'1': 1019, '0': 5}
-0.998046875 {'1': 1023, '0': 1}
-0.99609375 {'1': 1022, '0': 2}
-1.0 {'1': 1024}
-1.0 {'1': 1024}
-1.0 {'1': 1024}
-1.0 {'1': 1024}
-0.998046875 {'1': 1023, '0': 1}
-1.0 {'1': 1024}
-1.0 {'1': 1024}
-0.998046875 {'1': 1023, '0': 1}
-0.998046875 {'1': 1023, '0': 1}
-0.998046875 {'1': 1023, '0': 1}
-1.0 {'1': 1024}
-0.99609375 {'1': 1022, '0': 2}
-1.0 {'1': 1024}
-0.99609375 {'1': 1022, '0': 2}
-0.998046875 {'1': 1023, '0': 1}
-0.998046875 {'1': 1023, '0': 1}
-0.99609375 {'1': 1022, '0': 2}
-0.998046875 {'1': 1023, '0': 1}
-1.0 {'1': 1024}
-0.998046875 {'1': 1023, '0': 1}
-0.998046875 {'1': 1023, '0': 1}
-0.99609375 {'1': 1022, '0': 2}
-1.0 {'1': 1024}
-0.998046875 {'1': 1023, '0': 1}
-1.0 {'1': 1024}
-0.998046875 {'1': 1023, '0': 1}
-0.998046875 {'1': 1023, '0': 1}
-1.0 {'1': 1024}
-0.998046875 {'1': 1023, '0': 1}
-0.998046875 {'1': 1023, '0': 1}
-1.0 {'1': 1024}
-1.0 {'1': 1024}
-1.0 {'1': 1024}
-1.0 {'1': 1024}
-1.0 {'1': 1024}
-1.0 {'1': 1024}
-0.998046875 {'1': 1023, '0': 1}
-0.994140625 {'1': 1021, '0': 3}
-1.0 {'1': 1024}
-1.0 {'1': 1024}
-1.0 {'1': 1024}
-1.0 {'1': 1024}
-1.0 {'1': 1024}
-1.0 {'1': 1024}

 message: Optimization terminated successfully.
 success: True
  status: 1
     fun: -1.0
       x: [ 3.182e+00  1.338e+00  1.664e-01]
    nfev: 63
   maxcv: 0.0

# check counts of bitstrings with the optimal parameters
result = sampler.run([(qc, min_result.x)]).result()
result[0].data.c.get_counts()

Output:

{'0': 1, '1': 1023}

2.1 Übung

Berechnen Sie den minimalen Eigenwert von $Z \otimes Z$ mit VQE.

z2 = SparsePauliOp("ZZ")
print(z2)
print(z2.to_matrix())

Output:

SparsePauliOp(['ZZ'],
              coeffs=[1.+0.j])
[[ 1.+0.j  0.+0.j  0.+0.j  0.+0.j]
 [ 0.+0.j -1.+0.j  0.+0.j  0.+0.j]
 [ 0.+0.j  0.+0.j -1.+0.j  0.+0.j]
 [ 0.+0.j  0.+0.j  0.+0.j  1.+0.j]]

# compute eigenvalues with numpy

# define a variational form
# qc = ...

# compute counts of bitstrings with a random parameter values by Sampler
# result = sampler.run(...)
# result

# compute the expectation value of ZZ based on the counts

# verify the expectation value of ZZ with Estimator

# define a cost function to look for the minimum eigenvalue of ZZ
# def cost(x):
#    expval = ...
#    return expval

# minimize the cost function with scipy's minimize
# min_result = minimize(cost, [...], method="COBYLA", tol=1e-8)
# min_result

# check counts of bitstrings with the optimal parameter values
# result = sampler.run(qc, min_result.x).result()
# result

Lösungen zur Übung

Wir definieren den interessierenden Operator und betrachten ihn in Matrixform.

z2 = SparsePauliOp("ZZ")
print(z2)
print(z2.to_matrix())

Output:

SparsePauliOp(['ZZ'],
              coeffs=[1.+0.j])
[[ 1.+0.j  0.+0.j  0.+0.j  0.+0.j]
 [ 0.+0.j -1.+0.j  0.+0.j  0.+0.j]
 [ 0.+0.j  0.+0.j -1.+0.j  0.+0.j]
 [ 0.+0.j  0.+0.j  0.+0.j  1.+0.j]]

Um Eigenwerte mit Hilfe eines Quantenvariationsalgorithmus zu erhalten, konstruieren wir einen Schaltkreis mit Gattern, die Variationsparameter annehmen:

# define a variational form
param = ParameterVector("a", 6)
qc = QuantumCircuit(2, 2)
qc.u(param[0], param[1], param[2], 0)
qc.u(param[3], param[4], param[5], 1)
qc_estimator = qc.copy()
qc.measure([0, 1], [0, 1])
qc.draw("mpl")

Output:

Wenn wir den Erwartungswert eines Operators (z. B. $Z \otimes Z$ ) schätzen wollen, würden wir Estimator verwenden. Wenn wir uns die Zustände des Systems ansehen wollen, verwenden wir Sampler.

sampler = StatevectorSampler()
estimator = StatevectorEstimator()

# compute counts of bitstrings with random parameter values by Sampler
result = sampler.run([(qc, [1, 2, 3, 4, 5, 6])]).result()
counts = result[0].data.c.get_counts()
counts

Output:

{'10': 661, '11': 203, '01': 47, '00': 113}

# compute the expectation value of ZZ based on the counts
(
    counts.get("00", 0)
    - counts.get("01", 0)
    - counts.get("10", 0)
    + counts.get("11", 0)
) / sum(counts.values())

Output:

-0.3828125

Diese Schaltung funktionierte, aber die gewählten Parameterwerte entsprachen nicht einem Zustand mit sehr niedriger Energie (oder niedrigem Eigenwert). Der erhaltene Eigenwert ist um einiges höher als das Minimum. Das Ergebnis ist ähnlich, wenn der Schätzer verwendet wird.

# verify the expectation value of ZZ with Estimator
result = estimator.run([(qc_estimator, z2, [1, 2, 3, 4, 5, 6])]).result()
result[0].data.evs

Output:

array(-0.35316516)

Wir müssen die Parameter durchsuchen und diejenigen finden, die den niedrigsten Eigenwert ergeben.

# define a cost function to look for the minimum eigenvalue of ZZ
def cost(x):
    result = sampler.run([(qc, x)]).result()
    counts = result[0].data.c.get_counts()
    expval = (
        counts.get("00", 0)
        - counts.get("01", 0)
        - counts.get("10", 0)
        + counts.get("11", 0)
    ) / sum(counts.values())
    print(expval, counts)
    return expval

# minimize the cost function with scipy's minimize
min_result = minimize(cost, [0, 0, 0, 0, 0, 0], method="COBYLA", tol=1e-8)
min_result

Output:

1.0 {'00': 1024}
0.578125 {'00': 808, '01': 216}
0.5234375 {'00': 780, '01': 244}
0.548828125 {'00': 793, '01': 231}
0.3515625 {'00': 637, '10': 164, '11': 55, '01': 168}
0.3359375 {'00': 638, '11': 46, '10': 174, '01': 166}
0.283203125 {'00': 602, '10': 181, '01': 186, '11': 55}
-0.087890625 {'01': 414, '00': 184, '10': 143, '11': 283}
0.236328125 {'10': 27, '11': 623, '01': 364, '00': 10}
-0.0625 {'11': 261, '01': 403, '00': 219, '10': 141}
0.248046875 {'01': 366, '11': 628, '00': 11, '10': 19}
-0.0625 {'10': 145, '11': 254, '01': 399, '00': 226}
0.228515625 {'01': 373, '11': 609, '00': 20, '10': 22}
0.0546875 {'11': 376, '10': 273, '01': 211, '00': 164}
-0.447265625 {'01': 731, '10': 10, '11': 267, '00': 16}
-0.71484375 {'01': 871, '11': 99, '00': 47, '10': 7}
-0.46484375 {'01': 741, '00': 253, '10': 9, '11': 21}
-0.87890625 {'01': 962, '00': 39, '11': 23}
-0.640625 {'00': 176, '01': 837, '11': 8, '10': 3}
-0.88671875 {'01': 966, '00': 41, '11': 17}
-0.994140625 {'01': 1021, '11': 3}
-0.91796875 {'01': 982, '11': 35, '00': 7}
-0.994140625 {'01': 1021, '11': 2, '00': 1}
-0.939453125 {'01': 993, '00': 31}
-0.990234375 {'01': 1019, '11': 5}
-0.90234375 {'01': 974, '00': 21, '11': 29}
-0.98046875 {'01': 1014, '11': 10}
-0.994140625 {'01': 1021, '00': 3}
-0.990234375 {'01': 1019, '11': 4, '00': 1}
-0.98828125 {'01': 1018, '11': 6}
-0.990234375 {'01': 1019, '11': 4, '00': 1}
-0.994140625 {'01': 1021, '11': 2, '00': 1}
-0.99609375 {'01': 1022, '11': 2}
-0.998046875 {'01': 1023, '00': 1}
-0.99609375 {'01': 1022, '00': 2}
-1.0 {'01': 1024}
-1.0 {'01': 1024}
-1.0 {'01': 1024}
-0.998046875 {'01': 1023, '11': 1}
-1.0 {'01': 1024}
-1.0 {'01': 1024}
-1.0 {'01': 1024}
-1.0 {'01': 1024}
-1.0 {'01': 1024}
-1.0 {'01': 1024}
-1.0 {'01': 1024}
-0.998046875 {'01': 1023, '00': 1}
-0.998046875 {'01': 1023, '11': 1}
-0.998046875 {'01': 1023, '00': 1}
-1.0 {'01': 1024}
-1.0 {'01': 1024}
-1.0 {'01': 1024}
-1.0 {'01': 1024}
-1.0 {'01': 1024}
-1.0 {'01': 1024}
-0.998046875 {'01': 1023, '11': 1}
-0.998046875 {'01': 1023, '11': 1}
-1.0 {'01': 1024}
-1.0 {'01': 1024}
-0.998046875 {'01': 1023, '11': 1}
-0.998046875 {'01': 1023, '11': 1}
-0.998046875 {'01': 1023, '00': 1}
-1.0 {'01': 1024}
-1.0 {'01': 1024}
-0.998046875 {'01': 1023, '00': 1}
-1.0 {'01': 1024}
-1.0 {'01': 1024}
-1.0 {'01': 1024}
-1.0 {'01': 1024}
-0.998046875 {'01': 1023, '11': 1}
-0.998046875 {'01': 1023, '11': 1}
-1.0 {'01': 1024}
-0.998046875 {'01': 1023, '11': 1}
-1.0 {'01': 1024}
-1.0 {'01': 1024}
-1.0 {'01': 1024}
-0.998046875 {'01': 1023, '11': 1}
-0.998046875 {'01': 1023, '11': 1}
-1.0 {'01': 1024}
-1.0 {'01': 1024}
-0.998046875 {'01': 1023, '11': 1}
-0.998046875 {'01': 1023, '11': 1}
-0.998046875 {'01': 1023, '00': 1}
-1.0 {'01': 1024}
-1.0 {'01': 1024}
-1.0 {'01': 1024}
-0.998046875 {'01': 1023, '11': 1}
-1.0 {'01': 1024}
-0.99609375 {'01': 1022, '00': 1, '11': 1}
-0.998046875 {'01': 1023, '11': 1}
-0.998046875 {'01': 1023, '00': 1}
-0.998046875 {'01': 1023, '11': 1}
-1.0 {'01': 1024}
-0.99609375 {'01': 1022, '11': 1, '00': 1}
-1.0 {'01': 1024}
-0.998046875 {'01': 1023, '00': 1}
-0.994140625 {'01': 1021, '00': 3}
-0.998046875 {'01': 1023, '00': 1}
-0.99609375 {'01': 1022, '11': 2}
-1.0 {'01': 1024}
-1.0 {'01': 1024}
-0.998046875 {'01': 1023, '11': 1}
-1.0 {'01': 1024}
-1.0 {'01': 1024}
-1.0 {'01': 1024}
-1.0 {'01': 1024}
-1.0 {'01': 1024}
-1.0 {'01': 1024}
-0.998046875 {'01': 1023, '11': 1}
-1.0 {'01': 1024}
-1.0 {'01': 1024}
-1.0 {'01': 1024}
-0.998046875 {'01': 1023, '11': 1}
-0.998046875 {'01': 1023, '11': 1}
-1.0 {'01': 1024}
-0.998046875 {'01': 1023, '00': 1}
-1.0 {'01': 1024}
-1.0 {'01': 1024}
-1.0 {'01': 1024}
-1.0 {'01': 1024}
-1.0 {'01': 1024}
-1.0 {'01': 1024}
-0.998046875 {'01': 1023, '11': 1}
-0.998046875 {'01': 1023, '11': 1}
-0.998046875 {'01': 1023, '11': 1}
-0.99609375 {'01': 1022, '11': 2}
-1.0 {'01': 1024}
-0.998046875 {'01': 1023, '11': 1}

 message: Optimization terminated successfully.
 success: True
  status: 1
     fun: -0.998046875
       x: [ 3.167e+00  6.940e-01  1.033e+00 -2.894e-02  8.933e-01
            1.885e+00]
    nfev: 128
   maxcv: 0.0

 message: Optimization terminated successfully.
 success: True
  status: 1
     fun: -0.99609375
       x: [ 3.098e+00 -5.402e-01  1.091e+00 -1.004e-02  3.615e-01
            6.913e-01]
    nfev: 115
   maxcv: 0.0

Wir haben einen Eigenwert erhalten, der sehr nahe am Minimum liegt, das uns von numpy vorgegeben wurde.

# check counts of bitstrings with the optimal parameters
result = sampler.run([(qc, min_result.x)]).result()
result[0].data.c.get_counts()

Output:

{'01': 1024}

3. Quantenoptimierung mit Qiskit-Mustern

In dieser Anleitung lernen wir etwas über Qiskit-Muster und Quantenapproximation. Ein Qiskit-Muster ist eine intuitive, wiederholbare Reihe von Schritten zur Implementierung eines Quantencomputer-Workflows:

Wir werden die Muster auf den Kontext der kombinatorischen Optimierung anwenden und zeigen, wie sich das Max-Cut- Problem mithilfe des Quantum Approximate Optimization Algorithm (QAOA), einer hybriden (quanten-klassischen) iterativen Methode, lösen lässt.

Beachten Sie, dass dieser QAOA-Teil auf "Teil 1: QAOA im kleinen Maßstab" des Tutorials über den Algorithmus für die Quantenapproximation basiert. In der Anleitung erfahren Sie, wie Sie es vergrößern können.

3.1 (Kleines) Qiskit-Muster zur Optimierung

In diesem Abschnitt werden anhand eines kleinen Max-Cut-Problems die Schritte erläutert, die zur Lösung eines Optimierungsproblems mit einem Quantencomputer erforderlich sind.

Das Max-Cut-Problem ist ein schwer zu lösendes Optimierungsproblem (genauer gesagt handelt es sich um ein NP-schweres Problem), das in den Bereichen Clusteranalyse, Netzwerkwissenschaft und statistische Physik vielfältige Anwendungsmöglichkeiten findet. In diesem Tutorial wird ein Graph betrachtet, der aus durch Kanten verbundenen Knoten besteht, und es geht darum, die Knoten durch das „Durchtrennen“ von Kanten in zwei Mengen zu unterteilen, sodass die Anzahl der durchtrennten Kanten maximiert wird.

Um vor der Umsetzung dieses Problems in einen Quantenalgorithmus etwas Hintergrundwissen zu vermitteln, lässt sich besser nachvollziehen, wie das Max-Cut-Problem zu einem klassischen kombinatorischen Optimierungsproblem wird, wenn man zunächst die Minimierung einer Funktion betrachtet: $f(x)$

wobei die Eingabe $x$ ein Vektor ist, dessen Komponenten den einzelnen Knoten eines Graphen entsprechen. Dann beschränken Sie jede dieser Komponenten entweder auf $0$ oder $1$ (was bedeutet, dass sie in den Schnitt einbezogen oder nicht einbezogen werden). In diesem kleinen Beispielfall wird ein Graph mit $n=5$ Knoten verwendet.

Sie könnten eine Funktion für ein Knotenpaar $i,j$ schreiben, die anzeigt, ob die entsprechende Kante $(i,j)$ im Schnitt liegt. Zum Beispiel ist die Funktion $x_i + x_j - 2 x_i x_j$ nur dann 1, wenn entweder $x_i$ oder $x_j$ den Wert 1 haben (was bedeutet, dass die Kante im Schnitt liegt), und ansonsten Null. Das Problem der Maximierung der Kanten im Schnitt kann wie folgt formuliert werden

die als Minimierung der folgenden Form umgeschrieben werden kann

Das Minimum von $f(x)$ ist in diesem Fall, wenn die Anzahl der vom Schnitt durchquerten Kanten maximal ist. Wie Sie sehen können, gibt es noch nichts, was mit Quantencomputing zu tun hat. Sie müssen dieses Problem in etwas umformulieren, das ein Quantencomputer verstehen kann.

Initialisieren Sie Ihr Problem, indem Sie einen Graphen mit $n=5$ Knoten erstellen.

import matplotlib
import matplotlib.pyplot as plt
import numpy as np
import rustworkx as rx
from rustworkx.visualization import mpl_draw

n = 5

graph = rx.PyGraph()
graph.add_nodes_from(range(1, n + 1))
edge_list = [
    (0, 1, 1.0),
    (0, 2, 1.0),
    (1, 2, 1.0),
    (1, 3, 1.0),
    (2, 4, 1.0),
    (3, 4, 1.0),
]
graph.add_edges_from(edge_list)
pos = rx.spring_layout(graph, seed=2)
mpl_draw(graph, node_size=600, pos=pos, with_labels=True, labels=str)

Output:

3.2 Schritt 1. Klassische Eingaben einem Quantenproblem zuordnen

Der erste Schritt des Musters besteht darin, das klassische Problem (Graph) in Quantenschaltungen und -operatoren abzubilden. Dazu sind im Wesentlichen drei Schritte erforderlich:

1. Verwenden Sie eine Reihe von mathematischen Umformulierungen, um dieses Problem unter Verwendung der Notation für quadratische, uneingeschränkte binäre Optimierungsprobleme (QUBO) darzustellen.
2. Schreiben Sie das Optimierungsproblem in einen Hamiltonian um, bei dem der Grundzustand der Lösung entspricht, die die Kostenfunktion minimiert.
3. Erstellen Sie einen Quantenschaltkreis, der den Grundzustand dieses Hamiltonianers durch einen Prozess ähnlich dem Quantenglühen vorbereitet.

Hinweis: Bei der QAOA-Methode benötigen Sie letztendlich einen Operator (Hamiltonian ), der die Kostenfunktion unseres hybriden Algorithmus darstellt, sowie einen parametrisierten Schaltkreis (Ansatz ), der Quantenzustände mit Lösungsansätzen für das Problem darstellt. Sie können aus diesen Kandidatenzuständen eine Stichprobe ziehen und sie dann anhand der Kostenfunktion bewerten.

Grafik → Optimierungsproblem

Der erste Schritt der Zuordnung ist eine Änderung der Notation. Im Folgenden wird das Problem in QUBO-Notation ausgedrückt:

dabei ist $Q$ eine $n\times n$ Matrix aus reellen Zahlen, $n$ entspricht der Anzahl der Knoten in Ihrem Graphen, $x$ ist der oben eingeführte Vektor der binären Variablen und $x^T$ gibt die Transponierung des Vektors $x$ an.

Problem name: maxcut

Minimize
  2*x_1*x_2 + 2*x_1*x_3 + 2*x_2*x_3 + 2*x_2*x_4 + 2*x_3*x_5 + 2*x_4*x_5 - 2*x_1
  - 3*x_2 - 3*x_3 - 2*x_4 - 2*x_5

Subject to
  No constraints

  Binary variables (5)
    x_1 x_2 x_3 x_4 x_5

Optimierungsproblem → Hamilton-Funktion

Sie können dann das QUBO-Problem als Hamiltonian umformulieren (hier eine Matrix, die die Energie eines Systems darstellt):

Reformulierungsschritte vom QAOA-Problem zum Hamiltonian

Um zu demonstrieren, wie das QAOA-Problem auf diese Weise umgeschrieben werden kann, ersetzen Sie zunächst die binären Variablen $x_i$ durch einen neuen Satz von Variablen $z_i\in\{-1, 1\}$ über

Hier sieht man, dass, wenn $x_i$ gleich $0$ ist, $z_i$ gleich $1$ sein muss. Wenn die $x_i$ 's für die $z_i$ 's in der Optimierung Problem ( $x^TQx$ ) ersetzt werden, kann eine äquivalente Formulierung erhalten werden.

Wenn wir nun $b_i=-\sum_{j}(Q_{ij}+Q_{ji})$ definieren, den Vorfaktor und den konstanten Term $n^2$ entfernen, erhalten wir die beiden äquivalenten Formulierungen desselben Optimierungsproblems.

Hier hängt $b$ von $Q$ ab. Um $z^TQz + b^Tz$ zu erhalten, haben wir den Faktor 1/4 und einen konstanten Offset von $n^2$ weggelassen, die bei der Optimierung keine Rolle spielen.

Um nun eine Quantenformulierung des Problems zu erhalten, muss man die $z_i$ -Variablen zu einer Pauli-Matrix $Z$ befördern, beispielsweise zu einer $2\times 2$ -Matrix der Form

Wenn Sie diese Matrizen in das obige Optimierungsproblem einsetzen, erhalten Sie die folgende Hamiltonsche Funktion

Erinnern Sie sich auch daran, dass die $Z$ Matrizen in den Rechenraum des Quantencomputers eingebettet sind, d.h. in einen Hilbert-Raum der Größe $2^n\times 2^n$. Daher sollten Sie Begriffe wie $Z_iZ_j$ als das Tensorprodukt $Z_i\otimes Z_j$ verstehen, das in den Hilbert-Raum $2^n\times 2^n$ eingebettet ist. Bei einem Problem mit fünf Entscheidungsvariablen ist beispielsweise der Begriff $Z_1Z_3$ gleichbedeutend mit $I\otimes Z_3\otimes I\otimes Z_1\otimes I$, wobei $I$ die Identitätsmatrix $2\times 2$ ist.

Dieser Hamiltonian wird kostenfunktion Hamiltonian genannt. Er hat die Eigenschaft, dass sein Grundzustand der Lösung von minimiert die Kostenfunktion $f(x)$ entspricht. Um Ihr Optimierungsproblem zu lösen, müssen Sie nun den Grundzustand von $H_C$ (oder einen Zustand mit einer hohen Überlappung) auf dem Quantencomputer vorbereiten. Dann wird eine Stichprobe aus diesem Zustand mit hoher Wahrscheinlichkeit die Lösung von $\min~f(x)$ ergeben.

def build_max_cut_operator(graph: rx.PyGraph) -> tuple[SparsePauliOp, float]:
    sp_list = []
    constant = 0
    for s, t in graph.edge_list():
        w = graph.get_edge_data(s, t)
        sp_list.append(("ZZ", [s, t], w / 2))
        constant -= 1 / 2
    return SparsePauliOp.from_sparse_list(
        sp_list, num_qubits=graph.num_nodes()
    ), constant

cost_hamiltonian, constant = build_max_cut_operator(graph)
print("Cost Function Hamiltonian:", cost_hamiltonian)
print("Constant:", constant)

Output:

Cost Function Hamiltonian: SparsePauliOp(['IIIZZ', 'IIZIZ', 'IIZZI', 'IZIZI', 'ZIZII', 'ZZIII'],
              coeffs=[0.5+0.j, 0.5+0.j, 0.5+0.j, 0.5+0.j, 0.5+0.j, 0.5+0.j])
Constant: -3.0

Hamiltonian → Quantenschaltung

Der Hamiltonian $H_C$ enthält die Quantendefinition Ihres Problems. Jetzt können Sie einen Quantenschaltkreis erstellen, mit dessen Hilfe Sie gute Lösungen vom Quantencomputer abfragen können. Die QAOA ist von der Quantenglühung inspiriert und wendet alternierende Schichten von Operatoren in der Quantenschaltung an.

Die allgemeine Idee ist, im Grundzustand eines bekannten Systems zu beginnen ( $H^{\otimes n}|0\rangle$ ) und dann das System in den Grundzustand des Kostenoperators zu lenken, an dem man interessiert ist. Dies geschieht durch die Anwendung der Operatoren $\exp\{-i\gamma_k H_C\}$ und $\exp\{-i\beta_k H_m\}$ mit den Winkeln $\gamma_1,...,\gamma_p$ und $\beta_1,...,\beta_p~$.

Der von Ihnen erzeugte Quantenschaltkreis wird durch $\gamma_i$ und $\beta_i$** parametrisiert**, so dass Sie verschiedene Werte von $\gamma_i$ und $\beta_i$ ausprobieren und aus dem resultierenden Zustand eine Probe entnehmen können.

In diesem Fall werden wir ein Beispiel mit einer QAOA-Ebene ausprobieren, die zwei Parameter enthält: $\gamma_1$ und $\beta_1$.

from qiskit.circuit.library import QAOAAnsatz

circuit = QAOAAnsatz(cost_operator=cost_hamiltonian, reps=1)
circuit.measure_all()
circuit.draw("mpl")

Output:

circuit.decompose(reps=3).draw("mpl", fold=-1)

Output:

circuit.parameters

Output:

ParameterView([ParameterVectorElement(β[0]), ParameterVectorElement(γ[0])])

3.3 Schritt 2. Optimierung von Schaltkreisen für die Ausführung auf Quantenhardware

Die obige Schaltung enthält eine Reihe von Abstraktionen, die nützlich sind, um über Quantenalgorithmen nachzudenken, aber nicht auf der Hardware ausgeführt werden können. Um auf einer QPU ausgeführt werden zu können, muss die Schaltung eine Reihe von Operationen durchlaufen, die den Schritt der Transpilierung oder der Schaltungsoptimierung des Musters bilden.

Die Qiskit-Bibliothek bietet eine Reihe von Transpilationspässen, die für eine Vielzahl von Schaltkreistransformationen geeignet sind. Sie müssen sicherstellen, dass Ihre Schaltung für Ihren Zweck optimiert ist.

Die Umsetzung kann mehrere Schritte umfassen, wie zum Beispiel:

• Anfängliche Zuordnung der Qubits im Schaltkreis (z. B. Entscheidungsvariablen) zu den physischen Qubits auf dem Gerät.
• Abrollen der Befehle im Quantenschaltkreis zu den hardware-nativen Befehlen, die das Backend versteht.
• Weiterleitung aller Qubits im Schaltkreis, die miteinander interagieren, an benachbarte physikalische Qubits.
• Fehlerunterdrückung durch Hinzufügen von Ein-Qubit-Gattern zur Rauschunterdrückung mit dynamischer Entkopplung.

Weitere Informationen zur Transpilierung finden Sie in unserer Dokumentation.

Der folgende Code wandelt und optimiert die abstrakte Schaltung in ein Format, das für die Ausführung auf einem der Geräte bereit ist, die über die Cloud mit dem Qiskit IBM® Runtime Service zugänglich sind.

Beachten Sie, dass Sie Ihre Programme im „lokalen Testmodus“ lokal testen können, bevor Sie sie an echte Quantencomputer senden. Weitere Informationen zum lokalen Testmodus finden Sie in der Dokumentation.

from qiskit_ibm_runtime import QiskitRuntimeService
from qiskit.transpiler.preset_passmanagers import generate_preset_pass_manager

# Use a quantum device
service = QiskitRuntimeService()
backend = service.least_busy(min_num_qubits=127)
# backend = service.backend("ibm_kingston")

# You can test your programs locally with a fake backend (local testing mode)
# backend = FakeBrisbane()

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)

Output:

  service = QiskitRuntimeService(channel="ibm_quantum_platform")

<IBMBackend('ibm_strasbourg')>

3.4 Schritt 3. Mit den Primitiven von „ IBM Quantum “ ausführen

Im QAOA-Arbeitsablauf werden die optimalen QAOA-Parameter in einer iterativen Optimierungsschleife gefunden, die eine Reihe von Schaltkreisbewertungen durchführt und einen klassischen Optimierer verwendet, um die optimalen $\beta_k$ und $\gamma_k$ Parameter zu finden. Diese Ausführungsschleife wird in den folgenden Schritten ausgeführt:

1. Definieren Sie die Ausgangsparameter
2. Instanziieren Sie eine neue Session , die die Optimierungsschleife und das Primitiv enthält, das für die Abtastung der Schaltung verwendet wird
3. Sobald ein optimaler Parametersatz gefunden wurde, führen Sie den Schaltkreis ein letztes Mal aus, um eine endgültige Verteilung zu erhalten, die im Nachbearbeitungsschritt verwendet wird.

Schaltung mit Anfangsparametern definieren

Wir beginnen mit willkürlich gewählten Parametern.

initial_gamma = np.pi
initial_beta = np.pi / 2
init_params = [initial_gamma, initial_beta]

Backend und Ausführungsprimitive definieren

Verwenden Sie die „ IBM Quantum “-Primitive, um mit den „ IBM® “-Backends zu interagieren. Die beiden Primitive heißen „Sampler“ und „Estimator“, und die Wahl des Primitives hängt davon ab, welche Art von Messung Sie auf dem Quantencomputer durchführen möchten. Zur Minimierung von $H_C$ verwende den Schätzer, da die Messgröße der Kostenfunktion einfach der Erwartungswert von $\langle H_C \rangle$ ist.

Ausführen

Die Primitive bieten eine Vielzahl von Ausführungsmodi, um Arbeitslasten auf Quantengeräten zu planen, und ein QAOA-Workflow wird iterativ in einer Sitzung ausgeführt.

Sie können die samplerbasierte Kostenfunktion in die Minimierungsroutine SciPy einfügen, um die optimalen Parameter zu finden.

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

from qiskit_ibm_runtime import Session, EstimatorV2
from scipy.optimize import minimize

objective_func_vals = []  # Global variable
with Session(backend=backend) as session:
    # If using qiskit-ibm-runtime<0.24.0, change `mode=` to `session=`
    estimator = EstimatorV2(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)

Output:

 message: Optimization terminated successfully.
 success: True
  status: 1
     fun: -0.6557925874481715
       x: [ 2.873e+00  9.414e-01]
    nfev: 21
   maxcv: 0.0

Der Optimierer konnte die Kosten senken und bessere Parameter für die Schaltung finden.

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

Output:

Sobald Sie die optimalen Parameter für die Schaltung gefunden haben, können Sie diese Parameter zuweisen und die mit den optimierten Parametern erhaltene endgültige Verteilung stichprobenartig überprüfen. Hier sollte das Sampler-Primitiv verwendet werden, da es die Wahrscheinlichkeitsverteilung der Bitstring-Messungen ist, die dem optimalen Schnitt des Graphen entspricht.

Anmerkung: Dies bedeutet, dass ein Quantenzustand $\psi$ im Computer vorbereitet und dann gemessen wird. Durch eine Messung wird der Zustand in einen einzigen Berechnungsbasiszustand - zum Beispiel 010101110000... - umgewandelt, der einer Kandidatenlösung $x$ für unser anfängliches Optimierungsproblem entspricht ( $\max f(x)$ oder $\min f(x)$, je nach Aufgabe).

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

Output:

from qiskit_ibm_runtime import SamplerV2

# If using qiskit-ibm-runtime<0.24.0, change `mode=` to `backend=`
sampler = SamplerV2(mode=backend)

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

Output:

{12: 0.0652, 31: 0.0089, 4: 0.0085, 13: 0.0731, 26: 0.0256, 28: 0.0246, 17: 0.0405, 25: 0.0591, 20: 0.031, 15: 0.0221, 8: 0.017, 21: 0.0371, 14: 0.0461, 16: 0.0229, 19: 0.0723, 23: 0.0199, 22: 0.0478, 18: 0.0708, 24: 0.0165, 6: 0.0525, 7: 0.0155, 5: 0.0245, 3: 0.0231, 29: 0.0121, 30: 0.0062, 10: 0.0363, 1: 0.0097, 9: 0.042, 27: 0.0094, 11: 0.0349, 0: 0.0129, 2: 0.0119}

3.5 Schritt 4. Nachbearbeitung, Ergebnis im klassischen Format zurückgeben

Der Nachbearbeitungsschritt interpretiert die Stichprobenergebnisse, um eine Lösung für Ihr ursprüngliches Problem zu finden. In diesem Fall sind Sie an der Bitfolge mit der höchsten Wahrscheinlichkeit interessiert, da diese den optimalen Schnitt bestimmt. Die Symmetrien in dem Problem lassen vier mögliche Lösungen zu, und der Stichprobenprozess wird eine davon mit einer etwas höheren Wahrscheinlichkeit zurückgeben, aber Sie können in der unten aufgezeichneten Verteilung sehen, dass vier der Bitstrings deutlich wahrscheinlicher sind als der Rest.

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

Output:

Result bitstring: [1, 0, 1, 1, 0]

import matplotlib.pyplot as plt

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()[p[0].item()].set_color("tab:purple")
plt.show()

Output:

Den besten Schnitt visualisieren

Ausgehend von der optimalen Bitfolge können Sie dann diesen Schnitt auf dem ursprünglichen Graphen visualisieren.

colors = ["tab:grey" if i == 0 else "tab:purple" for i in most_likely_bitstring]
mpl_draw(graph, node_size=600, pos=pos, with_labels=True, labels=str, node_color=colors)

Output:

Und berechnen Sie den Wert des Schnitts. Die Lösung ist aufgrund von Rauschen nicht optimal (der Grenzwert der optimalen Lösung ist 5).

from typing import Sequence


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)

Output:

The value of the cut is: 5

Damit ist der kleine QAOA-Lehrgang abgeschlossen. In "Teil 2: Skalieren!" erfahren Sie, wie Sie QAOA an die Anforderungen von Versorgungsunternehmen anpassen können des Tutorials über den Algorithmus für die Quantenapproximation.

# Check Qiskit version
import qiskit

qiskit.__version__

Output:

'2.0.2'