פונקציות עלות

במהלך שיעור זה, נלמד כיצד להעריך פונקציית עלות:

• ראשית, נלמד על פרימיטיבים של Qiskit Runtime
• נגדיר פונקציית עלות $C(\vec\theta)$. זו פונקציה ספציפית לבעיה שמגדירה את מטרת הבעיה עבור המייעל למזעור (או למיקסום)
• הגדרת אסטרטגיית מדידה עם פרימיטיבים של Qiskit Runtime לאיזון בין מהירות לדיוק

 

פרימיטיבים

כל המערכות הפיזיות, בין אם קלאסיות ובין אם קוונטיות, יכולות להתקיים במצבים שונים. לדוגמה, מכונית על כביש יכולה להיות בעלת מסה, מיקום, מהירות או תאוצה מסוימים המאפיינים את מצבה. באופן דומה, מערכות קוונטיות יכולות גם הן להיות במצורות או מצבים שונים, אך הן שונות ממערכות קלאסיות באופן שבו אנו מתמודדים עם מדידות ואבולוציה של מצב. זה מוביל לתכונות ייחודיות כמו סופרפוזיציה ושזירה שהן בלעדיות למכניקת הקוונטים. ממש כמו שאנו יכולים לתאר את מצב המכונית באמצעות תכונות פיזיקליות כמו מהירות או תאוצה, אנו יכולים גם לתאר את מצב המערכת הקוונטית באמצעות אובזרוובלים, שהם אובייקטים מתמטיים.

במכניקת הקוונטים, מצבים מיוצגים על ידי וקטורי עמודה מרוכבים מנורמלים, או קטים ($|\psi\rangle$), ואובזרוובלים הם אופרטורים לינאריים הרמיטיים ($\hat{H}=\hat{H}^{\dagger}$) הפועלים על הקטים. וקטור עצמי ($|\lambda\rangle$) של אובזרוובל מכונה מצב עצמי. מדידת אובזרוובל עבור אחד ממצביו העצמיים ($|\lambda\rangle$) תיתן לנו את הערך העצמי המתאים ($\lambda$) כפלט.

אם אתם תוהים כיצד למדוד מערכת קוונטית ומה ניתן למדוד, Qiskit מציע שני פרימיטיבים שיכולים לעזור:

• Sampler: בהינתן מצב קוונטי $|\psi\rangle$, פרימיטיב זה מחשב את ההסתברות של כל מצב בסיס חישובי אפשרי.
• Estimator: בהינתן אובזרוובל קוונטי $\hat{H}$ ומצב $|\psi\rangle$, פרימיטיב זה מחשב את הערך הצפוי של $\hat{H}$.

הפרימיטיב Sampler

הפרימיטיב Sampler מחשב את ההסתברות לקבל כל מצב אפשרי $|k\rangle$ מהבסיס החישובי, בהינתן מעגל קוונטי המכין את המצב $|\psi\rangle$. הוא מחשב

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

כאשר $n$ הוא מספר ה-Qubitים, ו-$k$ הוא הייצוג השלם של כל מחרוזת בינארית אפשרית $\{0,1\}^n$ (כלומר, מספרים שלמים בבסיס $2$).

ה-Sampler של Qiskit Runtime מריץ את המעגל מספר פעמים על מכשיר קוונטי, מבצע מדידות בכל ריצה, ומשחזר את התפלגות ההסתברות ממחרוזות הביטים שהתקבלו. ככל שיבצע יותר ריצות (או shots), התוצאות יהיו מדויקות יותר, אך הדבר דורש יותר זמן ומשאבים קוונטיים.

עם זאת, מכיוון שמספר הפלטים האפשריים גדל בצורה אקספוננציאלית עם מספר ה-Qubitים $n$ (כלומר, $2^n$), גם מספר ה-shots יצטרך לגדול בצורה אקספוננציאלית כדי לתפוס התפלגות הסתברות צפופה. לכן, Sampler יעיל רק עבור התפלגויות הסתברות דלילות; שבהן המצב היעד $|\psi\rangle$ חייב להיות ניתן לביטוי כצירוף לינארי של מצבי בסיס חישובי, כאשר מספר האיברים גדל לכל היותר פולינומית עם מספר ה-Qubitים:

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

ניתן גם להגדיר את Sampler לשליפת הסתברויות מחלק מהמעגל, המייצג תת-קבוצה של המצבים האפשריים הכוללים.

הפרימיטיב Estimator

הפרימיטיב Estimator מחשב את הערך הצפוי של אובזרוובל $\hat{H}$ עבור מצב קוונטי $|\psi\rangle$; שבו הסתברויות האובזרוובל ניתנות לביטוי כ-$p_\lambda = |\langle\lambda|\psi\rangle|^2$, כאשר $|\lambda\rangle$ הם המצבים העצמיים של האובזרוובל $\hat{H}$. הערך הצפוי מוגדר אז כממוצע של כל התוצאות האפשריות $\lambda$ (כלומר, הערכים העצמיים של האובזרוובל) של מדידת המצב $|\psi\rangle$, משוקלל לפי ההסתברויות המתאימות:

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

עם זאת, חישוב הערך הצפוי של אובזרוובל אינו תמיד אפשרי, שכן לעתים קרובות אנו לא יודעים את הבסיס העצמי שלו. ה-Estimator של Qiskit Runtime משתמש בתהליך אלגברי מורכב כדי להעריך את הערך הצפוי על מכשיר קוונטי אמיתי על ידי פירוק האובזרוובל לצירוף של אובזרוובלים אחרים שהבסיס העצמי שלהם ידוע לנו.

במילים פשוטות יותר, Estimator מפרק כל אובזרוובל שאינו יודע כיצד למדוד לאובזרוובלים פשוטים וניתנים למדידה הנקראים אופרטורי פאולי. כל אופרטור ניתן לביטוי כצירוף של $4^n$ אופרטורי פאולי.

$$\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\}, \\$$

כך ש-

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

כאשר $n$ הוא מספר ה-Qubitים, $k \equiv k_{n-1} \cdots k_0$ עבור $k_l \in \mathbb{Z}_4 \equiv \{0, 1, 2, 3\}$ (כלומר, מספרים שלמים בבסיס $4$), ו-$(\sigma_0, \sigma_1, \sigma_2, \sigma_3) := (I, X, Y, Z)$.

לאחר ביצוע פירוק זה, Estimator גוזר מעגל חדש $V_k|\psi\rangle$ עבור כל אובזרוובל $\hat{P}_k$ (מהמעגל המקורי), כדי לאלכסן בפועל את אובזרוובל פאולי בבסיס החישובי ולמדוד אותו. אנו יכולים למדוד בקלות אובזרוובלי פאולי מכיוון שאנו יודעים את $V_k$ מראש, מה שאינו המקרה בדרך כלל עבור אובזרוובלים אחרים.

עבור כל $\hat{P}_{k}$, ה-Estimator מריץ את המעגל המתאים על מכשיר קוונטי מספר פעמים, מודד את מצב הפלט בבסיס החישובי, ומחשב את ההסתברות $p_{kj}$ לקבל כל פלט אפשרי $j$. לאחר מכן הוא מחפש את הערך העצמי $\lambda_{kj}$ של $P_k$ המתאים לכל פלט $j$, מכפיל ב-$w_k$, ומסכם את כל התוצאות כדי לקבל את הערך הצפוי של האובזרוובל $\hat{H}$ עבור המצב הנתון $|\psi\rangle$.

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

מכיוון שחישוב הערך הצפוי של $4^n$ פאוליים אינו מעשי (כלומר, גדל אקספוננציאלית), Estimator יכול להיות יעיל רק כאשר כמות גדולה של $w_k$ שווה לאפס (כלומר, פירוק פאולי דליל במקום צפוף). באופן פורמלי אנו אומרים שכדי שחישוב זה יהיה ניתן לפתרון ביעילות, מספר האיברים השונים מאפס חייב לגדול לכל היותר פולינומית עם מספר ה-Qubitים $n$: $\hat{H} = \sum^{\text{Poly}(n)}_k w_k \hat{P}_k.$

הקורא עשוי להבחין בהנחה הגלומה שגם דגימת ההסתברות צריכה להיות יעילה כפי שהוסבר עבור Sampler, כלומר

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

דוגמה מודרכת לחישוב ערכים צפויים

נניח את המצב החד-קיוביטי $|+\rangle := H|0\rangle = \frac{1}{\sqrt{2}}(|0\rangle + |1\rangle)$, ואובזרוובל

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

עם הערך הצפוי התיאורטי הבא $\langle\hat{H}\rangle_+ = \langle+|\hat{H}|+\rangle = 2.$

מכיוון שאנו לא יודעים כיצד למדוד אובזרוובל זה, לא ניתן לחשב את ערכו הצפוי ישירות, ועלינו לבטא אותו מחדש כ-$\langle\hat{H}\rangle_+ = 2\langle X \rangle_+ - \langle Z \rangle_+$. ניתן להראות שזה מוביל לאותה תוצאה בשל העובדה ש-$\langle+|X|+\rangle = 1$, ו-$\langle+|Z|+\rangle = 0$.

נראה כיצד לחשב את $\langle X \rangle_+$ ו-$\langle Z \rangle_+$ ישירות. מכיוון ש-$X$ ו-$Z$ אינם מתחלפים (כלומר, אין להם אותו בסיס עצמי), לא ניתן למדוד אותם בו-זמנית, ולכן אנחנו צריכים את המעגלים העזר:

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

כעת נוכל לבצע את החישוב ידנית באמצעות Sampler ולאמת את התוצאות עם Estimator:

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

דיוק מתמטי (אופציונלי)

כאשר מבטאים את $|\psi\rangle$ ביחס לבסיס המצבים העצמיים של $\hat{H}$, $|\psi\rangle = \sum_\lambda a_\lambda |\lambda\rangle$, נובע:

$$\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}$$

מכיוון שאיננו יודעים את הערכים העצמיים או המצבים העצמיים של האובזרוובל היעד $\hat{H}$, ראשית עלינו לשקול את האלכסון שלו. בהינתן ש-$\hat{H}$ הוא הרמיטי, קיים טרנספורמציה אוניטרית $V$ כך ש-$\hat{H}=V^\dagger \Lambda V,$ כאשר $\Lambda$ היא מטריצת הערכים העצמיים האלכסונית, כך ש-$\langle j | \Lambda | k \rangle = 0$ אם $j\neq k$, ו-$\langle j | \Lambda | j \rangle = \lambda_j$.

משמע שהערך הצפוי ניתן לכתיבה מחדש כ:

$$\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}$$

בהינתן שאם מערכת נמצאת במצב $|\phi\rangle = V |\psi\rangle$ ההסתברות למדוד $| j\rangle$ היא $p_j = |\langle j|\phi \rangle|^2$, ניתן לבטא את הערך הצפוי הנ"ל כ:

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

חשוב מאוד לציין שההסתברויות נלקחות מהמצב $V |\psi\rangle$ ולא מ-$|\psi\rangle$. זו הסיבה שהמטריצה $V$ הכרחית לחלוטין. אולי תשאלו כיצד לקבל את המטריצה $V$ ואת הערכים העצמיים $\Lambda$. אם כבר היו לכם הערכים העצמיים, לא היה צורך להשתמש במחשב קוונטי שכן מטרת האלגוריתמים הווריאציוניים היא למצוא את הערכים העצמיים הללו של $\hat{H}$.

למרבה המזל, קיים פתרון: כל מטריצה $2^n \times 2^n$ ניתנת לכתיבה כצירוף לינארי של $4^n$ מכפלות טנזוריות של $n$ מטריצות פאולי וזהויות, כולן הרמיטיות ואוניטריות עם $V$ ו-$\Lambda$ ידועים. זה מה שה-Estimator של Runtime עושה פנימית על ידי פירוק כל אובייקט Operator ל-SparsePauliOp.

הנה האופרטורים שניתן להשתמש בהם:

$$\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}$$

כעת נכתוב מחדש את $\hat{H}$ ביחס לפאוליים ולזהויות:

$$\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,$$

כאשר $k = \sum_{l=0}^{n-1} 4^l k_l \equiv k_{n-1}...k_0$ עבור $k_{n-1},...,k_0\in \{0,1,2,3\}$ (כלומר, בסיס $4$), ו-$\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}$$

כאשר $V_k := V_{k_{n-1}}\otimes ... \otimes V_{k_0}$ ו-$\Lambda_k := \Lambda_{k_{n-1}}\otimes ... \otimes \Lambda_{k_0}$, כך ש-$\hat{P_k}=V_k^\dagger \Lambda_k V_k.$

פונקציות עלות

באופן כללי, פונקציות עלות משמשות לתיאור המטרה של בעיה ועד כמה מצב ניסיון מסוים מתפקד ביחס למטרה זו. ניתן ליישם הגדרה זו על דוגמאות שונות בכימיה, למידת מכונה, פיננסים, אופטימיזציה, וכן הלאה.

בוא נבחן דוגמה פשוטה של מציאת מצב הבסיס של מערכת. המטרה שלנו היא למזער את ערך הציפייה של האובייקט המייצג אנרגיה (ה-Hamiltonian $\hat{\mathcal{H}}$):

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

נוכל להשתמש ב-Estimator כדי להעריך את ערך הציפייה ולהעביר ערך זה לאופטימייזר לצורך מינימיזציה. אם האופטימיזציה תצליח, היא תחזיר קבוצה של ערכי פרמטרים אופטימליים $\vec\theta^*$, שמהם נוכל לבנות את מצב הפתרון המוצע $|\psi(\vec\theta^*)\rangle$ ולחשב את ערך הציפייה הנצפה כ-$C(\vec\theta^*)$.

שים לב שנוכל למזער את פונקציית העלות רק עבור קבוצת המצבים המוגבלת שאנחנו שוקלים. זה מוביל אותנו לשתי אפשרויות נפרדות:

• ה-ansatz שלנו לא מגדיר את מצב הפתרון על פני מרחב החיפוש: אם זה המצב, האופטימייזר שלנו לעולם לא ימצא את הפתרון, ועלינו לנסות ansatzes אחרים שעשויים לייצג את מרחב החיפוש שלנו בצורה מדויקת יותר.
• האופטימייזר שלנו אינו מסוגל למצוא את הפתרון הזה: ניתן להגדיר אופטימיזציה בצורה גלובלית ובצורה מקומית. נחקור את המשמעות של זה בסעיף מאוחר יותר.

בסך הכל, נבצע לולאת אופטימיזציה קלאסית אך נסמוך על הערכת פונקציית העלות על ידי מחשב קוונטי. מנקודת מבט זו, ניתן לחשוב על האופטימיזציה כמאמץ קלאסי לחלוטין שבו אנחנו קוראים לאורקל קוונטי של קופסה שחורה בכל פעם שהאופטימייזר צריך להעריך את פונקציית העלות.

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

נבצע זאת תחילה באמצעות סימולטור: ה-StatevectorEstimator. זה בדרך כלל מומלץ לצורכי דיבאגינג, אך נעקוב מיד אחרי ריצת הדיבאגינג עם חישוב על חומרת קוונטי אמיתית. יותר ויותר, בעיות שמעניינות אותנו אינן ניתנות עוד לסימולציה קלאסית ללא מתקני מחשוב-על ברמה גבוהה.

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

[-0.58744589]

כעת נמשיך עם הרצה על מחשב קוונטי אמיתי. שים לב לשינויים בתחביר. הצעדים הכוללים את ה-pass_manager יידונו בהמשך בדוגמה הבאה. צעד בעל חשיבות מיוחדת באלגוריתמים וריאציוניים הוא השימוש ב-Qiskit Runtime Session. פתיחת Session מאפשרת לך להריץ איטרציות מרובות של אלגוריתם וריאציוני מבלי להמתין בתור חדש בכל פעם שהפרמטרים מתעדכנים. זה חשוב אם זמני התור ארוכים ו/או נדרשות איטרציות רבות. רק שותפים ברשת IBM Quantum® יכולים להשתמש ב-Runtime Sessions. אם אין לך גישה ל-Sessions, תוכל להפחית את מספר האיטרציות שאתה שולח בכל פעם, ולשמור את הפרמטרים האחרונים לשימוש בהרצות עתידיות. אם תשלח יותר מדי איטרציות או תיתקל בזמני תור ארוכים מדי, ייתכן שתיתקל בקוד שגיאה 1217, המתייחס לעיכובים ארוכים בין שליחות עבודות.

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

שים לב שהערכים שהתקבלו משני החישובים לעיל דומים מאוד. טכניקות לשיפור התוצאות יידונו בהמשך.

דוגמה: מיפוי למערכות לא-פיזיקליות

בעיית החיתוך המקסימלי (Max-Cut) היא בעיית אופטימיזציה קומבינטורית שכוללת חלוקת הצמתים של גרף לשתי קבוצות זרות כך שמספר הקשתות בין שתי הקבוצות יהיה מקסימלי. באופן פורמלי יותר, בהינתן גרף לא מכוון $G=(V,E)$, כאשר $V$ היא קבוצת הצמתים ו-$E$ היא קבוצת הקשתות, בעיית ה-Max-Cut מבקשת לחלק את הצמתים לשתי תתי-קבוצות זרות, $S$ ו-$T$, כך שמספר הקשתות עם נקודת קצה אחת ב-$S$ והשנייה ב-$T$ יהיה מקסימלי.

ניתן ליישם Max-Cut לפתרון בעיות שונות כולל: קיבוץ, עיצוב רשתות, מעברי פאזה, וכן הלאה. נתחיל ביצירת גרף הבעיה:

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

ניתן לבטא את הבעיה הזו כבעיית אופטימיזציה בינארית. עבור כל צומת $0 \leq i < n$, כאשר $n$ הוא מספר הצמתים בגרף (במקרה זה $n=4$), נשקול את המשתנה הבינארי $x_i$. למשתנה זה יהיה הערך $1$ אם צומת $i$ נמצאת באחת מהקבוצות שנסמן כ-$1$ ו-$0$ אם היא בקבוצה האחרת, שנסמנה כ-$0$. נסמן גם כ-$w_{ij}$ (האיבר $(i,j)$ של מטריצת הסמיכות $w$) את משקל הקשת שהולכת מצומת $i$ לצומת $j$. מכיוון שהגרף אינו מכוון, $w_{ij}=w_{ji}$. אז נוכל לנסח את הבעיה שלנו כמקסום של פונקציית העלות הבאה:

$$\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}$$

כדי לפתור את הבעיה הזו עם מחשב קוונטי, אנחנו הולכים לבטא את פונקציית העלות כערך הצפוי של אובייקט. עם זאת, האובייקטים שקיסקיט תומך בהם באופן מובנה מורכבים מאופרטורי Pauli, שיש להם ערכים עצמיים $1$ ו-$-1$ במקום $0$ ו-$1$. לכן נבצע את שינוי המשתנה הבא:

כאשר $\vec{x}=(x_0,x_1,\cdots ,x_{n-1})$. נוכל להשתמש במטריצת הסמיכות $w$ כדי לגשת בנוחות למשקלי כל הקשתות. זה ישמש לקבלת פונקציית העלות שלנו:

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

זה גורר:

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

אז פונקציית העלות החדשה שאנחנו רוצים למקסם היא:

$$\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}$$

יתרה מכך, הנטייה הטבעית של מחשב קוונטי היא למצוא מינימום (בדרך כלל האנרגיה הנמוכה ביותר) במקום מקסימום, כך שבמקום למקסם את $C(\vec{z})$ אנחנו הולכים למזער:

$$-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}$$

עכשיו שיש לנו פונקציית עלות למזעור שהמשתנים שלה יכולים לקבל ערכים $-1$ ו-$1$, נוכל לעשות את האנלוגיה הבאה עם ה-Pauli $Z$:

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

במילים אחרות, המשתנה $z_i$ יהיה שקול ל-Gate $Z$ הפועל על Qubit $i$. יתרה מכך:

$$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$$

אז האובייקט שאנחנו הולכים לשקול הוא:

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

שאליו נצטרך להוסיף את האיבר הבלתי תלוי לאחר מכן:

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

האופרטור הוא צירוף לינארי של איברים עם אופרטורי Z על צמתים המחוברים על ידי קשת (זכור שה-Qubit ה-0 הוא הכי ימני): $IIZZ + IZIZ + IZZI + ZIIZ + ZZII$. לאחר שהאופרטור נבנה, ניתן לבנות בקלות את ה-ansatz לאלגוריתם QAOA באמצעות Circuit QAOAAnsatz מספריית ה-Circuit של Qiskit.

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

כאשר Estimator של Runtime מקבל ישירות Hamiltonian ו-ansatz עם פרמטרים ומחזיר את האנרגיה הדרושה, פונקציית העלות עבור מופע QAOA היא פשוטה למדי:

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

נחזור לדוגמה זו ב-Applications כדי לחקור כיצד ממנפים אופטימייזר כדי לעבור על מרחב החיפוש. באופן כללי, זה כולל:

• מינוף אופטימייזר למציאת פרמטרים אופטימליים
• קישור פרמטרים אופטימליים ל-ansatz כדי למצוא ערכים עצמיים
• תרגום הערכים העצמיים להגדרת הבעיה שלנו

אסטרטגיית מדידה: מהירות מול דיוק

כפי שצוין, אנחנו משתמשים במחשב קוונטי רועש כ-אורקל קופסה שחורה, כאשר הרעש עלול להפוך את הערכים המתקבלים ללא-דטרמיניסטיים, מה שמוביל לתנודות אקראיות שבתורן יפגעו — או אפילו ימנעו לחלוטין — בהתכנסות של מיטביים מסוימים לפתרון מוצע. זהו בעיה כללית שעלינו לטפל בה בעת שאנו חוקרים בהדרגה את השימושיות הקוונטית ומתקדמים לעבר יתרון קוונטי:

אנחנו יכולים להשתמש באפשרויות דיכוי שגיאות והפחתת שגיאות של Qiskit Runtime Primitive כדי לטפל ברעש ולמקסם את שימושיות המחשבים הקוונטיים של ימינו.

דיכוי שגיאות

דיכוי שגיאות מתייחס לטכניקות המשמשות לאופטימיזציה ושינוי Circuit במהלך הקומפילציה כדי למזער שגיאות. זוהי טכניקת טיפול בסיסית בשגיאות שבדרך כלל גורמת לקצת עיבוד מקדים קלאסי תקורה לזמן הריצה הכולל. התקורה כוללת טרנספילציה של Circuits להרצה על חומרת קוונטום על ידי:

• ביטוי ה-Circuit באמצעות ה-Gates הנייטיביים הזמינים במערכת קוונטית
• מיפוי ה-Qubits הוירטואליים ל-Qubits פיזיים
• הוספת SWAPs בהתאם לדרישות הקישוריות
• אופטימיזציה של Gates חד-Qubit ודו-Qubit
• הוספת decoupling דינמי ל-Qubits סרלים כדי למנוע השפעות של דה-קוהרנטיות.

Primitives מאפשרים שימוש בטכניקות דיכוי שגיאות על ידי הגדרת האפשרות optimization_level ובחירת אפשרויות טרנספילציה מתקדמות. בקורס מאוחר יותר, נעמיק בשיטות בנייה שונות של Circuit כדי לשפר תוצאות, אבל ברוב המקרים, אנחנו ממליצים להגדיר optimization_level=3.

נמחיש את הערך של הגברת האופטימיזציה בתהליך הטרנספילציה על ידי בחינת Circuit לדוגמה עם התנהגות אידיאלית פשוטה.

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

ה-Circuit לעיל יכול לתת ערכי ציפייה סינוסואידיים של ה-observable הנתון, בתנאי שנכניס פאזות שמכסות טווח מתאים, כגון $[0,2\pi]$.

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

אנחנו יכולים להשתמש בסימולטור כדי להראות את השימושיות של טרנספילציה מאופטמת. נחזור למטה לשימוש בחומרה אמיתית כדי להדגים את שימושיות הפחתת השגיאות. נשתמש ב-QiskitRuntimeService כדי לקבל Backend אמיתי (במקרה זה, ibm_brisbane), ונשתמש ב-AerSimulator לסימולציה של אותו Backend, כולל התנהגות הרעש שלו.

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)

כעת אנחנו יכולים להשתמש ב-pass manager כדי לטרנספל את ה-Circuit ל-"instruction set architecture" או ISA של ה-Backend. זוהי דרישה חדשה ב-Qiskit Runtime: כל ה-Circuits המוגשים ל-Backend חייבים לעמוד במגבלות של יעד ה-Backend, כלומר הם חייבים להיות כתובים במונחים של ISA של ה-Backend — כלומר, מערכת ההוראות שהמכשיר יכול להבין ולבצע. אילוצי היעד הללו מוגדרים על ידי גורמים כמו ה-Gates הנייטיביים של המכשיר, קישוריות ה-Qubit שלו, ו — כשרלוונטי — מפרטי תזמון הדחף וההוראות האחרות שלו.

שים לב שבמקרה הנוכחי, נעשה זאת פעמיים: פעם אחת עם optimization_level = 0, ופעם אחת כשהוא מוגדר ל-3. בכל פעם נשתמש ב-Estimator Primitive כדי לאמוד את ערכי הציפייה של ה-observable בערכי פאזה שונים.

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

לבסוף, אנחנו יכולים לשרטט את התוצאות, ואנחנו רואים שדיוק החישוב היה טוב למדי אפילו ללא אופטימיזציה, אבל הוא בהחלט השתפר עם הגברת האופטימיזציה לרמה 3. שים לב שב-Circuits עמוקים ומורכבים יותר, ההבדל בין רמות אופטימיזציה 0 ו-3 צפוי להיות משמעותי יותר. זהו Circuit פשוט מאוד המשמש כמודל צעצוע.

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

הפחתת שגיאות

הפחתת שגיאות מתייחסת לטכניקות המאפשרות למשתמשים להפחית שגיאות Circuit על ידי מידול רעש המכשיר בזמן הביצוע. בדרך כלל, זה גורם לתקורת עיבוד מקדים קוונטי הקשורה לאימון המודל ותקורת עיבוד לאחר קלאסי להפחתת שגיאות בתוצאות הגולמיות באמצעות המודל שנוצר.

האפשרות resilience_level של Qiskit Runtime Primitive מציינת את כמות החוסן שיש לבנות נגד שגיאות. רמות גבוהות יותר מייצרות תוצאות מדויקות יותר במחיר זמני עיבוד ארוכים יותר עקב תקורת דגימה קוונטית. ניתן להשתמש ברמות החוסן כדי להגדיר את הפשרה בין עלות לדיוק בעת החלת הפחתת שגיאות על השאילתה ל-Primitive שלך.

כאשר מיישמים כל טכניקת הפחתת שגיאות, אנחנו מצפים שה-הטיה בתוצאות שלנו תצטמצם ביחס להטיה הקודמת, הבלתי-מופחתת. במקרים מסוימים, ההטיה עשויה אפילו להיעלם. אולם, זה כרוך בעלות. ככל שאנחנו מפחיתים את ההטיה בכמויות המשוערות שלנו, השונות הסטטיסטית תגדל (כלומר, שונות), שאותה אנחנו יכולים לפצות על ידי הגדלה נוספת של מספר ה-shots לכל Circuit בתהליך הדגימה שלנו. זה יכניס תקורה מעבר לנדרש להפחתת ההטיה, לכן זה לא נעשה כברירת מחדל. אנחנו יכולים להצטרף להתנהגות זו בקלות על ידי התאמת מספר ה-shots לכל Circuit ב-options.executions.shots, כפי שמוצג בדוגמה להלן.

לקורס זה, נחקור מודלי הפחתת שגיאות אלו ברמה גבוהה כדי להמחיש את הפחתת השגיאות שה-Primitives של Qiskit Runtime יכולים לבצע מבלי לדרוש פרטי מימוש מלאים.

כיבוי שגיאות קריאה מסולסל (T-REx)

כיבוי שגיאות קריאה מסולסל (T-REx) משתמש בטכניקה הידועה כ-Pauli twirling כדי להפחית את הרעש המוכנס במהלך תהליך המדידה הקוונטית. טכניקה זו אינה מניחה שום צורה ספציפית של רעש, מה שהופך אותה לכללית ויעילה מאוד.

זרימת עבודה כללית:

1. רכישת נתונים עבור מצב האפס עם היפוכי ביט אקראיים (Pauli X לפני המדידה)
2. רכישת נתונים עבור המצב הרצוי (הרועש) עם היפוכי ביט אקראיים (Pauli X לפני המדידה)
3. חישוב הפונקציה המיוחדת עבור כל קבוצת נתונים, וחלוקה.

 

אנחנו יכולים להגדיר זאת עם options.resilience_level = 1, כפי שמוצג בדוגמה להלן.

אקסטרפולציה של אפס רעש

אקסטרפולציה של אפס רעש (ZNE) פועלת על ידי הגברה תחילה של הרעש ב-Circuit המכין את המצב הקוונטי הרצוי, קבלת מדידות עבור מספר רמות שונות של רעש, ושימוש במדידות אלו להסקת התוצאה נטולת הרעש.

זרימת עבודה כללית:

1. הגבר את רעש ה-Circuit עבור מספר גורמי רעש
2. הרץ כל Circuit מוגבר רעש
3. אקסטרפל בחזרה לגבול אפס הרעש

 

אנחנו יכולים להגדיר זאת עם options.resilience_level = 2. נוכל לאופטמז זאת עוד יותר על ידי חקירה של מגוון noise_factors, noise_amplifiers ו-extrapolators, אך זה מחוץ לתחום הקורס הזה. אנחנו מעודדים אותך להתנסות עם אפשרויות אלו כמתואר כאן.

כל שיטה מגיעה עם תקורה משלה: פשרה בין מספר החישובים הקוונטיים הנדרשים (זמן) לבין דיוק התוצאות שלנו:

$$\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}$$

שימוש באפשרויות הפחתה ודיכוי של Qiskit Runtime

כך מחשבים ערך ציפייה תוך שימוש בהפחתת שגיאות ודיכוי ב-Qiskit Runtime. אנחנו יכולים להשתמש בדיוק באותו Circuit ואותו observable כמקודם, אבל הפעם שומרים על רמת האופטימיזציה קבועה ברמה 2, ועכשיו מכווננים את ה-resilience או טכניקות הפחתת השגיאות בשימוש. תהליך הפחתת השגיאות הזה מתרחש מספר פעמים לאורך לולאת האופטימיזציה.

אנחנו מבצעים חלק זה על חומרה אמיתית, שכן הפחתת שגיאות אינה זמינה על סימולטורים.

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

כמקודם, אנחנו יכולים לשרטט את ערכי הציפייה המתקבלים כפונקציה של זווית הפאזה עבור שלוש רמות הפחתת השגיאות שבהן השתמשנו. בקושי רב ניתן לראות שהפחתת השגיאות משפרת את התוצאות מעט. שוב, אפקט זה הרבה יותר בולט ב-Circuits עמוקים ומורכבים יותר.

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

סיכום

בשיעור זה, למדת כיצד ליצור פונקציית עלות:

• יצירת פונקציית עלות
• כיצד למנף את ה-Primitives של Qiskit Runtime להפחתת ודיכוי רעש
• כיצד להגדיר אסטרטגיית מדידה לאופטימיזציה של מהירות מול דיוק

הנה עומס העבודה הווריאציוני ברמה גבוהה שלנו:

פונקציית העלות שלנו רצה בכל איטרציה של לולאת האופטימיזציה. השיעור הבא יחקור כיצד המיטב הקלאסי משתמש בהערכת פונקציית העלות שלנו כדי לבחור פרמטרים חדשים.

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