Kostenfunktionen

In dieser Lektion lernst du, wie du eine Kostenfunktion auswertest:

• Zunächst lernst du die Qiskit Runtime Primitives kennen
• Du definierst eine Kostenfunktion $C(\vec\theta)$. Das ist eine problemspezifische Funktion, die das Ziel des Problems vorgibt, das der Optimierer minimieren (oder maximieren) soll
• Du legst eine Messstrategie mit den Qiskit Runtime Primitives fest, um Geschwindigkeit und Genauigkeit gegeneinander abzuwägen

 

Primitives

Alle physikalischen Systeme – ob klassisch oder quantenmechanisch – können in verschiedenen Zuständen vorliegen. Ein Auto auf einer Straße hat zum Beispiel eine bestimmte Masse, Position, Geschwindigkeit oder Beschleunigung, die seinen Zustand beschreiben. Quantensysteme können ebenfalls verschiedene Konfigurationen oder Zustände annehmen, unterscheiden sich aber von klassischen Systemen darin, wie Messungen und Zustandsentwicklung behandelt werden. Das führt zu einzigartigen Eigenschaften wie Superposition und Verschränkung, die ausschließlich in der Quantenmechanik auftreten. Genauso wie sich der Zustand eines Autos durch physikalische Größen wie Geschwindigkeit oder Beschleunigung beschreiben lässt, kann auch der Zustand eines Quantensystems durch Observablen beschrieben werden – mathematische Objekte, die messbare Größen repräsentieren.

In der Quantenmechanik werden Zustände durch normierte komplexe Spaltenvektoren, sogenannte Kets ($|\psi\rangle$), dargestellt, und Observablen sind hermitesche lineare Operatoren ($\hat{H}=\hat{H}^{\dagger}$), die auf die Kets wirken. Ein Eigenvektor ($|\lambda\rangle$) einer Observablen wird als Eigenzustand bezeichnet. Misst man eine Observable in einem ihrer Eigenzustände ($|\lambda\rangle$), erhält man den zugehörigen Eigenwert ($\lambda$) als Messergebnis.

Wenn du dich fragst, wie und was man an einem Quantensystem messen kann: Qiskit bietet zwei Primitives an:

• Sampler: Gegeben einen Quantenzustand $|\psi\rangle$, ermittelt dieses Primitive die Wahrscheinlichkeit jedes möglichen Rechenbasisszustands.
• Estimator: Gegeben eine quantenmechanische Observable $\hat{H}$ und einen Zustand $|\psi\rangle$, berechnet dieses Primitive den Erwartungswert von $\hat{H}$.

Das Sampler-Primitive

Das Sampler-Primitive berechnet die Wahrscheinlichkeit, jeden möglichen Zustand $|k\rangle$ aus der Rechenbasis zu erhalten, wenn ein Quantenkreis den Zustand $|\psi\rangle$ präpariert. Es berechnet

$$p_k = |\langle k | \psi \rangle|^2 \quad \forall k \in \mathbb{Z}_2^n \equiv \{0,1,\cdots,2^n-1\},$$

wobei $n$ die Anzahl der Qubits ist und $k$ die ganzzahlige Darstellung einer möglichen binären Ausgabezeichenkette $\{0,1\}^n$ (also ganze Zahlen zur Basis $2$).

Der Qiskit Runtime Sampler führt den Circuit mehrfach auf einem Quantengerät aus, nimmt bei jedem Durchlauf Messungen vor und rekonstruiert die Wahrscheinlichkeitsverteilung aus den erhaltenen Bitstrings. Je mehr Durchläufe (sogenannte Shots) durchgeführt werden, desto genauer sind die Ergebnisse – allerdings auf Kosten von mehr Zeit und Quantenressourcen.

Da die Anzahl möglicher Ausgaben jedoch exponentiell mit der Qubit-Anzahl $n$ wächst (nämlich $2^n$), muss auch die Anzahl der Shots exponentiell wachsen, um eine dichte Wahrscheinlichkeitsverteilung abzubilden. Deshalb ist Sampler nur für sparse Wahrscheinlichkeitsverteilungen effizient; der Zielzustand $|\psi\rangle$ muss als Linearkombination von Rechenbasisszuständen ausdrückbar sein, wobei die Anzahl der Terme höchstens polynomiell mit der Qubit-Anzahl wächst:

$$|\psi\rangle = \sum^{\text{Poly}(n)}_k w_k |k\rangle.$$

Der Sampler kann auch so konfiguriert werden, dass er Wahrscheinlichkeiten nur für einen Teilbereich des Circuits abruft, was einem Teilmenge aller möglichen Zustände entspricht.

Das Estimator-Primitive

Das Estimator-Primitive berechnet den Erwartungswert einer Observablen $\hat{H}$ für einen Quantenzustand $|\psi\rangle$; die Wahrscheinlichkeiten der Observablen können als $p_\lambda = |\langle\lambda|\psi\rangle|^2$ ausgedrückt werden, wobei $|\lambda\rangle$ die Eigenzustände der Observablen $\hat{H}$ sind. Der Erwartungswert ist dann als Durchschnitt aller möglichen Messergebnisse $\lambda$ (also der Eigenwerte der Observablen) des Zustands $|\psi\rangle$ definiert, gewichtet mit den entsprechenden Wahrscheinlichkeiten:

$$\langle\hat{H}\rangle_\psi := \sum_\lambda p_\lambda \lambda = \langle \psi | \hat{H} | \psi \rangle$$

Den Erwartungswert einer Observablen zu berechnen ist jedoch nicht immer möglich, da die Eigenbasis oft unbekannt ist. Der Qiskit Runtime Estimator verwendet einen komplexen algebraischen Prozess, um den Erwartungswert auf einem echten Quantengerät zu schätzen, indem er die Observable in eine Kombination anderer Observablen zerlegt, deren Eigenbasis bekannt ist.

Vereinfacht gesagt zerlegt Estimator jede Observable, die er nicht direkt messen kann, in einfachere, messbare Observablen – sogenannte Pauli-Operatoren. Jeder Operator lässt sich als Kombination von $4^n$ Pauli-Operatoren ausdrücken.

$$\hat{P}_k := \sigma_{k_{n-1}}\otimes \cdots \otimes \sigma_{k_0} \quad \forall k \in \mathbb{Z}_4^n \equiv \{0,1,\cdots,4^n-1\}, \\$$

sodass

$$\hat{H} = \sum^{4^n-1}_{k=0} w_k \hat{P}_k$$

wobei $n$ die Anzahl der Qubits ist, $k \equiv k_{n-1} \cdots k_0$ für $k_l \in \mathbb{Z}_4 \equiv \{0, 1, 2, 3\}$ (also ganze Zahlen zur Basis $4$) und $(\sigma_0, \sigma_1, \sigma_2, \sigma_3) := (I, X, Y, Z)$.

Nach dieser Zerlegung leitet Estimator für jede Observable $\hat{P}_k$ einen neuen Circuit $V_k|\psi\rangle$ aus dem ursprünglichen Circuit ab, um die Pauli-Observable in der Rechenbasis effektiv zu diagonalisieren und zu messen. Pauli-Observablen lassen sich leicht messen, weil $V_k$ im Voraus bekannt ist – was bei anderen Observablen im Allgemeinen nicht der Fall ist.

Für jedes $\hat{P}_{k}$ führt Estimator den entsprechenden Circuit mehrfach auf einem Quantengerät aus, misst den Ausgangszustand in der Rechenbasis und berechnet die Wahrscheinlichkeit $p_{kj}$, jede mögliche Ausgabe $j$ zu erhalten. Anschließend sucht es den Eigenwert $\lambda_{kj}$ von $P_k$ für jede Ausgabe $j$, multipliziert ihn mit $w_k$ und addiert alle Ergebnisse, um den Erwartungswert der Observablen $\hat{H}$ für den gegebenen Zustand $|\psi\rangle$ zu erhalten.

$$\langle\hat{H}\rangle_\psi = \sum_{k=0}^{4^n-1} w_k \sum_{j=0}^{2^n-1}p_{kj} \lambda_{kj},$$

Da die Berechnung des Erwartungswerts von $4^n$ Paulis unpraktisch wäre (da sie exponentiell wächst), kann Estimator nur dann effizient sein, wenn ein Großteil der $w_k$ gleich null ist (also eine sparse statt dichte Pauli-Zerlegung). Formal gesagt muss für eine effizient lösbare Berechnung die Anzahl der von null verschiedenen Terme höchstens polynomiell mit der Qubit-Anzahl $n$ wachsen: $\hat{H} = \sum^{\text{Poly}(n)}_k w_k \hat{P}_k.$

Der aufmerksame Leser bemerkt die implizite Annahme, dass auch das Wahrscheinlichkeits-Sampling effizient sein muss, wie es für Sampler erläutert wurde, was bedeutet:

$$\langle\hat{H}\rangle_\psi = \sum_{k}^{\text{Poly}(n)} w_k \sum_{j}^{\text{Poly}(n)}p_{kj} \lambda_{kj}.$$

Geführtes Beispiel zur Berechnung von Erwartungswerten

Nehmen wir den Ein-Qubit-Zustand $|+\rangle := H|0\rangle = \frac{1}{\sqrt{2}}(|0\rangle + |1\rangle)$ und die Observable

$$\begin{aligned} \hat{H} & = \begin{pmatrix} -1 & 2 \\ 2 & 1 \\ \end{pmatrix}\\[1mm] & = 2X - Z \end{aligned}$$

mit dem theoretischen Erwartungswert $\langle\hat{H}\rangle_+ = \langle+|\hat{H}|+\rangle = 2.$

Da wir nicht wissen, wie wir diese Observable direkt messen können, müssen wir sie umschreiben zu $\langle\hat{H}\rangle_+ = 2\langle X \rangle_+ - \langle Z \rangle_+$. Dass dies zum gleichen Ergebnis führt, ergibt sich daraus, dass $\langle+|X|+\rangle = 1$ und $\langle+|Z|+\rangle = 0$.

Schauen wir uns an, wie man $\langle X \rangle_+$ und $\langle Z \rangle_+$ direkt berechnet. Da $X$ und $Z$ nicht kommutieren (sie teilen also keine gemeinsame Eigenbasis), können sie nicht gleichzeitig gemessen werden. Deshalb brauchen wir Hilfsschaltkreise:

# Added by doQumentation — required packages for this notebook
!pip install -q matplotlib numpy qiskit qiskit-aer qiskit-ibm-runtime rustworkx

from qiskit import QuantumCircuit
from qiskit.quantum_info import SparsePauliOp

# The following code will work for any other initial single-qubit state and observable
original_circuit = QuantumCircuit(1)
original_circuit.h(0)

H = SparsePauliOp(["X", "Z"], [2, -1])

aux_circuits = []
for pauli in H.paulis:
    aux_circ = original_circuit.copy()
    aux_circ.barrier()
    if str(pauli) == "X":
        aux_circ.h(0)
    elif str(pauli) == "Y":
        aux_circ.sdg(0)
        aux_circ.h(0)
    else:
        aux_circ.id(0)
    aux_circ.measure_all()
    aux_circuits.append(aux_circ)

original_circuit.draw("mpl")

# Auxiliary circuit for X
aux_circuits[0].draw("mpl")

# Auxiliary circuit for Z
aux_circuits[1].draw("mpl")

Jetzt können wir die Berechnung manuell mit Sampler durchführen und die Ergebnisse mit Estimator überprüfen:

from qiskit.primitives import StatevectorSampler, StatevectorEstimator
from qiskit.result import QuasiDistribution
import numpy as np

## SAMPLER
shots = 10000
sampler = StatevectorSampler()
job = sampler.run(aux_circuits, shots=shots)

# Run the sampler job and step through results
expvals = []
for index, pauli in enumerate(H.paulis):
    data_pub = job.result()[index].data
    bitstrings = data_pub.meas.get_bitstrings()
    counts = data_pub.meas.get_counts()
    quasi_dist = QuasiDistribution(
        {outcome: freq / shots for outcome, freq in counts.items()}
    )

    # Use the probabilities and known eigenvalues of Pauli operators to estimate the expectation value.
    val = 0

    if str(pauli) == "X":
        val += -1 * quasi_dist.get(1, 0)
        val += 1 * quasi_dist.get(0, 0)

    if str(pauli) == "Y":
        val += -1 * quasi_dist.get(1, 0)
        val += 1 * quasi_dist.get(0, 0)

    if str(pauli) == "Z":
        val += 1 * quasi_dist.get(0, 0)
        val += -1 * quasi_dist.get(1, 0)

    expvals.append(val)

# Print expectation values

print("Sampler results:")
for pauli, expval in zip(H.paulis, expvals):
    print(f"  >> Expected value of {str(pauli)}: {expval:.5f}")

total_expval = np.sum(H.coeffs * expvals).real
print(f"  >> Total expected value: {total_expval:.5f}")

# Use estimator for comparison
observables = [
    *H.paulis,
    H,
]  # Note: run for individual Paulis as well as full observable H

estimator = StatevectorEstimator()
job = estimator.run([(original_circuit, observables)])
estimator_expvals = job.result()[0].data.evs

# Print results
print("Estimator results:")
for obs, expval in zip(observables, estimator_expvals):
    if obs is not H:
        print(f"  >> Expected value of {str(obs)}: {expval:.5f}")
    else:
        print(f"  >> Total expected value: {expval:.5f}")

Sampler results:
  >> Expected value of X: 1.00000
  >> Expected value of Z: 0.00420
  >> Total expected value: 1.99580
Estimator results:
  >> Expected value of X: 1.00000
  >> Expected value of Z: 0.00000
  >> Total expected value: 2.00000

Mathematische Grundlagen (optional)

Wenn man $|\psi\rangle$ bezüglich der Eigenzustandsbasis von $\hat{H}$ ausdrückt, $|\psi\rangle = \sum_\lambda a_\lambda |\lambda\rangle$, ergibt sich:

$$\begin{aligned} \langle \psi | \hat{H} | \psi \rangle & = \bigg(\sum_{\lambda'}a^*_{\lambda'} \langle \lambda'|\bigg) \hat{H} \bigg(\sum_{\lambda} a_\lambda | \lambda\rangle\bigg)\\[1mm] & = \sum_{\lambda}\sum_{\lambda'} a^*_{\lambda'}a_{\lambda} \langle \lambda'|\hat{H}| \lambda\rangle\\[1mm] & = \sum_{\lambda}\sum_{\lambda'} a^*_{\lambda'}a_{\lambda} \lambda \langle \lambda'| \lambda\rangle\\[1mm] & = \sum_{\lambda}\sum_{\lambda'} a^*_{\lambda'}a_{\lambda} \lambda \cdot \delta_{\lambda, \lambda'}\\[1mm] & = \sum_\lambda |a_\lambda|^2 \lambda\\[1mm] & = \sum_\lambda p_\lambda \lambda\\[1mm] \end{aligned}$$

Da wir die Eigenwerte oder Eigenzustände der Ziel-Observablen $\hat{H}$ nicht kennen, müssen wir zunächst ihre Diagonalisierung betrachten. Da $\hat{H}$ hermitesch ist, gibt es eine unitäre Transformation $V$ mit $\hat{H}=V^\dagger \Lambda V$, wobei $\Lambda$ die diagonale Eigenwertmatrix ist, sodass $\langle j | \Lambda | k \rangle = 0$ für $j\neq k$ und $\langle j | \Lambda | j \rangle = \lambda_j$.

Daraus folgt, dass sich der Erwartungswert umschreiben lässt zu:

$$\begin{aligned} \langle\psi|\hat{H}|\psi\rangle & = \langle\psi|V^\dagger \Lambda V|\psi\rangle\\[1mm] & = \langle\psi|V^\dagger \bigg(\sum_{j=0}^{2^n-1} |j\rangle \langle j|\bigg) \Lambda \bigg(\sum_{k=0}^{2^n-1} |k\rangle \langle k|\bigg) V|\psi\rangle\\[1mm] & = \sum_{j=0}^{2^n-1} \sum_{k=0}^{2^n-1}\langle\psi|V^\dagger |j\rangle \langle j| \Lambda |k\rangle \langle k| V|\psi\rangle\\[1mm] & = \sum_{j=0}^{2^n-1}\langle\psi|V^\dagger |j\rangle \langle j| \Lambda |j\rangle \langle j| V|\psi\rangle\\[1mm] & = \sum_{j=0}^{2^n-1}|\langle j| V|\psi\rangle|^2 \lambda_j\\[1mm] \end{aligned}$$

Da die Wahrscheinlichkeit, $| j\rangle$ zu messen, im Zustand $|\phi\rangle = V |\psi\rangle$ gleich $p_j = |\langle j|\phi \rangle|^2$ ist, lässt sich der obige Erwartungswert ausdrücken als:

$$\langle\psi|\hat{H}|\psi\rangle = \sum_{j=0}^{2^n-1} p_j \lambda_j.$$

Es ist sehr wichtig zu beachten, dass die Wahrscheinlichkeiten aus dem Zustand $V |\psi\rangle$ und nicht aus $|\psi\rangle$ stammen. Deshalb ist die Matrix $V$ unbedingt erforderlich. Vielleicht fragst du dich, wie man $V$ und die Eigenwerte $\Lambda$ erhält. Wenn du die Eigenwerte bereits kennen würdest, bräuchtest du keinen Quantencomputer mehr – denn das Ziel variationeller Algorithmen ist es ja gerade, diese Eigenwerte von $\hat{H}$ zu finden.

Zum Glück gibt es einen Ausweg: Jede $2^n \times 2^n$-Matrix lässt sich als Linearkombination von $4^n$ Tensorprodukten von $n$ Pauli-Matrizen und Identitäten schreiben, die alle sowohl hermitesch als auch unitär sind und für die $V$ und $\Lambda$ bekannt sind. Genau das macht der Runtime Estimator intern, indem er jedes Operator-Objekt in ein SparsePauliOp zerlegt.

Die verfügbaren Operatoren sind:

$$\begin{array}{c|c|c|c} \text{Operator} & \sigma & V & \Lambda \\[1mm] \hline I & \sigma_0 = \begin{pmatrix} 1 & 0 \\ 0 & 1 \end{pmatrix} & V_0 = I & \Lambda_0 = I = \begin{pmatrix} 1 & 0 \\ 0 & 1 \end{pmatrix} \\[4mm] X & \sigma_1 = \begin{pmatrix} 0 & 1 \\ 1 & 0 \end{pmatrix} & V_1 = H =\frac{1}{\sqrt{2}} \begin{pmatrix} 1 & 1 \\ 1 & -1 \end{pmatrix} & \Lambda_1 = \sigma_3 = \begin{pmatrix} 1 & 0 \\ 0 & -1 \end{pmatrix} \\[4mm] Y & \sigma_2 = \begin{pmatrix} 0 & -i \\ i & 0 \end{pmatrix} & V_2 = HS^\dagger =\frac{1}{\sqrt{2}} \begin{pmatrix} 1 & 1 \\ 1 & -1 \end{pmatrix}\cdot \begin{pmatrix} 1 & 0 \\ 0 & -i \end{pmatrix} = \frac{1}{\sqrt{2}} \begin{pmatrix} 1 & -i \\ 1 & i \end{pmatrix}\quad & \Lambda_2 = \sigma_3 = \begin{pmatrix} 1 & 0 \\ 0 & -1 \end{pmatrix} \\[4mm] Z & \sigma_3 = \begin{pmatrix} 1 & 0 \\ 0 & -1 \end{pmatrix} & V_3 = I & \Lambda_3 = \sigma_3 = \begin{pmatrix} 1 & 0 \\ 0 & -1 \end{pmatrix} \end{array}$$

Schreiben wir also $\hat{H}$ bezüglich der Paulis und Identitäten um:

$$\hat{H} = \sum_{k_{n-1}=0}^3... \sum_{k_0=0}^3 w_{k_{n-1}...k_0} \sigma_{k_{n-1}}\otimes ... \otimes \sigma_{k_0} = \sum_{k=0}^{4^n-1} w_k \hat{P}_k,$$

wobei $k = \sum_{l=0}^{n-1} 4^l k_l \equiv k_{n-1}...k_0$ für $k_{n-1},...,k_0\in \{0,1,2,3\}$ (also Basis $4$) und $\hat{P}_{k} := \sigma_{k_{n-1}}\otimes ... \otimes \sigma_{k_0}$:

$$\begin{aligned} \langle\psi|\hat{H}|\psi\rangle & = \sum_{k=0}^{4^n-1} w_k \sum_{j=0}^{2^n-1}|\langle j| V_k|\psi\rangle|^2 \langle j| \Lambda_k |j\rangle \\[1mm] & = \sum_{k=0}^{4^n-1} w_k \sum_{j=0}^{2^n-1}p_{kj} \lambda_{kj}, \\[1mm] \end{aligned}$$

wobei $V_k := V_{k_{n-1}}\otimes ... \otimes V_{k_0}$ und $\Lambda_k := \Lambda_{k_{n-1}}\otimes ... \otimes \Lambda_{k_0}$, sodass $\hat{P_k}=V_k^\dagger \Lambda_k V_k$ gilt.

Kostenfunktionen

Kostenfunktionen beschreiben allgemein das Ziel eines Problems und geben an, wie gut ein Testzustand dieses Ziel erfüllt. Diese Definition lässt sich auf verschiedene Anwendungsgebiete übertragen: Chemie, maschinelles Lernen, Finanzen, Optimierung und viele weitere.

Betrachten wir ein einfaches Beispiel: die Suche nach dem Grundzustand eines Systems. Unser Ziel ist es, den Erwartungswert der Observablen, die die Energie darstellt (Hamiltonian $\hat{\mathcal{H}}$), zu minimieren:

$$\min_{\vec\theta} \langle\psi(\vec\theta)|\hat{\mathcal{H}}|\psi(\vec\theta)\rangle$$

Wir können den Estimator verwenden, um den Erwartungswert auszuwerten, und diesen Wert an einen Optimierer übergeben. Wenn die Optimierung erfolgreich ist, liefert sie einen optimalen Parametersatz $\vec\theta^*$, aus dem sich der vorgeschlagene Lösungszustand $|\psi(\vec\theta^*)\rangle$ konstruieren und der beobachtete Erwartungswert als $C(\vec\theta^*)$ berechnen lässt.

Beachte, dass wir die Kostenfunktion nur innerhalb der eingeschränkten Menge von Zuständen minimieren können, die wir betrachten. Das führt zu zwei möglichen Szenarien:

• Unser Ansatz deckt den Lösungszustand nicht im Suchraum ab: In diesem Fall wird der Optimierer die Lösung nie finden, und wir müssen andere Ansätze ausprobieren, die den Suchraum besser abbilden können.
• Der Optimierer findet die gültige Lösung nicht: Optimierung kann global oder lokal definiert sein. Was das bedeutet, werden wir im späteren Abschnitt näher beleuchten.

Insgesamt führen wir eine klassische Optimierungsschleife durch, stützen uns aber bei der Auswertung der Kostenfunktion auf einen Quantencomputer. Aus dieser Perspektive könnte man die Optimierung als rein klassisches Unterfangen betrachten, bei dem wir jedes Mal, wenn der Optimierer die Kostenfunktion auswerten muss, ein quantenmechanisches Black-Box-Orakel aufrufen.

def cost_func_vqe(params, circuit, hamiltonian, estimator):
    """Return estimate of energy from estimator

    Parameters:
        params (ndarray): Array of ansatz parameters
        ansatz (QuantumCircuit): Parameterized ansatz circuit
        hamiltonian (SparsePauliOp): Operator representation of Hamiltonian
        estimator (Estimator): Estimator primitive instance

    Returns:
        float: Energy estimate
    """
    pub = (circuit, hamiltonian, params)
    cost = estimator.run([pub]).result()[0].data.evs
    return cost

from qiskit.circuit.library import TwoLocal

observable = SparsePauliOp.from_list([("XX", 1), ("YY", -3)])

reference_circuit = QuantumCircuit(2)
reference_circuit.x(0)

variational_form = TwoLocal(
    2,
    rotation_blocks=["rz", "ry"],
    entanglement_blocks="cx",
    entanglement="linear",
    reps=1,
)
ansatz = reference_circuit.compose(variational_form)

theta_list = (2 * np.pi * np.random.rand(1, 8)).tolist()
ansatz.decompose().draw("mpl")

Zunächst führen wir dies mit einem Simulator durch: dem StatevectorEstimator. Das empfiehlt sich in der Regel zum Debuggen, aber direkt im Anschluss werden wir die Berechnung auf echter Quantenhardware wiederholen. Immer häufiger sind interessante Probleme klassisch nicht mehr simulierbar – zumindest nicht ohne modernste Supercomputer-Infrastruktur.

estimator = StatevectorEstimator()
cost = cost_func_vqe(theta_list, ansatz, observable, estimator)
print(cost)

[-0.58744589]

Jetzt führen wir die Berechnung auf einem echten Quantencomputer durch. Beachte die geänderte Syntax. Die Schritte mit dem pass_manager werden im nächsten Beispiel näher erläutert. Ein besonders wichtiger Schritt bei variationellen Algorithmen ist die Verwendung einer Qiskit Runtime Session. Wenn du eine Session startest, kannst du mehrere Iterationen eines variationellen Algorithmus ausführen, ohne jedes Mal neu in die Warteschlange eingereiht zu werden, wenn Parameter aktualisiert werden. Das ist wichtig, wenn die Wartezeiten lang sind und/oder viele Iterationen benötigt werden. Runtime Sessions können nur von Partnern im IBM Quantum® Network genutzt werden. Falls du keinen Zugang zu Sessions hast, kannst du die Anzahl der Iterationen reduzieren und die zuletzt verwendeten Parameter für spätere Durchläufe speichern. Bei zu vielen Iterationen oder zu langen Wartezeiten kann Fehlercode 1217 auftreten, der auf zu große Verzögerungen zwischen Job-Einreichungen hinweist.

# Estimated usage: < 1 min. Benchmarked at 7 seconds on an Eagle processor
# Load necessary packages:

from qiskit_ibm_runtime import (
    QiskitRuntimeService,
    Session,
    EstimatorOptions,
    EstimatorV2 as Estimator,
)
from qiskit.transpiler.preset_passmanagers import generate_preset_pass_manager

# Select the least busy backend:

service = QiskitRuntimeService()
backend = service.least_busy(
    operational=True, min_num_qubits=ansatz.num_qubits, simulator=False
)
# Or get a specific backend:
# backend = service.backend("ibm_brisbane")

# Use a pass manager to transpile the circuit and observable for the specific backend being used:

pm = generate_preset_pass_manager(backend=backend, optimization_level=1)
isa_ansatz = pm.run(ansatz)
isa_observable = observable.apply_layout(layout=isa_ansatz.layout)

# Set estimator options
estimator_options = EstimatorOptions(resilience_level=1, default_shots=10_000)

# Open a Runtime session:

with Session(backend=backend) as session:
    estimator = Estimator(mode=session, options=estimator_options)
    cost = cost_func_vqe(theta_list, isa_ansatz, isa_observable, estimator)

session.close()
print(cost)

Die Werte aus den beiden obigen Berechnungen sind sehr ähnlich. Techniken zur Verbesserung der Ergebnisse werden weiter unten besprochen.

Beispiel: Abbildung auf nicht-physikalische Systeme

Das Maximum-Cut-Problem (Max-Cut) ist ein kombinatorisches Optimierungsproblem, bei dem die Knoten eines Graphen in zwei disjunkte Mengen aufgeteilt werden sollen, sodass die Anzahl der Kanten zwischen den beiden Mengen maximiert wird. Formal gesagt: Gegeben ein ungerichteter Graph $G=(V,E)$, wobei $V$ die Menge der Knoten und $E$ die Menge der Kanten ist, fragt das Max-Cut-Problem danach, die Knoten in zwei disjunkte Teilmengen $S$ und $T$ zu partitionieren, sodass die Anzahl der Kanten mit einem Endpunkt in $S$ und dem anderen in $T$ maximiert wird.

Max-Cut lässt sich auf verschiedene Probleme anwenden, darunter Clustering, Netzwerkdesign, Phasenübergänge und mehr. Beginnen wir mit der Erstellung eines Problemgraphen:

import rustworkx as rx
from rustworkx.visualization import mpl_draw

n = 4
G = rx.PyGraph()
G.add_nodes_from(range(n))
# The edge syntax is (start, end, weight)
edges = [(0, 1, 1.0), (0, 2, 1.0), (0, 3, 1.0), (1, 2, 1.0), (2, 3, 1.0)]
G.add_edges_from(edges)

mpl_draw(
    G, pos=rx.shell_layout(G), with_labels=True, edge_labels=str, node_color="#1192E8"
)

Dieses Problem lässt sich als binäres Optimierungsproblem formulieren. Für jeden Knoten $0 \leq i < n$, wobei $n$ die Anzahl der Knoten des Graphen ist (hier $n=4$), betrachten wir die binäre Variable $x_i$. Diese Variable hat den Wert $1$, wenn Knoten $i$ zu der Gruppe gehört, die wir mit $1$ bezeichnen, und $0$, wenn er zur anderen Gruppe mit der Bezeichnung $0$ gehört. Außerdem bezeichnen wir mit $w_{ij}$ (Element $(i,j)$ der Adjazenzmatrix $w$) das Gewicht der Kante von Knoten $i$ zu Knoten $j$. Da der Graph ungerichtet ist, gilt $w_{ij}=w_{ji}$. Wir können unser Problem nun als Maximierung der folgenden Kostenfunktion formulieren:

$$\begin{aligned} C(\vec{x}) & =\sum_{i,j=0}^n w_{ij} x_i(1-x_j)\\[1mm] & = \sum_{i,j=0}^n w_{ij} x_i - \sum_{i,j=0}^n w_{ij} x_ix_j\\[1mm] & = \sum_{i,j=0}^n w_{ij} x_i - \sum_{i=0}^n \sum_{j=0}^i 2w_{ij} x_ix_j \end{aligned}$$

Um dieses Problem mit einem Quantencomputer zu lösen, drücken wir die Kostenfunktion als Erwartungswert einer Observablen aus. Die von Qiskit nativ unterstützten Observablen bestehen jedoch aus Pauli-Operatoren mit den Eigenwerten $1$ und $-1$ statt $0$ und $1$. Deshalb nehmen wir folgende Variablensubstitution vor:

Dabei ist $\vec{x}=(x_0,x_1,\cdots ,x_{n-1})$. Wir können die Adjazenzmatrix $w$ verwenden, um komfortabel auf die Gewichte aller Kanten zuzugreifen. Damit ergibt sich unsere Kostenfunktion:

$$z_i = 1-2x_i \rightarrow x_i = \frac{1-z_i}{2}$$

Daraus folgt:

$$\begin{array}{lcl} x_i=0 & \rightarrow & z_i=1 \\ x_i=1 & \rightarrow & z_i=-1.\end{array}$$

Die neue Kostenfunktion, die wir maximieren wollen, lautet also:

$$\begin{aligned} C(\vec{z}) & = \sum_{i,j=0}^n w_{ij} \bigg(\frac{1-z_i}{2}\bigg)\bigg(1-\frac{1-z_j}{2}\bigg)\\[1mm] & = \sum_{i,j=0}^n \frac{w_{ij}}{4} - \sum_{i,j=0}^n \frac{w_{ij}}{4} z_iz_j\\[1mm] & = \sum_{i=0}^n \sum_{j=0}^i \frac{w_{ij}}{2} - \sum_{i=0}^n \sum_{j=0}^i \frac{w_{ij}}{2} z_iz_j \end{aligned}$$

Da ein Quantencomputer von Natur aus Minima sucht (in der Regel den niedrigsten Energiezustand) statt Maxima, werden wir nicht $C(\vec{z})$ maximieren, sondern stattdessen minimieren:

$$-C(\vec{z}) = \sum_{i=0}^n \sum_{j=0}^i \frac{w_{ij}}{2} z_iz_j - \sum_{i=0}^n \sum_{j=0}^i \frac{w_{ij}}{2}$$

Da wir nun eine Kostenfunktion haben, deren Variablen die Werte $-1$ und $1$ annehmen können, können wir folgende Analogie zum Pauli-$Z$ herstellen:

$$z_i \equiv Z_i = \overbrace{I}^{n-1}\otimes ... \otimes \overbrace{Z}^{i} \otimes ... \otimes \overbrace{I}^{0}$$

Mit anderen Worten: Die Variable $z_i$ entspricht einem $Z$-Gate, das auf Qubit $i$ wirkt. Außerdem gilt:

$$Z_i|x_{n-1}\cdots x_0\rangle = z_i|x_{n-1}\cdots x_0\rangle \rightarrow \langle x_{n-1}\cdots x_0 |Z_i|x_{n-1}\cdots x_0\rangle = z_i$$

Die Observable, die wir betrachten, ist daher:

$$\hat{H} = \sum_{i=0}^n \sum_{j=0}^i \frac{w_{ij}}{2} Z_iZ_j$$

zu der wir anschließend den konstanten Term addieren müssen:

$$\texttt{offset} = - \sum_{i=0}^n \sum_{j=0}^i \frac{w_{ij}}{2}$$

Der Operator ist eine Linearkombination von Termen mit Z-Operatoren auf durch eine Kante verbundenen Knoten (beachte, dass das 0. Qubit ganz rechts steht): $IIZZ + IZIZ + IZZI + ZIIZ + ZZII$. Sobald der Operator konstruiert ist, lässt sich der Ansatz für den QAOA-Algorithmus einfach mit dem QAOAAnsatz-Circuit aus der Qiskit Circuit Library erstellen.

from qiskit.circuit.library import QAOAAnsatz
from qiskit.quantum_info import SparsePauliOp

hamiltonian = SparsePauliOp.from_list(
    [("IIZZ", 1), ("IZIZ", 1), ("IZZI", 1), ("ZIIZ", 1), ("ZZII", 1)]
)

ansatz = QAOAAnsatz(hamiltonian, reps=2)
# Draw
ansatz.decompose(reps=3).draw("mpl")

# Sum the weights, and divide by 2

offset = -sum(edge[2] for edge in edges) / 2
print(f"""Offset: {offset}""")

Offset: -2.5

Da der Runtime Estimator direkt einen Hamiltonian und einen parametrisierten Ansatz entgegennimmt und die benötigte Energie zurückgibt, ist die Kostenfunktion für eine QAOA-Instanz recht einfach:

def cost_func(params, ansatz, hamiltonian, estimator):
    """Return estimate of energy from estimator

    Parameters:
        params (ndarray): Array of ansatz parameters
        ansatz (QuantumCircuit): Parameterized ansatz circuit
        hamiltonian (SparsePauliOp): Operator representation of Hamiltonian
        estimator (Estimator): Estimator primitive instance

    Returns:
        float: Energy estimate
    """
    pub = (ansatz, hamiltonian, params)
    cost = estimator.run([pub]).result()[0].data.evs
    #    cost = estimator.run(ansatz, hamiltonian, parameter_values=params).result().values[0]
    return cost

import numpy as np

x0 = 2 * np.pi * np.random.rand(ansatz.num_parameters)

estimator = StatevectorEstimator()
cost = cost_func_vqe(x0, ansatz, hamiltonian, estimator)
print(cost)

1.473098768180865

# Estimated usage: < 1 min, benchmarked at 6 seconds on ibm_osaka, 5-23-24
# Load some necessary packages:

from qiskit_ibm_runtime import QiskitRuntimeService
from qiskit_ibm_runtime import Session, EstimatorV2 as Estimator

# Select the least busy backend:

backend = service.least_busy(
    operational=True, min_num_qubits=ansatz.num_qubits, simulator=False
)

# Or get a specific backend:
# backend = service.backend("ibm_brisbane")

# Use a pass manager to transpile the circuit and observable for the specific backend being used:

pm = generate_preset_pass_manager(backend=backend, optimization_level=1)
isa_ansatz = pm.run(ansatz)
isa_hamiltonian = hamiltonian.apply_layout(layout=isa_ansatz.layout)

# Set estimator options
estimator_options = EstimatorOptions(resilience_level=1, default_shots=10_000)

# Open a Runtime session:

with Session(backend=backend) as session:
    estimator = Estimator(mode=session, options=estimator_options)
    cost = cost_func_vqe(x0, isa_ansatz, isa_hamiltonian, estimator)

# Close session after done
session.close()
print(cost)

1.1120776913677988

Wir werden auf dieses Beispiel in den Anwendungen zurückkommen, um zu zeigen, wie ein Optimierer genutzt werden kann, um den Suchraum zu durchlaufen. Im Allgemeinen umfasst das:

• Einsatz eines Optimierers zur Bestimmung optimaler Parameter
• Einbinden der optimalen Parameter in den Ansatz zur Berechnung der Eigenwerte
• Übersetzen der Eigenwerte in unsere Problemdefinition

Messstrategie: Geschwindigkeit versus Genauigkeit

Wie bereits erwähnt, verwenden wir einen verrauschten Quantencomputer als Black-Box-Orakel, wobei Rauschen die abgerufenen Werte nicht-deterministisch machen kann. Das führt zu zufälligen Schwankungen, die die Konvergenz bestimmter Optimierer auf eine vorgeschlagene Lösung erschweren oder sogar vollständig verhindern können. Das ist ein allgemeines Problem, das wir angehen müssen, während wir schrittweise den Quantennutzen erkunden und auf einen Quantenvorteil hinarbeiten:

Wir können die Fehlerunterdrückungs- und Fehlerminderungsoptionen der Qiskit Runtime Primitives nutzen, um Rauschen zu begegnen und den Nutzen heutiger Quantencomputer zu maximieren.

Fehlerunterdrückung

Fehlerunterdrückung bezeichnet Techniken zur Optimierung und Transformation eines Circuits während der Kompilierung, um Fehler zu minimieren. Das ist eine grundlegende Fehlerbehandlungstechnik, die in der Regel zu einem gewissen klassischen Vorverarbeitungs-Overhead führt. Dieser Overhead umfasst das Transpilieren von Circuits für die Ausführung auf Quantenhardware durch:

• Ausdrücken des Circuits mit den nativen Gates des Quantensystems
• Abbilden der virtuellen Qubits auf physikalische Qubits
• Hinzufügen von SWAPs gemäß den Konnektivitätsanforderungen
• Optimieren von 1Q- und 2Q-Gates
• Hinzufügen von Dynamical Decoupling für inaktive Qubits, um Dekohärenzeffekte zu verhindern.

Primitives ermöglichen den Einsatz von Fehlerunterdrückungstechniken durch Setzen der Option optimization_level und die Auswahl erweiterter Transpilierungsoptionen. In einem späteren Kurs werden wir verschiedene Circuit-Konstruktionsmethoden zur Ergebnisverbesserung kennenlernen, aber für die meisten Fälle empfehlen wir optimization_level=3.

Den Wert einer höheren Optimierungsstufe im Transpilierungsprozess veranschaulichen wir anhand eines Beispielcircuits mit einem einfachen idealen Verhalten.

from qiskit.circuit import Parameter, QuantumCircuit
from qiskit.quantum_info import SparsePauliOp

theta = Parameter("theta")

qc = QuantumCircuit(2)
qc.x(1)
qc.h(0)
qc.cp(theta, 0, 1)
qc.h(0)
observables = SparsePauliOp.from_list([("ZZ", 1)])

qc.draw("mpl")

Der obige Circuit kann sinusförmige Erwartungswerte der gegebenen Observablen liefern, wenn wir Phasen über ein geeignetes Intervall wie $[0,2\pi]$ einfügen.

## Setup phases
import numpy as np

phases = np.linspace(0, 2 * np.pi, 50)

# phases need to be expressed as a list of lists in order to work
individual_phases = [[phase] for phase in phases]

Wir nutzen einen Simulator, um den Nutzen einer optimierten Transpilierung zu zeigen. Im Anschluss kehren wir zur echten Hardware zurück, um den Nutzen der Fehlerminderung zu demonstrieren. Wir verwenden QiskitRuntimeService, um ein echtes Backend zu erhalten (hier ibm_brisbane), und nutzen AerSimulator, um dieses Backend einschließlich seines Rauschverhaltens zu simulieren.

from qiskit_ibm_runtime import QiskitRuntimeService
from qiskit_aer import AerSimulator

# get a real backend from the runtime service
service = QiskitRuntimeService()
backend = service.backend("ibm_brisbane")

# generate a simulator that mimics the real quantum system with the latest calibration results
backend_sim = AerSimulator.from_backend(backend)

Jetzt können wir einen Pass Manager verwenden, um den Circuit in die „Instruction Set Architecture" (ISA) des Backends zu transpilieren. Das ist eine neue Anforderung in Qiskit Runtime: Alle an ein Backend übermittelten Circuits müssen den Einschränkungen des Backend-Targets entsprechen – sie müssen also in der ISA des Geräts geschrieben sein, also der Menge der Instruktionen, die das Gerät verstehen und ausführen kann. Diese Zieleinschränkungen werden durch Faktoren wie die nativen Basistransformationen des Geräts, seine Qubit-Konnektivität und – wenn relevant – seine Puls- und andere Instruktions-Timing-Spezifikationen definiert.

Beachte, dass wir dies im vorliegenden Fall zweimal tun: einmal mit optimization_level = 0 und einmal mit dem Wert 3. Jedes Mal verwenden wir das Estimator-Primitive, um die Erwartungswerte der Observablen bei verschiedenen Phasenwerten zu schätzen.

# Import estimator and specify that we are using the simulated backend:

from qiskit_ibm_runtime import EstimatorV2 as Estimator

estimator = Estimator(mode=backend_sim)

circuit = qc

# Use a pass manager to transpile the circuit and observable for the backend being simulated.
# Start with no optimization:

from qiskit.transpiler.preset_passmanagers import generate_preset_pass_manager

pm = generate_preset_pass_manager(backend=backend_sim, optimization_level=0)
isa_circuit = pm.run(circuit)
isa_observables = observables.apply_layout(layout=isa_circuit.layout)

noisy_exp_values = []
pub = (isa_circuit, isa_observables, [individual_phases])
cost = estimator.run([pub]).result()[0].data.evs
noisy_exp_values = cost[0]

# Repeat above steps, but now with optimization = 3:

exp_values_with_opt_es = []
pm = generate_preset_pass_manager(backend=backend_sim, optimization_level=3)
isa_circuit = pm.run(circuit)
isa_observables = observables.apply_layout(layout=isa_circuit.layout)

pub = (isa_circuit, isa_observables, [individual_phases])
cost = estimator.run([pub]).result()[0].data.evs
exp_values_with_opt_es = cost[0]

Schließlich können wir die Ergebnisse darstellen. Wir sehen, dass die Genauigkeit der Berechnung schon ohne Optimierung recht gut war, aber durch Erhöhen der Optimierungsstufe auf 3 deutlich verbessert wurde. Beachte, dass bei tieferen, komplexeren Circuits der Unterschied zwischen Optimierungsstufe 0 und 3 deutlich ausgeprägter sein dürfte. Hier handelt es sich um einen sehr einfachen Circuit als Spielzeugmodell.

import matplotlib.pyplot as plt

plt.plot(phases, noisy_exp_values, "o", label="opt=0")
plt.plot(phases, exp_values_with_opt_es, "o", label="opt=3")
plt.plot(phases, 2 * np.sin(phases / 2) ** 2 - 1, label="ideal")
plt.ylabel("Expectation")
plt.legend()
plt.show()

Fehlerminderung

Fehlerminderung bezeichnet Techniken, die es ermöglichen, Circuit-Fehler zu reduzieren, indem das Geräterauschen zum Ausführungszeitpunkt modelliert wird. In der Regel entsteht dabei ein Quantum-Vorverarbeitungs-Overhead für das Modelltraining sowie ein klassischer Nachverarbeitungs-Overhead, um Fehler in den Rohergebnissen mithilfe des erstellten Modells zu korrigieren.

Die Option resilience_level der Qiskit Runtime Primitives gibt an, wie viel Widerstandsfähigkeit gegenüber Fehlern aufgebaut werden soll. Höhere Stufen erzeugen genauere Ergebnisse, erfordern aber aufgrund des Quantum-Sampling-Overheads mehr Verarbeitungszeit. Mit den Resilience-Stufen lässt sich der Kompromiss zwischen Kosten und Genauigkeit beim Einsatz von Fehlerminderung steuern.

Bei jeder Fehlerminderungstechnik erwarten wir, dass der Bias in unseren Ergebnissen im Vergleich zur ungeminderten Version reduziert wird. In manchen Fällen kann der Bias sogar vollständig verschwinden. Das hat jedoch seinen Preis: Wenn wir den Bias in unseren geschätzten Größen reduzieren, nimmt die statistische Variabilität zu (also die Varianz). Dem können wir begegnen, indem wir die Anzahl der Shots pro Circuit im Sampling-Prozess erhöhen. Da das über den Aufwand zur Bias-Reduktion hinausgeht, ist es standardmäßig nicht aktiviert. Dieses Verhalten lässt sich einfach durch Anpassen der Anzahl der Shots pro Circuit in options.executions.shots aktivieren, wie im folgenden Beispiel gezeigt.

In diesem Kurs betrachten wir diese Fehlerminderungsmodelle auf einem hohen Abstraktionsniveau, um die Fehlerminderung zu veranschaulichen, die Qiskit Runtime Primitives durchführen können, ohne alle Implementierungsdetails vorauszusetzen.

Twirled Readout Error Extinction (T-REx)

Twirled Readout Error Extinction (T-REx) verwendet eine Technik namens Pauli Twirling, um das während des Quantenmessvorgangs eingeführte Rauschen zu reduzieren. Diese Technik setzt keine spezifische Rauschform voraus, was sie sehr allgemein und effektiv macht.

Gesamtablauf:

1. Daten für den Nullzustand mit zufälligen Bit-Flips erfassen (Pauli X vor der Messung)
2. Daten für den gewünschten (verrauschten) Zustand mit zufälligen Bit-Flips erfassen (Pauli X vor der Messung)
3. Die spezielle Funktion für jeden Datensatz berechnen und dividieren.

 

Das lässt sich mit options.resilience_level = 1 einstellen, wie im folgenden Beispiel gezeigt.

Zero Noise Extrapolation

Zero Noise Extrapolation (ZNE) funktioniert, indem zunächst das Rauschen im Circuit, der den gewünschten Quantenzustand präpariert, verstärkt wird. Dann werden Messungen für verschiedene Rauschstärken durchgeführt, und aus diesen Messungen wird das rauschfreie Ergebnis extrapoliert.

Gesamtablauf:

1. Circuit-Rauschen für verschiedene Rauschfaktoren verstärken
2. Jeden rauschverstärkten Circuit ausführen
3. Zum Null-Rausch-Grenzwert zurückextrapolieren

 

Das lässt sich mit options.resilience_level = 2 einstellen. Eine weitere Optimierung ist durch das Erkunden verschiedener noise_factors, noise_amplifiers und extrapolators möglich, aber das würde den Rahmen dieses Kurses sprengen. Wir empfehlen dir, mit diesen Optionen wie hier beschrieben zu experimentieren.

Jede Methode bringt ihren eigenen Overhead mit sich – einen Kompromiss zwischen der Anzahl der benötigten Quantenberechnungen (Zeit) und der Genauigkeit unserer Ergebnisse:

$$\begin{array}{c|c|c|c} \text{Methods} & R=1 \text{, T-REx} & R=2 \text{, ZNE} \\[1mm] \hline \text{Assumptions} & \text{None} & \text{Ability to scale noise} \\[1mm] \text{Qubit overhead} & 1 & 1 \\[1mm] \text{Sampling overhead} & 2 & N_{\text{noise-factors}} \\[1mm] \text{Bias} & 0 & \mathcal{O}(\lambda^{N_{\text{noise-factors}}}) \\[1mm] \end{array}$$

Nutzung der Fehlerminderungs- und Fehlerunterdrückungsoptionen von Qiskit Runtime

So berechnet man einen Erwartungswert unter Verwendung von Fehlerminderung und Fehlerunterdrückung in Qiskit Runtime. Wir können genau denselben Circuit und dieselbe Observable wie zuvor verwenden, diesmal jedoch die Optimierungsstufe auf Stufe 2 festlegen und die Resilience-Stufe bzw. die verwendeten Fehlerminderungstechniken anpassen. Dieser Fehlerminderungsprozess findet mehrfach während einer Optimierungsschleife statt.

Diesen Teil führen wir auf echter Hardware durch, da Fehlerminderung auf Simulatoren nicht verfügbar ist.

# Estimated usage: 8 minutes, benchmarked on an Eagle processor, 5-23-24

from qiskit_ibm_runtime import QiskitRuntimeService
from qiskit_ibm_runtime import (
    Session,
    EstimatorOptions,
    EstimatorV2 as Estimator,
)

# We select the least busy backend

# Select the least busy backend
# backend = service.least_busy(
#    operational=True, min_num_qubits=ansatz.num_qubits, simulator=False
# )

# Or use a specific backend
backend = service.backend("ibm_brisbane")

# Initialize some variables to save the results from different runs:

exp_values_with_em0_es = []
exp_values_with_em1_es = []
exp_values_with_em2_es = []

# Use a pass manager to optimize the circuit and observables for the backend chosen:

pm = generate_preset_pass_manager(backend=backend, optimization_level=2)
isa_circuit = pm.run(circuit)
isa_observables = observables.apply_layout(layout=isa_circuit.layout)

# Open a session and run with no error mitigation:

estimator_options = EstimatorOptions(resilience_level=0, default_shots=10_000)

with Session(backend=backend) as session:
    estimator = Estimator(mode=session, options=estimator_options)

    pub = (isa_circuit, isa_observables, [individual_phases])
    cost = estimator.run([pub]).result()[0].data.evs

session.close()

exp_values_with_em0_es = cost[0]

# Open a session and run with resilience = 1:

estimator_options = EstimatorOptions(resilience_level=1, default_shots=10_000)

with Session(backend=backend) as session:
    estimator = Estimator(mode=session, options=estimator_options)

    pub = (isa_circuit, isa_observables, [individual_phases])
    cost = estimator.run([pub]).result()[0].data.evs

session.close()

exp_values_with_em1_es = cost[0]

# Open a session and run with resilience = 2:

estimator_options = EstimatorOptions(resilience_level=2, default_shots=10_000)

with Session(backend=backend) as session:
    estimator = Estimator(mode=session, options=estimator_options)

    pub = (isa_circuit, isa_observables, [individual_phases])
    cost = estimator.run([pub]).result()[0].data.evs

session.close()

exp_values_with_em2_es = cost[0]

Wie zuvor können wir die resultierenden Erwartungswerte als Funktion des Phasenwinkels für die drei verwendeten Fehlerminderungsstufen darstellen. Mit etwas Mühe lässt sich erkennen, dass die Fehlerminderung die Ergebnisse leicht verbessert. Auch dieser Effekt ist bei tieferen, komplexeren Circuits deutlich ausgeprägter.

import matplotlib.pyplot as plt

plt.plot(phases, exp_values_with_em0_es, "o", label="unmitigated")
plt.plot(phases, exp_values_with_em1_es, "o", label="resil = 1")
plt.plot(phases, exp_values_with_em2_es, "o", label="resil = 2")
plt.plot(phases, 2 * np.sin(phases / 2) ** 2 - 1, label="ideal")
plt.ylabel("Expectation")
plt.legend()
plt.show()

Zusammenfassung

In dieser Lektion hast du gelernt, wie du eine Kostenfunktion erstellst:

• Eine Kostenfunktion erstellen
• Qiskit Runtime Primitives zur Rauschminderung und -unterdrückung einsetzen
• Eine Messstrategie zur Optimierung von Geschwindigkeit und Genauigkeit definieren

Hier ist unser übergeordneter variationeller Arbeitsablauf:

Unsere Kostenfunktion wird in jeder Iteration der Optimierungsschleife ausgeführt. Die nächste Lektion zeigt, wie der klassische Optimierer die Auswertung unserer Kostenfunktion nutzt, um neue Parameter auszuwählen.

import qiskit
import qiskit_ibm_runtime

print(qiskit.version.get_version_info())
print(qiskit_ibm_runtime.version.get_version_info())

1.1.0
0.23.0