אומדן אנרגיית מצב-היסוד של שרשרת היזנברג עם VQE

אומדן שימוש: 37 דקות על מעבד Heron (הערה: זהו אומדן בלבד. זמן הריצה עשוי להשתנות.)

תוצאות למידה

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

• כיצד להשתמש ב-optimizer SPSA לאומדן אנרגיית מצב-היסוד של מערכת קוונטית

• כיצד להריץ תהליכי עבודה וריאציוניים על חומרת קוונטום של IBM® באמצעות IBM Quantum primitives וסשנים

דרישות מוקדמות

• יסודות מידע קוונטי

• מבוא לתבניות Qiskit

• עיצוב אלגוריתם וריאציוני

רקע

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

$H = \sum_{\langle i,j \rangle} \left( J_x X_i X_j + J_y Y_i Y_j + J_z Z_i Z_j \right) + \sum_{i} h_i Z_i,$

כאשר $X_i$, $Y_i$, ו-$Z_i$ הם אופרטורי פאולי הפועלים על אתר $i$, הסכום $\langle i,j \rangle$ רץ על פני זוגות שכנים קרובים, $J_x = J_y = J_z = 0.5$ הם קבועי הצימוד החילופין (איזוטרופי במדריך זה), ו-$h_i$ מייצג שדה מגנטי חיצוני תלוי-אתר. במדריך זה, ערכי השדה המגנטי נדגמים אקראית מהקטע $[-1, 1]$. שימו לב שבמימוש שלהלן, קבוצת זוגות "שכנים קרובים" נקבעת על פי הצימוד הטבעי של ה-backend לחומרה בין $N$ ה-qubits הראשונים, שאינם בהכרח מהווים שרשרת לינארית קפדנית בהתאם לטופולוגיית המכשיר.

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

ה-Variational Quantum Eigensolver (VQE) הוא אלגוריתם היברידי קוונטי-קלאסי המיועד לאומדן אנרגיית מצב-היסוד של המילטוניאן. הוא פועל על ידי הכנת מצב קוונטי פרמטרי $|\psi(\theta)\rangle$ (הנקרא ansatz) על מחשב קוונטי ומדידת ערך הציפייה $\langle \psi(\theta) | H | \psi(\theta) \rangle$. optimizer קלאסי מכוונן אז איטרטיבית את הפרמטרים $\theta$ כדי למזעֵר אנרגיה זו, מנצל את עקרון הוריאציה המבטיח שהאנרגיה הנמדדת היא תמיד חסם עליון לאנרגיית מצב-היסוד האמיתית.

במדריך זה, אנו משתמשים ב-ansatz ‏efficient_su2 מספריית המעגלים של Qiskit, אשר בונה שכבות של סיבובי qubit יחיד ושערי שזירה. האופטימיזציה מתבצעת באמצעות אלגוריתם Simultaneous Perturbation Stochastic Approximation (SPSA), המתאים במיוחד לחומרה קוונטית רועשת מפני שהוא מאמד גרדיאנטים תוך שימוש בשתי הערכות פונקציה בלבד לכל איטרציה ללא תלות במספר הפרמטרים.

דרישות

לפני תחילת מדריך זה, ודאו שיש לכם את הדברים הבאים מותקנים:

• Qiskit SDK v2.0 או מאוחר יותר, עם תמיכה בהדמייה

• Qiskit Runtime v0.44 או מאוחר יותר (pip install qiskit-ibm-runtime)

הגדרה

# Added by doQumentation — installs the packages this notebook needs if they are missing
import importlib.util

_needed = {"matplotlib": "matplotlib", "numpy": "numpy", "qiskit": "qiskit", "qiskit_ibm_runtime": "qiskit-ibm-runtime"}
_missing = [pip for module, pip in _needed.items()
            if importlib.util.find_spec(module) is None]
# One at a time, so a package that fails to install does not block the others
for _pip in _missing:
    %pip install -q {_pip}
if not _missing:
    print("\u2713 All required packages are installed")

import numpy as np
import matplotlib.pyplot as plt
from typing import Sequence

from qiskit import QuantumCircuit
from qiskit.quantum_info import SparsePauliOp
from qiskit.primitives import BaseEstimatorV2
from qiskit.circuit.library import XGate
from qiskit.circuit.library import efficient_su2
from qiskit.transpiler import PassManager
from qiskit.transpiler.preset_passmanagers import generate_preset_pass_manager
from qiskit.transpiler.passes.scheduling import (
    ALAPScheduleAnalysis,
    PadDynamicalDecoupling,
)
from qiskit_ibm_runtime import QiskitRuntimeService, Session, EstimatorV2

def visualize_results(results):
    plt.plot(results["cost_history"], lw=2)
    plt.xlabel("Number of function evaluations")
    plt.ylabel("Energy")
    plt.show()

דוגמה בקנה מידה קטן

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

שלב 1: מיפוי קלטים קלאסיים לבעיה קוונטית

• קלט: מספר הספינים

• פלט: Ansatz והמילטוניאן הממדלים את שרשרת היזנברג

נבנה ansatz והמילטוניאן הממדלים שרשרת היזנברג בעלת 10 ספינים. בשלב זה, נבנה המילטוניאן של היזנברג בעל 10 ספינים על מפת הצימוד של ה-backend הפחות עמוס ונכין את ה-ansatz ‏efficient_su2.

num_spins = 10
ansatz = efficient_su2(num_qubits=num_spins, reps=2)

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

coupling = backend.target.build_coupling_map()
reduced_coupling = coupling.reduce(list(range(num_spins)))

edge_list = reduced_coupling.graph.edge_list()
ham_list = []

for edge in edge_list:
    ham_list.append(("ZZ", edge, 0.5))
    ham_list.append(("YY", edge, 0.5))
    ham_list.append(("XX", edge, 0.5))

for qubit in reduced_coupling.physical_qubits:
    ham_list.append(("Z", [qubit], np.random.random() * 2 - 1))

hamiltonian = SparsePauliOp.from_sparse_list(ham_list, num_qubits=num_spins)

ansatz.draw("mpl", style="iqp")

שלב 2: אופטימיזציה של הבעיה להרצה על חומרה קוונטית

• קלט: מעגל מופשט, observable

• פלט: מעגל יעד ו-observable, מותאמים עבור ה-QPU שנבחר

השתמשו בפונקציה generate_preset_pass_manager מ-Qiskit כדי ליצור אוטומטית שגרת אופטימיזציה עבור המעגל שלנו ביחס ל-QPU שנבחר. אנחנו בוחרים optimization_level=3, המספק את הרמה הגבוהה ביותר של אופטימיזציה מבין מנהלי ההעברות הקבועים מראש. אנחנו גם כוללים העברות תזמון ALAPScheduleAnalysis ו-PadDynamicalDecoupling כדי לדכא שגיאות דה-קוהרנטיות.

target = backend.target
pm = generate_preset_pass_manager(optimization_level=3, target=target)
pm.scheduling = PassManager(
    [
        ALAPScheduleAnalysis(durations=target.durations()),
        PadDynamicalDecoupling(
            durations=target.durations(),
            dd_sequence=[XGate(), XGate()],
            pulse_alignment=target.pulse_alignment,
        ),
    ]
)
isa_ansatz = pm.run(ansatz)
isa_observable = hamiltonian.apply_layout(isa_ansatz.layout)
isa_ansatz.draw("mpl", scale=0.6, style="iqp", fold=-1, idle_wires=False)

שלב 3: הרצה באמצעות primitives של Qiskit

• קלט: מעגל יעד ו-observable

• פלט: תוצאות האופטימיזציה

מזעור אנרגיית מצב-היסוד המוערכת של המערכת על ידי אופטימיזציה של פרמטרי המעגל. השתמשו ב-IBM Quantum primitive Estimator כדי להעריך את פונקציית העלות במהלך האופטימיזציה.

מאחר שאופטמנו את המעגל עבור ה-backend בשלב 2, ניתן להימנע מטרנספילציה בשרת Quantum Compute על ידי הגדרת skip_transpilation=True והעברת המעגל המאופטם. להדגמה זו, נריץ על QPU באמצעות IBM Quantum primitives. כדי להריץ עם primitives מבוססי statevector של qiskit, החליפו את בלוק הקוד המשתמש ב-IBM Quantum primitives בבלוק שבהערה.

במדריך זה אנו משתמשים ב-Simultaneous Perturbation Stochastic Approximation (SPSA), שהוא optimizer מבוסס גרדיאנט. להלן הצגה קצרה שלו, וקוד לממש SPSA עם Qiskit v2.0.

הצגת SPSA

Simultaneous Perturbation Stochastic Approximation (SPSA) [1] הוא אלגוריתם אופטימיזציה המקרב את וקטור הגרדיאנט כולו בשימוש בשתי קריאות פונקציה בלבד בכל איטרציה. תהי $f:\mathbb{R}^p\rightarrow \mathbb{R}$ פונקציית העלות עם $p$ פרמטרים לאופטימיזציה, ו-$x_i\in \mathbb{R}^p$ וקטור הפרמטרים בשלב ה-$i^{th}$ של האיטרציה. לחישוב הגרדיאנט, נוצר וקטור אקראי $\Delta_i$ בגודל $p$, כאשר כל איבר $\Delta_{ij}$, $\forall$ $j\in \{1,2,...,p\}$, נדגם אחיד מ-$\{-1, 1\}$. לאחר מכן, כל איבר של הוקטור האקראי $\Delta_i$ מוכפל בערך קטן $c_i$ כדי ליצור הפרעה אקראית. הגרדיאנט מוערך אז כ-

$[\nabla f(x_i)]_j \approx \frac{f(x_i + c_i \Delta_i) - f(x_i - c_i \Delta_i)}{2c_i\Delta_{ij}}.$

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

במדריך זה, ה-hyperparameters עבור האיטרציה ה-$i^{th}$, $a_i$ ו-$c_i$, מחושבים כ-

$a_i = \frac{a}{(A + i + 1)^\alpha} \quad \text{and} \quad c_i = \frac{c}{(i+1)^\gamma},$

כאשר הערכים הקבועים הם $A = 30$, $\alpha = 0.9$, $a = 0.3$, $c = 0.1$, ו-$\gamma = 0.4$. ערכים אלה נבחרו מ-[2]. כיוון נאות של hyperparameters הכרחי להפקת ביצועים טובים מ-SPSA.

def spsa(
    fun, x0, args=(), A=30, alpha=0.9, a=0.3, c=0.1, gamma=0.4, maxiter=100
):
    nparams = len(x0)
    x = np.copy(x0)

    for i in range(maxiter):
        a_i = a / (A + i + 1) ** alpha
        c_i = c / (i + 1) ** gamma
        delta_i = np.random.choice([-1, 1], nparams)

        # two hardware calls
        eval_1 = fun(x + c_i * delta_i, *args)
        eval_2 = fun(x - c_i * delta_i, *args)

        # compute the gradient and update the parameters
        grad = (eval_1 - eval_2) / (2 * c_i) * np.reciprocal(delta_i)
        x = x - a_i * grad

    return x

def cost_func(
    params: Sequence,
    ansatz: QuantumCircuit,
    hamiltonian: SparsePauliOp,
    estimator: BaseEstimatorV2,
    cost_history_dict: dict,
) -> float:
    """Ground state energy evaluation."""
    energy = (
        estimator.run([(ansatz, hamiltonian, [params])]).result()[0].data.evs
    )

    cost_history_dict["iters"] += 1
    cost_history_dict["prev_vector"] = list(params)
    cost_history_dict["cost_history"].append(float(energy[0]))

    print(
        f"Fx Iters. done: {cost_history_dict['iters']} [Current cost: {round(energy[0], 5)}]",
        end="\r",
    )

    return energy

def solve(x0, isa_ansatz, isa_observable, maxiter=150):
    cost_history_dict = {
        "prev_vector": None,
        "iters": 0,
        "cost_history": [],
        "y_min": None,
    }

    # Evaluate the problem using a QPU via `qiskit-ibm-runtime`
    with Session(backend=backend) as session:
        estimator = EstimatorV2(mode=session)
        estimator.skip_transpilation = True
        estimator.options.environment.job_tags = ["TUT_HSVQE"]
        x_opt = spsa(
            cost_func,
            x0=x0,
            args=(isa_ansatz, isa_observable, estimator, cost_history_dict),
            maxiter=maxiter,
        )

        y_min = cost_func(
            x_opt, isa_ansatz, isa_observable, estimator, cost_history_dict
        )

    return y_min, cost_history_dict

np.random.seed(42)
num_params = ansatz.num_parameters
params = 2 * np.pi * np.random.random(num_params)

כאן אנו מגדירים את maxiter = 50. שימו לב שמאחר שכל איטרציה מצריכה שתי קריאות לפונקציה לחישוב הגרדיאנט, סך מספר קריאות הפונקציה יהיה $2 \times \text{maxiter}$. ניתן להגדיל את maxiter לכל ערך גבוה יותר לאומדן אנרגיה טוב יותר.

maxiter = 50
spsa_min, spsa_history = solve(
    params, isa_ansatz, isa_observable, maxiter=maxiter
)

Fx Iters. done: 101 [Current cost: -3.03843]

שלב 4: עיבוד לאחר והחזרת התוצאה בפורמט קלאסי רצוי

• קלט: אומדני אנרגיית מצב-היסוד במהלך האופטימיזציה

• פלט: אנרגיית מצב-יסוד מוערכת

print(f"Estimated ground state energy: {spsa_min}")

Estimated ground state energy: [-3.03842968]

results = {
    "spsa": spsa_history,
}

visualize_results(spsa_history)

דוגמת חומרה בקנה מידה גדול

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

אתגר

עכשיו שיש לכם מימוש VQE עובד עבור שרשרת היזנברג, נסו את הדברים הבאים:

1. ניסוי עם עומק ה-ansatz: שנו את פרמטר reps ב-efficient_su2 (לדוגמה, נסו reps=1 ו-reps=3). כיצד עומק ה-ansatz משפיע על אנרגיית מצב-היסוד המוערכת ומהירות ההתכנסות? באיזו נקודה אתם מבחינים בתשואות פוחתות או חוסר יציבות?

2. כוונון hyperparameters של SPSA: כוונו את פרמטרי לוח הזמנים של קצב הלמידה (a, c, alpha, gamma, A) והתבוננו כיצד הם משפיעים על ההתכנסות. האם תוכלו למצוא הגדרות שמתכנסות מהר יותר מהברירות המחדל שבשימוש כאן?

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

מקורות

[1] Spall, J. C. (2002). Implementation of the simultaneous perturbation algorithm for stochastic optimization. IEEE Transactions on Aerospace and Electronic Systems, 34(3), 817-823.

[2] Sahin, M. Emre, et al. (2025). Qiskit Machine Learning: an open-source library for quantum machine learning tasks at scale on quantum hardware and classical simulators. arXiv:2505.17756.

צעדים הבאים

המלצות

אם מצאתם עבודה זו מעניינת, ייתכן שתתעניינו בחומר הבא:

• נסו Sample-based Quantum Diagonalization (SQD): כפי שהודגם במדריך זה, VQE מתמודד עם אתגרים בקנה מידה גדול בשל מישורים בורים ועומס מדידות גבוה. IBM פיתחה Sample-based Quantum Diagonalization (SQD) כחלופה ניתנת להרחבה יותר. בשונה מ-VQE, SQD נמנעת כלל מאופטימיזציה וריאציונית; במקום זאת, מחשב קוונטי מייצר דגימות ומחשב קלאסי מקרין את ההמילטוניאן על תת-מרחב שנפרש על ידי הדגימות הללו ומבצע עליו דיאגונליזציה. זה מספק חסם עליון לאנרגיית מצב-היסוד עם הרבה פחות מדידות וללא רגישות למישורים בורים. עקבו אחר מדריך SQD כדי לראות גישה זו בפעולה.

• חקרו את קורס אלגוריתמי דיאגונליזציה קוונטית: העמיקו את ההבנה שלכם הן ב-VQE והן ב-SQD, כולל פשרותיהם, בקורס אלגוריתמי דיאגונליזציה קוונטית ב-IBM Quantum Learning.