Оцінка енергії основного стану ланцюга Гайзенберга за допомогою VQE

Орієнтовне використання: 37 хвилин на процесорі Heron (ПРИМІТКА: це лише орієнтовна оцінка. Твій час виконання може відрізнятися.)

Результати навчання

Після завершення цього посібника ти зможеш розуміти наступне:

• Як моделювати спіновий ланцюг Гайзенберга як квантовий гамільтоніан за допомогою Qiskit
• Як використовувати оптимізатор SPSA для оцінки енергії основного стану квантової системи
• Як виконувати варіаційні робочі процеси на квантовому обладнанні IBM® за допомогою примітивів та сесій Qiskit Runtime

Передумови

Рекомендується ознайомитися з такими темами:

• Основи квантової інформації
• Вступ до шаблонів 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]$. Зверни увагу, що в реалізації нижче набір «найближчих сусідів» визначається нативним з'єднанням апаратного бекенду між першими $N$ кубітами, яке може не утворювати суворий лінійний ланцюг залежно від топології пристрою.

Розуміння енергії основного стану цього гамільтоніана має фундаментальне значення у фізиці. Основний стан містить інформацію про квантові фазові переходи, структуру заплутаності та магнітний порядок. Класично обчислення точної енергії основного стану стає нездійсненним у міру зростання кількості спінів, оскільки розмірність простору Гільберта експоненційно збільшується як $2^N$ для $N$ спінів. Це робить задачу природним кандидатом для квантового моделювання.

Варіаційний квантовий власний вирішувач (VQE) — це гібридний квантово-класичний алгоритм, призначений для оцінки енергії основного стану гамільтоніана. Він працює, готуючи параметризований квантовий стан $|\psi(\theta)\rangle$ (що називається анзацем) на квантовому комп'ютері та вимірюючи очікуване значення $\langle \psi(\theta) | H | \psi(\theta) \rangle$. Класичний оптимізатор потім ітеративно коригує параметри $\theta$, щоб мінімізувати цю енергію, використовуючи варіаційний принцип, який гарантує, що виміряна енергія завжди є верхньою межею справжньої енергії основного стану.

У цьому посібнику ми використовуємо анзац efficient_su2 з бібліотеки схем Qiskit, який будує шари однокубітних обертань та заплутувальних вентилів. Оптимізація виконується за допомогою алгоритму Simultaneous Perturbation Stochastic Approximation (SPSA), який добре підходить для шумного квантового обладнання, оскільки оцінює градієнти, використовуючи лише два обчислення функції на ітерацію, незалежно від кількості параметрів.

Вимоги

Перед початком роботи з цим посібником переконайся, що у тебе встановлено наступне:

• Qiskit SDK v2.0 або новішої версії, з підтримкою візуалізації
• Qiskit Runtime v0.44 або новішої версії (pip install qiskit-ibm-runtime)

Налаштування

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

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: Перетворення класичних вхідних даних у квантову задачу

• Вхідні дані: кількість спінів
• Вихідні дані: анзац та гамільтоніан, що моделюють ланцюг Гайзенберга

Побудуй анзац та гамільтоніан, які моделюють ланцюг Гайзенберга з 10 спінами. На цьому кроці ми побудуємо гамільтоніан Гайзенберга з 10 спінами на карті зв'язків найменш зайнятого бекенду та підготуємо анзац 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: Оптимізація задачі для виконання на квантовому обладнанні

• Вхідні дані: абстрактна схема, спостережувана величина
• Вихідні дані: цільова схема та спостережувана величина, оптимізовані для обраного 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: Виконання за допомогою примітивів Qiskit

• Вхідні дані: цільова схема та спостережувана величина
• Вихідні дані: результати оптимізації

Мінімізуй оцінену енергію основного стану системи шляхом оптимізації параметрів схеми. Використовуй примітив Estimator з Qiskit Runtime для обчислення функції вартості під час оптимізації.

Оскільки ми оптимізували схему для бекенду на кроці 2, ми можемо уникнути транспіляції на сервері Runtime, встановивши skip_transpilation=True та передавши оптимізовану схему. Для цієї демонстрації ми будемо виконувати обчислення на QPU за допомогою примітивів qiskit-ibm-runtime. Щоб виконати обчислення з примітивами на основі вектора стану qiskit, заміни блок коду, що використовує примітиви Qiskit Runtime, на закоментований блок.

У цьому посібнику ми використовуємо метод Simultaneous Perturbation Stochastic Approximation (SPSA) — оптимізатор на основі градієнта. Нижче ми коротко його представимо та надамо код для реалізації 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$-му кроці ітерації. Для обчислення градієнта створюється випадковий вектор $\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 особливо відомий своєю стійкістю до шуму та вимагає лише двох апаратних викликів для кожної ітерації. Тому він є одним з найбільш переважних оптимізаторів для реалізації варіаційних алгоритмів.

У цьому посібнику гіперпараметри для $i$-ї ітерації, $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]. Відповідне налаштування гіперпараметрів необхідне для отримання хорошої продуктивності від 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)

Приклад на великому масштабі апаратного забезпечення

Приклад у великому масштабі апаратного забезпечення не включено до цього посібника. Зі збільшенням кількості кубітів VQE стикається зі значними труднощами через явище barren plateau: градієнт функції вартості зникає експоненційно зі збільшенням розміру системи, що робить оптимізацію практично нездійсненною для великих схем. У поєднанні з апаратними шумами це означає, що масштабування VQE до більших спінових ланцюгів не дає надійно відтворюваних результатів. Щодо підходів, які долають ці обмеження, дивися розділ «Наступні кроки» нижче.

Завдання

Тепер, коли ти маєш робочу реалізацію VQE для ланцюга Гайзенберга, спробуй наступне:

1. Експериментуй із глибиною анзацу: Змінюй параметр reps в efficient_su2 (наприклад, спробуй reps=1 та reps=3). Як глибина анзацу впливає на оцінену енергію основного стану та швидкість збіжності? При якому значенні ти спостерігаєш спадну віддачу або нестабільність?
2. Налаштовуй гіперпараметри SPSA: Коригуй параметри графіка швидкості навчання (a, c, alpha, gamma, A) та спостерігай, як вони впливають на збіжність. Чи зможеш ти знайти конфігурацію, яка збігається швидше, ніж значення за замовчуванням, використані тут?
3. Порівнюй топології з'єднань: Замість використання нативної карти зв'язків бекенду спробуй побудувати простий лінійний ланцюг найближчих сусідів і порівняй результати. Як зв'язність фізичного апаратного забезпечення впливає на глибину транспільованої схеми та кінцеву оцінку енергії?

Посилання

[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 стикається з труднощами при масштабуванні через barren plateau та великі витрати на вимірювання. IBM розробила Sample-based Quantum Diagonalization (SQD) як більш масштабовану альтернативу. На відміну від VQE, SQD повністю уникає варіаційної оптимізації; натомість квантовий комп'ютер генерує зразки, а класичний комп'ютер проектує гамільтоніан на підпростір, натягнутий на ці зразки, та діагоналізує його. Це забезпечує верхню межу енергії основного стану зі значно меншою кількістю вимірювань і без сприйнятливості до barren plateau. Дотримуйся посібника SQD, щоб побачити цей підхід у дії.
• Вивчи курс алгоритмів квантової діагоналізації: Поглиб своє розуміння VQE та SQD, включно з їхніми компромісами, у курсі Алгоритми квантової діагоналізації на IBM Quantum Learning.