Instance a rozšíření

Tato kapitola pokryje několik kvantových variačních algoritmů, včetně

• Variační kvantový řešič vlastních čísel (VQE)
• Subspace Search VQE (SSVQE)
• Variational Quantum Deflation (VQD)
• Quantum Sampling Regression (QSR)

Pomocí těchto algoritmů se seznámíme s několika návrhovými myšlenkami, které lze začlenit do vlastních variačních algoritmů, jako jsou váhy, penalizace, převzorkování a podvzorkování. Doporučujeme ti, abys s těmito koncepty experimentoval a podělil se o svá zjištění s komunitou.

Framework Qiskit patterns se vztahuje na všechny tyto algoritmy – jednotlivé kroky však explicitně vyjmenujeme pouze v prvním příkladu.

Variační kvantový řešič vlastních čísel (VQE)

VQE je jedním z nejpoužívanějších variačních kvantových algoritmů a tvoří šablonu, na které mohou stavět další algoritmy.

Krok 1: Mapování klasických vstupů na kvantový problém

Teoretické uspořádání

Uspořádání VQE je jednoduché:

• Příprava referenčních operátorů $U_R$
  • Začínáme ze stavu $|0\rangle$ a přecházíme do referenčního stavu $|\rho\rangle$
• Aplikace variační formy $U_V(\vec\theta_{i,j})$ k vytvoření ansatzu $U_A(\vec\theta_{i,j})$
  • Přecházíme ze stavu $|\rho\rangle$ do $U_V(\vec\theta_{i,j})|\rho\rangle = |\psi(\vec\theta_{i,j})\rangle$
• Bootstrap při $i=0$, pokud máme podobný problém (obvykle zjištěný pomocí klasické simulace nebo vzorkování)
  • Každý optimalizátor bude bootstrapován odlišně, což vede k počáteční sadě vektorů parametrů $\Theta_0 := \\{ {\vec\theta_{0,j} | j \in \mathcal{J}_\text{opt}^0} \\}$ (například z počátečního bodu $\vec\theta_0$).
• Vyhodnocení nákladové funkce $C(\vec\theta_{i,j}) := \langle \psi(\vec{\theta}) | \hat{H} | \psi(\vec{\theta})\rangle$ pro všechny připravené stavy na kvantovém počítači.
• Použití klasického optimalizátoru k výběru další sady parametrů $\Theta_{i+1}$.
• Opakování procesu, dokud není dosaženo konvergence.

Jedná se o jednoduchou klasickou optimalizační smyčku, ve které vyhodnocujeme nákladovou funkci. Některé optimalizátory mohou vyžadovat více vyhodnocení k výpočtu gradientu, určení další iterace nebo posouzení konvergence.

Zde je příklad pro následující pozorovatelnou:

$$\hat{O}_1 = 2 II - 2 XX + 3 YY - 3 ZZ,$$

Implementace

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

from qiskit import QuantumCircuit
from qiskit.quantum_info import SparsePauliOp
from qiskit.circuit.library import TwoLocal
import numpy as np

theta_list = (2 * np.pi * np.random.rand(1, 8)).tolist()
observable = SparsePauliOp.from_list([("II", 2), ("XX", -2), ("YY", 3), ("ZZ", -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)

ansatz.decompose().draw("mpl")

def cost_func_vqe(parameters, 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
    """

    estimator_job = estimator.run([(ansatz, hamiltonian, [parameters])])
    estimator_result = estimator_job.result()[0]

    cost = estimator_result.data.evs[0]
    return cost

from qiskit.primitives import StatevectorEstimator

estimator = StatevectorEstimator()

Tuto nákladovou funkci můžeme použít k výpočtu optimálních parametrů

# SciPy minimizer routine
from scipy.optimize import minimize

x0 = np.ones(8)

result = minimize(
    cost_func_vqe, x0, args=(ansatz, observable, estimator), method="COBYLA"
)

result

message: Optimization terminated successfully.
 success: True
  status: 1
     fun: -5.999999982445723
       x: [ 1.741e+00  9.606e-01  1.571e+00  2.115e-05  1.899e+00
            1.243e+00  6.063e-01  6.063e-01]
    nfev: 136
   maxcv: 0.0

Krok 2: Optimalizace problému pro kvantové spouštění

Vybereme nejméně vytížený backend a naimportujeme potřebné komponenty z qiskit_ibm_runtime.

from qiskit_ibm_runtime import SamplerV2 as Sampler
from qiskit_ibm_runtime import EstimatorV2 as Estimator
from qiskit_ibm_runtime import Session, EstimatorOptions
from qiskit_ibm_runtime import QiskitRuntimeService

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

<IBMBackend('ibm_brisbane')>

Circuit transpilujeme pomocí přednastaveného pass manageru s úrovní optimalizace 3 a na pozorovatelnou veličinu použijeme odpovídající rozložení.

from qiskit.transpiler.preset_passmanagers import generate_preset_pass_manager

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

Krok 3: Spuštění pomocí primitiv Qiskit Runtime

Nyní jsme připraveni spustit náš výpočet na hardwaru IBM Quantum®. Protože minimalizace nákladové funkce je vysoce iterativní, spustíme Runtime session. Díky tomu budeme muset čekat ve frontě pouze jednou. Jakmile úloha začne běžet, každá iterace s aktualizovanými parametry proběhne okamžitě.

x0 = np.ones(8)

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

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

    result = minimize(
        cost_func_vqe,
        x0,
        args=(isa_ansatz, isa_observable, estimator),
        method="COBYLA",
        options={"maxiter": 200, "disp": True},
    )
session.close()
print(result)

Krok 4: Zpracování výsledků do klasického formátu

Vidíme, že minimalizační rutina úspěšně skončila, což znamená, že jsme dosáhli výchozí tolerance klasického optimalizátoru COBYLA. Pokud potřebujeme přesnější výsledek, můžeme specifikovat menší toleranci. To může být skutečně nutné, protože výsledek se lišil o několik procent oproti výsledku získanému výše simulátorem.

Získaná hodnota x je aktuální nejlepší odhad parametrů, které minimalizují nákladovou funkci. Pokud iterujeme pro získání vyšší přesnosti, měly by se tyto hodnoty použít místo původně použitého x0 (vektoru jedniček).

Nakonec poznamenejme, že funkce byla v procesu optimalizace vyhodnocena 96krát. To se může lišit od počtu optimalizačních kroků, protože některé optimalizátory vyžadují více vyhodnocení funkce v jednom kroku, například při odhadu gradientu.

VQE hledající podprostory (SSVQE)

SSVQE je varianta VQE, která umožňuje získat prvních $k$ vlastních hodnot pozorovatelné veličiny $\hat{H}$ s vlastními hodnotami $\{\lambda_0, \lambda_1,...,\lambda_{N-1}\}$, kde $N\geq k$. Bez újmy na obecnosti předpokládáme, že $\lambda_0<\lambda_1<...<\lambda_{N-1}$. SSVQE zavádí novou myšlenku tím, že přidává váhy, které pomáhají upřednostnit optimalizaci členu s největší váhou.

Abychom tento algoritmus implementovali, potřebujeme $k$ vzájemně ortogonálních referenčních stavů $\{ |\rho_j\rangle \}_{j=0}^{k-1}$, což znamená $\langle \rho_j | \rho_l \rangle = \delta_{jl}$ pro $j,l<k$. Tyto stavy lze sestrojit pomocí Pauliho operátorů. Cenová funkce tohoto algoritmu pak je:

$$\begin{aligned} C(\vec{\theta}) & := \sum_{j=0}^{k-1} w_j \langle \rho_j | U_{V}^{\dagger}(\vec{\theta})\hat{H} U_{V}(\vec{\theta})|\rho_j \rangle \\[1mm] & := \sum_{j=0}^{k-1} w_j \langle \psi_{j}(\vec{\theta}) | \hat{H} | \psi_{j}(\vec{\theta}) \rangle \\[1mm] \end{aligned}$$

kde $w_j$ je libovolné kladné číslo takové, že pokud $j<l<k$, pak $w_j>w_l$, a $U_V(\vec{\theta})$ je uživatelem definovaná variační forma.

Algoritmus SSVQE se opírá o skutečnost, že vlastní stavy odpovídající různým vlastním hodnotám jsou vzájemně ortogonální. Konkrétně lze skalární součin $U_V(\vec{\theta})|\rho_j\rangle$ a $U_V(\vec{\theta})|\rho_l\rangle$ vyjádřit jako:

$$\begin{aligned} \langle \rho_j | U_{V}^{\dagger}(\vec{\theta})U_{V}(\vec{\theta})|\rho_l \rangle & = \langle \rho_j | I |\rho_l \rangle \\[1mm] & = \langle \rho_j | \rho_l \rangle \\[1mm] & = \delta_{jl} \end{aligned}$$

První rovnost platí, protože $U_{V}(\vec{\theta})$ je kvantový operátor, a je tedy unitární. Poslední rovnost platí díky ortogonalitě referenčních stavů $|\rho_j\rangle$. Skutečnost, že ortogonalita je zachována při unitárních transformacích, hluboce souvisí s principem zachování informace, jak je vyjádřen v kvantové informační vědě. Z tohoto pohledu neunitární transformace představují procesy, při nichž je informace buď ztracena, nebo vnesena.

Váhy $w_j$ pomáhají zajistit, aby všechny stavy byly vlastními stavy. Pokud se váhy dostatečně liší, bude během optimalizace upřednostněn člen s největší vahou (tj. $w_0$) před ostatními. V důsledku toho se výsledný stav $U_{V}(\vec{\theta})|\rho_0 \rangle$ stane vlastním stavem odpovídajícím $\lambda_0$. Protože $\{ U_{V}(\vec{\theta})|\rho_j\rangle \}_{j=0}^{k-1}$ jsou vzájemně ortogonální, zbývající stavy k němu budou ortogonální, a tedy budou obsaženy v podprostoru odpovídajícím vlastním hodnotám $\{\lambda_1,...,\lambda_{N-1}\}$.

Aplikací stejného argumentu na zbývající členy bude další prioritou člen s vahou $w_1$, takže $U_{V}(\vec{\theta})|\rho_1 \rangle$ bude vlastním stavem odpovídajícím $\lambda_1$ a ostatní členy budou obsaženy ve vlastním prostoru $\{\lambda_2,...,\lambda_{N-1}\}$.

Induktivním uvažováním odvodíme, že $U_{V}(\vec{\theta})|\rho_j \rangle$ bude přibližným vlastním stavem $\lambda_j$ pro $0\leq j < k.$

Teoretické uspořádání

SSVQE lze shrnout následovně:

• Připravte několik referenčních stavů aplikací unitárního operátoru U_R na k různých stavů výpočetní báze
  • Tento algoritmus vyžaduje použití $k$ vzájemně ortogonálních referenčních stavů $\{ |\rho_j\rangle \}_{j=0}^{k-1}$, takže $\langle \rho_j | \rho_l \rangle = \delta_{jl}$ pro $j,l<k$.
• Aplikujte variační formu $U_V(\vec\theta_{i,j})$ na každý referenční stav, což vede k následujícímu ansatzu $U_A(\vec\theta_{i,j})$.
• Použijte bootstrap při $i=0$, pokud je k dispozici podobný problém (obvykle nalezený prostřednictvím klasické simulace nebo vzorkování).
• Vyhodnoťte cenovou funkci $C(\vec\theta_{i,j}) := \sum_{j=0}^{k-1} w_j \langle \psi_{j}(\vec{\theta}) | \hat{H} | \psi_{j}(\vec{\theta}) \rangle$ pro všechny připravené stavy na kvantovém počítači.
  • To lze rozdělit na výpočet střední hodnoty pro pozorovatelnou veličinu $\langle \psi_{j}(\vec{\theta}) | \hat{H} | \psi_{j}(\vec{\theta}) \rangle$ a vynásobení tohoto výsledku hodnotou $w_j$.
  • Poté cenová funkce vrátí součet všech vážených středních hodnot.
• Použijte klasický optimalizátor k určení další sady parametrů $\Theta_{i+1}$.
• Opakujte výše uvedené kroky, dokud nebude dosaženo konvergence.

V testu budete rekonstruovat cenovou funkci SSVQE, ale máme následující úryvek, který vám pomůže s řešením:

import numpy as np

def cost_func_ssvqe(
    params, initialized_anastz_list, weights, ansatz, hamiltonian, estimator
):
    # """Return estimate of energy from estimator

    # Parameters:
    #     params (ndarray): Array of ansatz parameters
    #     initialized_anastz_list (list QuantumCircuit): Array of initialised ansatz with reference
    #     weights (list): List of weights
    #     ansatz (QuantumCircuit): Parameterized ansatz circuit
    #     hamiltonian (SparsePauliOp): Operator representation of Hamiltonian
    #     estimator (Estimator): Estimator primitive instance

    # Returns:
    #     float: Weighted energy estimate
    # """

    energies = []

    # Define SSVQE

    weighted_energy_sum = np.dot(energies, weights)
    return weighted_energy_sum

Variační kvantová deflace (VQD)

VQD je iterativní metoda, která rozšiřuje VQE tak, aby místo pouze nejnižší vlastní hodnoty získala prvních $k$ vlastních hodnot pozorovatelné $\hat{H}$ s vlastními hodnotami $\{\lambda_0, \lambda_1,...,\lambda_{N-1}\}$, kde $N\geq k$. Ve zbytku této sekce budeme bez ztráty obecnosti předpokládat, že $\lambda_0\leq\lambda_1\leq...\leq\lambda_{N-1}$. VQD zavádí pojem penalizační ceny, která vede optimalizační proces.

VQD zavádí penalizační člen, označený jako $\beta$, který vyvažuje příspěvek každého překryvového členu k ceně. Tento penalizační člen slouží k penalizaci optimalizačního procesu, pokud není dosaženo ortogonality. Tuto podmínku zavádíme proto, že vlastní stavy pozorovatelné, nebo hermitovského operátoru, odpovídající různým vlastním hodnotám jsou vždy vzájemně ortogonální, nebo je lze takovými učinit v případě degenerace či opakovaných vlastních hodnot. Vynucením ortogonality s vlastním stavem odpovídajícím $\lambda_0$ tedy efektivně optimalizujeme nad podprostorem, který odpovídá zbývajícím vlastním hodnotám $\{\lambda_1, \lambda_2,..., \lambda_{N-1}\}$. Zde je $\lambda_1$ nejnižší vlastní hodnotou ze zbývajících vlastních hodnot, a proto lze optimální řešení nového problému získat pomocí variačního teorému.

Obecná myšlenka VQD spočívá v tom, že se běžně použije VQE k získání nejnižší vlastní hodnoty $\lambda_0 := C_0(\vec\theta^0) \equiv C_\text{VQE}(\vec\theta^0)$ spolu s odpovídajícím (přibližným) vlastním stavem $|\psi(\vec{\theta^0})\rangle$ pro nějaký optimální vektor parametrů $\vec{\theta^0}$. Poté, abychom získali další vlastní hodnotu $\lambda_1 > \lambda_0$, místo minimalizace cenové funkce $C_0(\vec{\theta}) := \langle \psi(\vec{\theta}) | \hat{H} | \psi(\vec{\theta})\rangle$ optimalizujeme:

$$C_1(\vec{\theta}) := C_0(\vec{\theta})+ \beta_0 |\langle \psi(\vec{\theta})| \psi(\vec{\theta^0})\rangle |^2$$

Kladná hodnota $\beta_0$ by ideálně měla být větší než $\lambda_1-\lambda_0$.

Tím se zavádí nová cenová funkce, kterou lze chápat jako omezenou úlohu, v níž minimalizujeme $C_\text{VQE}(\vec{\theta}) = \langle \psi(\vec{\theta}) | \hat{H} | \psi(\vec{\theta})\rangle$ s podmínkou, že stav musí být ortogonální k dříve získanému $|\psi(\vec{\theta^0})\rangle$, přičemž $\beta_0$ funguje jako penalizační člen, pokud podmínka není splněna.

Alternativně lze tuto novou úlohu interpretovat jako spuštění VQE nad novou pozorovatelnou:

$$\hat{H_1} := \hat{H} + \beta_0 |\psi(\vec{\theta^0})\rangle \langle \psi(\vec{\theta^0})| \quad \Rightarrow \quad C_1(\vec{\theta}) = \langle \psi(\vec{\theta}) | \hat{H_1} | \psi(\vec{\theta})\rangle,$$

Za předpokladu, že řešením nové úlohy je $|\psi(\vec{\theta^1})\rangle$, střední hodnota $\hat{H}$ (nikoli $\hat{H_1}$) by měla být $\langle \psi(\vec{\theta^1}) | \hat{H} | \psi(\vec{\theta^1})\rangle = \lambda_1$. Pro získání třetí vlastní hodnoty $\lambda_2$ je cenovou funkcí, kterou je třeba optimalizovat:

$$C_2(\vec{\theta}) := C_1(\vec{\theta}) + \beta_1 |\langle \psi(\vec{\theta})| \psi(\vec{\theta^1})\rangle |^2$$

kde $\beta_1$ je kladná konstanta dostatečně velká na to, aby vynutila ortogonalitu řešícího stavu vůči $|\psi(\vec{\theta^0})\rangle$ i $|\psi(\vec{\theta^1})\rangle$. To penalizuje stavy v prohledávaném prostoru, které tento požadavek nesplňují, čímž se prohledávaný prostor efektivně zužuje. Optimálním řešením nové úlohy by tedy měl být vlastní stav odpovídající $\lambda_2$.

Stejně jako v předchozím případě lze i tuto novou úlohu interpretovat jako VQE s pozorovatelnou:

$$\hat{H_2} := \hat{H_1} + \beta_1 |\psi(\vec{\theta^1})\rangle \langle \psi(\vec{\theta^1})| \quad \Rightarrow \quad C_2(\vec{\theta}) = \langle \psi(\vec{\theta}) | \hat{H_2} | \psi(\vec{\theta})\rangle.$$

Pokud je řešením této nové úlohy $|\psi(\vec{\theta^2})\rangle$, střední hodnota $\hat{H}$ (nikoli $\hat{H_2}$) by měla být $\langle \psi(\vec{\theta^2}) | \hat{H} | \psi(\vec{\theta^2})\rangle = \lambda_2$. Analogicky pro získání $k$-té vlastní hodnoty $\lambda_{k-1}$ byste minimalizovali cenovou funkci:

$$C_{k-1}(\vec{\theta}) := C_{k-2}(\vec{\theta}) + \beta_{k-2} |\langle \psi(\vec{\theta})| \psi(\vec{\theta^{k-2}})\rangle |^2,$$

Připomeňme, že jsme definovali $\vec{\theta^j}$ tak, že $\langle \psi(\vec{\theta^j}) | \hat{H} | \psi(\vec{\theta^j})\rangle = \lambda_j, \forall j<k$. Tato úloha je ekvivalentní minimalizaci $C(\vec{\theta}) = \langle \psi(\vec{\theta}) | \hat{H} | \psi(\vec{\theta})\rangle$ s podmínkou, že stav musí být ortogonální k $|\psi(\vec{\theta^j})\rangle ; \forall j \in {0, \cdots, k-1}$, čímž se prohledávaný prostor omezuje na podprostor odpovídající vlastním hodnotám $\{\lambda_{k-1},\cdots,\lambda_{N-1}\}$.

Tato úloha je ekvivalentní VQE s pozorovatelnou:

$$\hat{H}_{k-1} := \hat{H}_{k-2} + \beta_{k-2} |\psi(\vec{\theta^{k-2}})\rangle \langle \psi(\vec{\theta^{k-2}})| \quad \Rightarrow \quad C_{k-1}(\vec{\theta}) = \langle \psi(\vec{\theta}) | \hat{H}_{k-1} | \psi(\vec{\theta})\rangle,$$

Jak je z postupu patrné, pro získání $k$-té vlastní hodnoty potřebujete (přibližné) vlastní stavy předchozích $k-1$ vlastních hodnot, takže VQE byste museli spustit celkem $k$-krát. Cenová funkce VQD tedy vypadá následovně:

$$C_k(\vec{\theta}) = \langle \psi(\vec{\theta}) | \hat{H} | \psi(\vec{\theta})\rangle + \sum_{j=0}^{k-1}\beta_j |\langle \psi(\vec{\theta})| \psi(\vec{\theta^j})\rangle |^2$$

Teoretické uspořádání

Uspořádání VQD lze shrnout takto:

• Připravte referenční operátor $U_R$
• Aplikujte variační formu $U_V(\vec\theta_{i,j})$ na referenční stav, čímž vytvoříte následující ansatze $U_A(\vec\theta_{i,j})$
• Pokud máme podobnou úlohu (typicky nalezenou pomocí klasické simulace nebo vzorkování), proveďte bootstrap při $i=0$.
• Vyhodnoťte cenovou funkci $C_k(\vec{\theta})$, která zahrnuje výpočet $k$ excitovaných stavů a pole hodnot $\beta$ definujících překryvovou penalizaci pro každý překryvový člen.
  • Vypočítejte střední hodnotu pozorovatelné $\langle \psi_{j}(\vec{\theta}) | \hat{H} | \psi_{j}(\vec{\theta}) \rangle$ pro každé $k$
  • Vypočítejte penalizaci $\sum_{j=0}^{k-1}\beta_j |\langle \psi(\vec{\theta})| \psi(\vec{\theta^j})\rangle |^2$.
  • Cenová funkce by pak měla vrátit součet těchto dvou členů
• Pomocí klasického optimalizátoru zvolte další sadu parametrů $\Theta_{i+1}$.
• Tento proces opakujte, dokud nedojde ke konvergenci.

Implementace

Pro tuto implementaci vytvoříme funkci pro penalizaci překryvu. Tato penalizace bude použita v účelové funkci v každé iteraci. Tento proces se zopakuje pro každý excitovaný stav.

from qiskit.circuit.library import TwoLocal

ansatz = TwoLocal(2, rotation_blocks=["ry", "rz"], entanglement_blocks="cz", reps=1)

ansatz.decompose().draw("mpl")

Nejprve si připravíme funkci, která spočítá fidelitu stavu – procentuální překryv mezi dvěma stavy, který použijeme jako penalizaci pro VQD:

import numpy as np

def calculate_overlaps(ansatz, prev_circuits, parameters, sampler):
    def create_fidelity_circuit(circuit_1, circuit_2):
        """
        Constructs the list of fidelity circuits to be evaluated.
        These circuits represent the state overlap between pairs of input circuits,
        and their construction depends on the fidelity method implementations.
        """

        if len(circuit_1.clbits) > 0:
            circuit_1.remove_final_measurements()
        if len(circuit_2.clbits) > 0:
            circuit_2.remove_final_measurements()

        circuit = circuit_1.compose(circuit_2.inverse())
        circuit.measure_all()
        return circuit

    overlaps = []

    for prev_circuit in prev_circuits:
        fidelity_circuit = create_fidelity_circuit(ansatz, prev_circuit)
        sampler_job = sampler.run([(fidelity_circuit, parameters)])
        meas_data = sampler_job.result()[0].data.meas

        counts_0 = meas_data.get_int_counts().get(0, 0)
        shots = meas_data.num_shots
        overlap = counts_0 / shots
        overlaps.append(overlap)

    return np.array(overlaps)

Je čas napsat účelovou funkci pro VQD. Stejně jako dříve, když jsme počítali pouze základní stav, určíme stav s nejnižší energií pomocí primitivy Estimator. Nyní však, jak bylo popsáno výše, přidáme penalizační člen, který zajistí ortogonalitu stavů s vyšší energií. To znamená, že pro každý nový excitovaný stav je přidána penalizace za jakýkoli překryv mezi aktuálním variačním stavem a již nalezenými vlastními stavy s nižší energií.

def cost_func_vqd(
    parameters, ansatz, prev_states, step, betas, estimator, sampler, hamiltonian
):
    estimator_job = estimator.run([(ansatz, hamiltonian, [parameters])])

    total_cost = 0

    if step > 1:
        overlaps = calculate_overlaps(ansatz, prev_states, parameters, sampler)
        total_cost = np.sum(
            [np.real(betas[state] * overlap) for state, overlap in enumerate(overlaps)]
        )

    estimator_result = estimator_job.result()[0]

    value = estimator_result.data.evs[0] + total_cost

    return value

Všimni si zejména, že účelová funkce výše odkazuje na funkci calculate_overlaps, která ve skutečnosti vytváří nový kvantový obvod. Pokud chceme spouštět na reálném hardwaru, musí být i tento nový obvod transpilován, nejlépe optimálním způsobem, aby mohl běžet na námi zvoleném backendu. Všimni si, že transpilace byla zabudována do funkcí calculate_overlaps nebo cost_func_vqd. Klidně si zkus sám upravit kód tak, abys do něj zabudoval tuto dodatečnou (a podmíněnou) transpilaci – ale také to za tebe bude uděláno v další lekci.

V této lekci spustíme algoritmus VQD s využitím Statevector Sampleru a Statevector Estimatoru:

from qiskit.primitives import StatevectorEstimator as Estimator

sampler = Sampler()
estimator = Estimator()

Zavedeme pozorovatelnou, kterou budeme odhadovat. V další lekci k tomu přidáme nějaký fyzikální kontext, například excitovaný stav molekuly. Může být užitečné si tuto pozorovatelnou představit jako Hamiltonián systému, který může mít excitované stavy, i když tato pozorovatelná nebyla zvolena tak, aby odpovídala nějaké konkrétní molekule nebo atomu.

from qiskit.quantum_info import SparsePauliOp

observable = SparsePauliOp.from_list([("II", 2), ("XX", -2), ("YY", 3), ("ZZ", -3)])

Zde nastavíme celkový počet stavů, které chceme spočítat (základní stav a excitované stavy, k), a penalizace (betas) za překryv mezi stavovými vektory, které by měly být ortogonální. Důsledky volby příliš vysokých nebo příliš nízkých hodnot bet budou trochu prozkoumány v další lekci. Prozatím jednoduše použijeme hodnoty uvedené níže. Začneme tím, že jako počáteční parametry použijeme samé nuly. Ve vlastních výpočtech možná budeš chtít použít chytřejší počáteční parametry založené na znalosti systému nebo na předchozích výpočtech.

k = 3
betas = [33, 33, 33]
x0 = np.zeros(8)

Nyní můžeme spustit výpočet:

from scipy.optimize import minimize

prev_states = []
prev_opt_parameters = []
eigenvalues = []

for step in range(1, k + 1):
    if step > 1:
        prev_states.append(ansatz.assign_parameters(prev_opt_parameters))

    result = minimize(
        cost_func_vqd,
        x0,
        args=(ansatz, prev_states, step, betas, estimator, sampler, observable),
        method="COBYLA",
        options={
            "maxiter": 200,
        },
    )
    print(result)

    prev_opt_parameters = result.x
    eigenvalues.append(result.fun)

message: Optimization terminated successfully.
 success: True
  status: 1
     fun: -5.999999979545955
       x: [-5.150e-01 -5.452e-02 -1.571e+00 -2.853e-05  2.671e-01
           -2.672e-01 -8.509e-01 -8.510e-01]
    nfev: 131
   maxcv: 0.0
 message: Optimization terminated successfully.
 success: True
  status: 1
     fun: 4.024550284767612
       x: [-3.745e-01  1.041e+00  8.637e-01  1.202e+00 -8.847e-02
            1.181e-02  7.611e-01 -3.006e-01]
    nfev: 110
   maxcv: 0.0
 message: Optimization terminated successfully.
 success: True
  status: 1
     fun: 5.608925562838559
       x: [-2.670e-01  1.280e+00  1.070e+00 -8.031e-01 -1.524e-01
           -6.956e-02  7.018e-01  1.514e+00]
    nfev: 90
   maxcv: 0.0

Hodnoty, které jsme získali z účelové funkce, jsou přibližně -6,00, 4,02 a 5,61. Důležité na těchto výsledcích je, že hodnoty funkce rostou. Pokud bychom získali první excitovaný stav s nižší energií, než byl náš prvotní neomezený výpočet základního stavu, naznačovalo by to chybu někde v našem kódu.

Hodnoty x jsou parametry, které daly stavový vektor odpovídající každé z těchto cen (energií).

Na závěr podotýkáme, že všechny tři minimalizace konvergovaly v rámci výchozí tolerance klasického optimalizátoru (zde COBYLA). Vyžadovaly 131, 110 a 90 vyhodnocení funkce.

Kvantová vzorkovací regrese (QSR)

Jedním z hlavních problémů VQE je potřeba vícenásobných volání kvantového počítače k získání parametrů pro každý krok, například $k$, $k-1$ a tak dále. To je obzvláště problematické, když je přístup ke kvantovým zařízením ve frontě. Zatímco Session lze použít ke seskupení více iterativních volání, alternativním přístupem je použití vzorkování. Využitím více klasických zdrojů můžeme dokončit celý optimalizační proces v jediném volání. Zde přichází na řadu kvantová vzorkovací regrese. Jelikož je přístup ke kvantovým počítačům stále komoditou s nízkou nabídkou a vysokou poptávkou, považujeme tento kompromis za možný a výhodný pro mnoho současných studií. Tento přístup využívá všechny dostupné klasické schopnosti a přitom zachycuje mnoho vnitřních mechanismů a vlastních vlastností kvantových výpočtů, které se v simulaci neobjevují.

Myšlenka za QSR spočívá v tom, že nákladovou funkci $C(\theta) := \langle \psi(\theta) | \hat{H} | \psi(\theta)\rangle$ lze vyjádřit jako Fourierovu řadu následujícím způsobem:

$$\begin{aligned} C(\vec{\theta}) & := \langle \psi(\theta) | \hat{H} | \psi(\theta)\rangle \\[1mm] & := a_0 + \sum_{k=1}^S[a_k\cos(k\theta)+ b_k\sin(k\theta)] \\[1mm] \end{aligned}$$

V závislosti na periodicitě a šířce pásma původní funkce může být množina $S$ konečná nebo nekonečná. Pro účely této diskuse budeme předpokládat, že je nekonečná. Dalším krokem je několikanásobné vzorkování nákladové funkce $C(\theta)$, abychom získali Fourierovy koeficienty $\{a_0, a_k, b_k\}_{k=1}^S$. Konkrétně, jelikož máme $2S+1$ neznámých, budeme muset nákladovou funkci vzorkovat $2S+1$-krát.

Pokud pak vzorkujeme nákladovou funkci pro $2S+1$ hodnot parametrů $\{\theta_1,...,\theta_{2S+1}\}$, můžeme získat následující soustavu:

$$\begin{pmatrix} 1 & \cos(\theta_1) & \sin(\theta_1) & \cos(2\theta_1) & ... & \sin(S\theta_1) \\ 1 & \cos(\theta_2) & \sin(\theta_2) & \cos(2\theta_2) & \cdots & \sin(S\theta_2)\\ \vdots & \vdots & \vdots & \vdots & \ddots & \vdots\\ 1 & \cos(\theta_{2S+1}) & \sin(\theta_{2S+1}) & \cos(2\theta_{2S+1}) & \cdots & \sin(S\theta_{2S+1}) \end{pmatrix} \begin{pmatrix} a_0 \\ a_1 \\ b_1 \\ a_2 \\ \vdots \\ b_S \end{pmatrix} = \begin{pmatrix} C(\theta_1) \\ C(\theta_2) \\ \vdots \\ C(\theta_{2S+1}) \end{pmatrix},$$

kterou přepíšeme jako

$$Fa=c.$$

V praxi tato soustava obvykle není konzistentní, protože hodnoty nákladové funkce $c$ nejsou přesné. Proto je obvykle dobrým nápadem je normalizovat vynásobením zleva $F^\dagger$, což vede k:

$$F^\dagger Fa = F^\dagger c.$$

Tato nová soustava je vždy konzistentní a jejím řešením je řešení původní úlohy metodou nejmenších čtverců. Pokud máme $k$ parametrů místo jednoho a každý parametr $\theta^i$ má své vlastní $S_i$ pro $i \in {1,...,k}$, pak je celkový počet požadovaných vzorků:

$$T=\prod_{i=1}^k(2S_i+1)\leq \prod_{i=1}^k(2S_{max}+1) = (2S_{max}+1)^n,$$

kde $S_{\max} = \max_i(S_i)$. Navíc nastavení $S_{\max}$ jako laditelného parametru (místo jeho odvozování) otevírá nové možnosti, jako například:

• Nadvzorkování pro zlepšení přesnosti.
• Podvzorkování pro zvýšení výkonu snížením režie doby běhu nebo odstraněním lokálních minim.

Teoretické uspořádání

Uspořádání QSR lze shrnout následovně:

• Připrav referenční operátory $U_R$.
  • Přejdeme ze stavu $|0\rangle$ do referenčního stavu $|\rho\rangle$
• Aplikuj variační formu $U_V(\vec\theta_{i,j})$ k vytvoření ansatzu $U_A(\vec\theta_{i,j})$.
  • Urči šířku pásma spojenou s každým parametrem v ansatzu. Horní mez postačuje.
• Proveď bootstrap při $i=0$, pokud máme podobný problém (typicky nalezený klasickou simulací nebo vzorkováním).
• Vzorkuj nákladovou funkci $C(\vec\theta) := a_0 + \sum_{k=1}^S[a_k\cos(k\theta)+ b_k\sin(k\theta)]$ alespoň $T$-krát.
  • $T=\prod_{i=1}^k(2S_i+1)\leq \prod_{i=1}^k(2S_{max}+1) = (2S_{max}+1)^n$
  • Rozhodni, zda nadvzorkovat/podvzorkovat pro vyvážení rychlosti vs. přesnosti úpravou $T$.
• Spočítej Fourierovy koeficienty ze vzorků (tj. vyřeš normalizovanou lineární soustavu rovnic).
• Vyřeš globální minimum výsledné regresní funkce na klasickém stroji.

Shrnutí

V této lekci ses seznámil/a s vícero dostupnými variačními instancemi:

• Obecné uspořádání
• Zavedení vah a penalizací pro úpravu nákladové funkce
• Prozkoumání podvzorkování vs. nadvzorkování pro kompromis mezi rychlostí a přesností

Tyto myšlenky lze upravit tak, aby vytvořily vlastní variační algoritmus, který odpovídá tvému problému. Povzbuzujeme tě, abys své výsledky sdílel/a s komunitou. Další lekce prozkoumá, jak použít variační algoritmus k řešení aplikace.