Algoritmul SqDRIFT pentru estimarea stării fundamentale

Python

Descarcă notebook-ul Python

Descarcă notebook-ul Python pentru a rula singur acest tutorial. Notebook-ul descarcă datele de intrare FCIDump dacă nu sunt deja prezente.

Estimare de utilizare: 180 de secunde pe un procesor Heron r3 (NOTĂ: Aceasta este doar o estimare. Timpul tău de execuție poate varia.)

Rezultate ale învățării

• Învață cum să creezi circuite cu adâncime mai mică în comparație cu Trotterizarea

• Parcurge un flux de lucru complet pentru estimarea stării fundamentale folosind qDRIFT și SQD

• Învață cum să folosești qiskit-fermions împreună cu alte addon-uri Qiskit pentru a implementa un astfel de flux de lucru

Cerințe preliminare

• Citește prezentarea generală a Diagonalizării cuantice bazate pe eșantioane (SQD)

• Citește lecția Diagonalizare cuantică Krylov bazată pe eșantioane (SKQD)

Context

SqDRIFT este o variantă a SKQD care elimină nevoia de a alege un ansatz din care să eșantionezi șiruri de biți, folosind în schimb un ansamblu de circuite de evoluție temporală construite direct din Hamiltonianul țintă. Acest lucru se realizează prin subeșantionarea unor operatori de evoluție temporală mai mici din Hamiltonian, pe baza coeficienților săi, cunoscută sub numele de metoda de Trotterizare qDRIFT.

Ambele implementări folosesc Qiskit Fermions pentru a reprezenta Hamiltonianul molecular și a construi circuite pentru algoritmul qDRIFT. Implementarea Python folosește pass-uri fermionice de Transpiler, în timp ce implementarea C++ eșantionează explicit grupuri și le mapează la operatori de qubiți.

Fie Hamiltonianul de forma:

$$H = \sum_{i=1}^{N} c_i h_i$$

unde, fără a pierde din generalitate, cerem ca $c_i > 0$ și ca cea mai mare valoare proprie a lui $h_i$ să fie egală, în valoare absolută, cu $1$. Orice prefactor cu semn sau complex este absorbit în $h_i$, deci coeficienții $c_i$ sunt ponderi strict pozitive, în timp ce $h_i$ poartă direcția fiecărui termen. Aici $N$ este numărul de termeni (sau, după grupare, numărul de grupuri) din hamiltonian; este o proprietate a hamiltonianului și este distinct de numărul de operatori eșantionați într-un singur circuit, notat $n$ mai jos.

Algoritmul qDRIFT realizează apoi, pentru timpul țintă $t$, un anumit operator $V_k$, unde $k$ merge de la $1 \cdots K$ și semnifică al $k_{th}$ circuit SqDRIFT, definit ca:

$$V_k = \prod_{j=1}^{n} e^{-i h_{k_j} \lambda t / n }$$

Aici $n$ este numărul de operatori eșantionați per circuit, iar $K$ este numărul de circuite din ansamblu. Produsul se face peste cele $n$ trageri, nu peste toți cei $N$ termeni ai hamiltonianului, iar deoarece termenii sunt trași cu revenire, același $h_i$ poate apărea de mai multe ori într-un singur $V_k$.

Cantitatea:

$$\lambda = \sum_{i=1}^{N} c_i$$

este norma $L_1$ a coeficienților, deci fiecare dintre cei $n$ pași evoluează pe aceeași durată $\lambda t / n$, indiferent de termenul care a fost tras. Uniformitatea unghiului de pas este caracteristica esențială a qDRIFT: un coeficient influențează rezultatul prin cât de des este tras termenul său, nu prin cât de mult este rotit acel termen. Indicii sunt eșantionați din distribuția:

$$P[k_i] = \frac{c_i}{\lambda}$$

deci seria $(k_1, \ldots, k_n)$ este o secvență aleatoare de indici de termeni trași din această distribuție. Deoarece $c_i$ sunt pozitivi și însumează $\lambda$, aceasta este o distribuție de probabilitate normalizată, iar valoarea așteptată a canalului rezultat, peste tragerile aleatoare, aproximează evoluția sub $H$, cu o eroare care scade pe măsură ce $n$ crește. Observă că eroarea de aproximare depinde de $\lambda$, nu de numărul de termeni $N$.

(Articolul SqDRIFT scrie numărul de termeni ca $\mathcal{N}$ și lungimea secvenței ca $N$; noi folosim $N$ și $n$ aici pentru a-i menține clar distincți.)

Acest tutorial arată cum să generezi un ansamblu de astfel de circuite randomizate. După ce am creat aceste circuite, similar modului în care creăm un subspațiu Krylov pentru operatori diferiți, eșantionăm șiruri de biți din mai mulți astfel de operatori cu parametri de timp diferiți. Aceasta asigură o suprapunere mai mare între vectorii stării fundamentale și șirurile de biți eșantionate.

Cerințe

Înainte de a începe acest tutorial, asigură-te că ai instalat

• Un mediu virtual Python (>=3.10)
• pip>=25.1
• qiskit ~= 2.5
• qiskit-fermions==0.1.0 (Reține că numele este la plural)
• numpy
• pyscf
• qiskit-aer
• qiskit-ibm-runtime
• qiskit-addon-sqd

Poți instala toate pachetele necesare cu:

pip install "qiskit~=2.5" "qiskit-fermions==0.1.0" qiskit-aer qiskit-ibm-runtime qiskit-addon-sqd pyscf numpy

Configurare

# Third-party scientific computing
import numpy as np

# PySCF
from pyscf import tools, ao2mo, fci

# Qiskit core
from qiskit import transpile
from qiskit.primitives import BitArray

# Qiskit Aer
from qiskit_aer import AerSimulator

# IBM Quantum Compute Service
from qiskit_ibm_runtime import QiskitRuntimeService, SamplerV2 as Sampler

# Qiskit Fermions
from qiskit_fermions.operators.library import FCIDump
from qiskit_fermions.operators import FermionOperator
from qiskit_fermions.operators.terms.filtering import filter_diagonal_terms
from qiskit_fermions.operators.terms.grouping import (
    group_terms_by_electronic_structure,
)
from qiskit_fermions.operators.terms.ordering import canonical_order
from qiskit_fermions.circuit import FermionicCircuit
from qiskit_fermions.circuit.library import Evolution
from qiskit_fermions.transpiler import FermionicPassManager
from qiskit_fermions.transpiler.presets import generate_preset_jw_pass_manager
from qiskit_fermions.transpiler.passes import QDriftTrotterization
from qiskit_fermions.circuit.library import InitializeModes

# Qiskit addon SQD
from qiskit_addon_sqd.fermion import (
    diagonalize_fermionic_hamiltonian,
    SCIResult,
)

Implementare

Pasul 1: Mapează intrările clasice la o problemă cuantică

Pașii Python de mai jos folosesc simulatorul Aer. Exemplul pe hardware urmează după aceea.

Citirea și pregătirea FCIDump

Pentru acest tutorial, vom încărca Hamiltonianul de structură electronică pentru azot (N2). Există și alte moduri de a crea operatori fermionici. Consultă documentația la qiskit_fermions.operators.library.

Despre acest FCIDump. Fișierul N2_sto_3g descrie o moleculă de azot ($N_2$) în baza minimală STO-3G, la o separare interatomică de 1,09 $\AA$, lungimea experimentală de echilibru a legăturii. Antetul său declară NORB=10, NELEC=14 și MS2=0: 10 orbitali spațiali (deci 20 de orbitali de spin și 20 de qubiți sub transformarea Jordan-Wigner), 14 electroni într-un singlet de spin, deci șapte electroni $\alpha$ și șapte $\beta$. Toți orbitalii primesc eticheta de simetrie 1, adică nu se exploatează nicio simetrie de grup punctual. Fiind un dump STO-3G pe întregul spațiu, niciun orbital nu este înghețat, iar spațiul de corelare este suficient de mic încât o energie de referință FCI exactă poate fi calculată clasic pentru comparație, așa cum se arată în blocul de cod următor.

Un fișier echivalent poate fi regenerat cu PySCF:

from pyscf import gto, scf, tools

mol = gto.M(atom="N 0 0 0; N 0 0 1.09", basis="sto-3g", symmetry=False)
mf = scf.RHF(mol).run()
tools.fcidump.from_scf(mf, "N2_sto_3g")

Deoarece integralele depind de orbitalii SCF convergenți, un fișier regenerat poate diferi de cel livrat în privința fazei sau ordinii orbitalilor; energiile totale nu sunt afectate.

Obținerea fișierului. Găsește FCIDump în acest repozitoriu GitHub. Poți rula codul de mai jos pentru a-l prelua în locația pe care restul tutorialului o așteaptă.

Mai întâi folosim cisolver oferit de pyscf pentru a obține energia de referință. Aceasta este energia reală a stării fundamentale a moleculei cu care lucrăm. Pentru asta, vom declara mai întâi norb și nelec, care reprezintă numărul de orbitali, respectiv numărul de electroni. Apoi declarăm h1e și h2e, care reprezintă integralele monoelectronice, respectiv bielectronice. Toate acestea vor fi folosite mai târziu și pentru SQD.

import os
from urllib.request import urlopen

# The FCIDump is stored with this tutorial in the Qiskit documentation repository.
FCIDUMP_URL = "https://raw.githubusercontent.com/Qiskit/documentation/main/docs/tutorials/assets/sqdrift/fcidump_files/N2_sto_3g"
FCIDUMP_PATH = "fcidump_files/N2_sto_3g"

if not os.path.exists(FCIDUMP_PATH):
    os.makedirs(os.path.dirname(FCIDUMP_PATH), exist_ok=True)
    with urlopen(FCIDUMP_URL) as response:
        contents = response.read()
    with open(FCIDUMP_PATH, "wb") as f:
        f.write(contents)
    print(f"Downloaded FCIDump to {FCIDUMP_PATH}")
else:
    print(f"Using existing FCIDump at {FCIDUMP_PATH}")

name = "fcidump_files/N2_sto_3g"

fcidump = tools.fcidump.read(name)

# Extract metadata from the FCIDump header
norb = fcidump["NORB"]  # number of spatial orbitals
nelec = fcidump["NELEC"]  # total number of electrons
e_nuc = fcidump["ECORE"]  # nuclear repulsion / core energy
ms2 = fcidump["MS2"]  # 2S (spin)

num_elec_a = (nelec + ms2) // 2  # alpha electrons
num_elec_b = (nelec - ms2) // 2  # beta  electrons

# Reconstruct full 4-index ERIs from the FCIDump (stored in 8-fold symmetry)
h1e = fcidump["H1"]  # shape (norb, norb)
h2e = ao2mo.restore(  # shape (norb, norb, norb, norb)
    1, fcidump["H2"], norb
)

cisolver = fci.direct_spin1.FCI()
cisolver.max_cycle = 200
cisolver.conv_tol = 1e-12

e_fci, _ = cisolver.kernel(
    h1e,
    h2e,
    norb,
    (num_elec_a, num_elec_b),
    ecore=e_nuc,  # adds nuclear repulsion to the final energy
)

reference_energy = e_fci

print(f"Reference FCI Energy  = {reference_energy:.10f} Ha")

nuclear_repulsion_energy = fcidump["ECORE"]
print(f"Nuclear Repulsion Energy = {nuclear_repulsion_energy:.10f} Ha")

Reference FCI Energy  = -107.6481842917 Ha
Nuclear Repulsion Energy = 23.7887003074 Ha

Încărcarea hamiltonianului

Cu datele necesare pregătite, citim hamiltonianul din fișierul FCI într-un format compatibil cu qiskit-fermions

fcidump = FCIDump.from_file(name)
hamiltonian = FermionOperator.from_fcidump(fcidump)
num_modes = 2 * fcidump.norb

Fluxuri de lucru fermionice cu qiskit-fermions

Vom mapa mai întâi hamiltonianul într-un model de circuit fermionic folosind qiskit-fermions, care oferă pași de transpiler și porți specifice circuitelor fermionice. Aceștia vor fi folosiți mai târziu, înainte de pașii tradiționali de transpiler ai Qiskit pentru acest flux de lucru.

Gruparea termenilor

Pentru a asigura reproductibilitatea rezultatelor, folosim mai întâi canonical_order pentru a sorta termenii doar pe baza structurii lor. Ordinea operatorilor din lista canon este astfel fixă. Acest lucru asigură reproductibilitatea operatorilor creați deoarece pasul QDriftTrotterization, pe care îl vom folosi mai târziu, eșantionează indici aleatorii pentru a crea operatorii qDRIFT.

În acest pas, exploatăm multiplele simetrii prezente în hamiltonianul structurii electronice, grupând termenii înrudiți cu coeficienți identici. Deși acest lucru schimbă distribuția coeficienților operatorului din care eșantionează protocolul qDRIFT, acest fapt nu afectează garanțiile sale de convergență. Esențial, gruparea termenilor înrudiți prin simetrie duce la o anulare favorabilă a termenilor Pauli și la o adâncime globală mai mică a circuitului atunci când o stare este evoluată în timp sub acțiunea lor.

qiskit-fermions oferă funcția group_terms_by_electronic_structure care realizează această grupare pentru noi.

Reține că group_terms_by_electronic_structure presupune termeni în ordine normală.

Filtrarea termenilor diagonali

Eliminăm termenii diagonali din hamiltonianul folosit pentru a genera circuitele, astfel încât cele $n$ sloturi de eșantionare qDRIFT să fie folosite pentru termenii care mută populația între configurații. Astfel de termeni sunt cel mai bine filtrați din hamiltonian în acest punct, înainte ca poarta Evolution să fie construită în pasul următor.

Termenii în discuție sunt cei diagonali în baza numărului de ocupare, adică produsele operatorilor de număr $a^\dagger_i a_i$. Trei tipuri de termeni se încadrează în această descriere:

• decalajul constant de energie, un produs de zero operatori de număr, a cărui evoluție în timp contribuie doar cu o fază globală;

• operatorii individuali de număr $n_i$, a căror evoluție în timp se reduce la rotații $Z$ pe un singur qubit;

• produsele de ordin superior, precum $n_i n_j$.

Luate separat, niciunul dintre acestea nu mută populația între configurațiile de număr de ocupare; ele acționează doar asupra fazelor configurațiilor deja prezente. Totuși, nu sunt inerte: aceste faze relative alimentează interferența generată de termenii de excitație mai târziu în circuit, astfel încât filtrarea lor schimbă evoluția generată efectiv și poate schimba distribuția de eșantionare. Aceasta este o aproximare deliberată în pasul de generare a circuitului, făcută pentru a concentra eșantionarea pe termenii de excitație, și nu un pas care lasă neatinsă distribuția eșantionată. Spre deosebire de gruparea prin simetrie de mai sus, care lasă intacte garanțiile de convergență ale qDRIFT, acest filtru schimbă operatorul care este evoluat. Prin urmare, circuitele nu mai aproximează evoluția sub hamiltonianul complet, iar limitele de eroare ale qDRIFT se aplică operatorului filtrat, nu celui original. Acest lucru este acceptabil aici deoarece circuitele sunt doar o euristică de eșantionare folosită pentru a propune configurații: niciun termen nu se pierde din estimarea energiei în sine, deoarece filtrul se aplică doar hamiltonianului folosit pentru a construi circuitele, în timp ce diagonalizarea clasică folosește ulterior hamiltonianul complet, inclusiv termenii diagonali. Precizia SQD depinde de acel pas clasic, care rămâne variațional în subspațiul eșantionat, indiferent de modul în care au fost propuse configurațiile.

Funcția filter_diagonal_terms() elimină astfel de termeni dintr-un operator direct (in place). Ea îi identifică pe baza structurii lor în ordine normală — multisetul modurilor de creare care se potrivește cu multisetul modurilor de anihilare — astfel încât este validă doar pe un operator care este deja în ordine normală. Această presupunere nu este verificată la momentul execuției.

# Apply automatic grouping
canon = canonical_order(hamiltonian.normal_ordered().simplify(atol=1e-16))
exit_code = group_terms_by_electronic_structure(
    canon, num_modes, two_body_physicist_order=False
)
filter_diagonal_terms(canon)

print(len(canon.groups))

5060

Acum că am grupat termenii din hamiltonian, vom stabili următorii parametri pentru a genera ansamblul de circuite:

• Numărul de circuite de generat: num_circuits
• Lungimea fiecărui circuit în termeni de grupuri de excitație: num_exc
• Factorul pentru diferiții timpi de evoluție: times

Crearea circuitelor fermionice

Vom crea acum circuite fermionice pentru fiecare dintre pașii de timp. Fiecare circuit va consta dintr-o singură poartă de evoluție, cu timpul de evoluție declarat anterior. Operatorul de evoluție este hamiltonianul. Mai târziu rulăm pași de transpiler pe aceste circuite pentru a crea circuite qDRIFT.

Pregătirea ansatzului

Pregătim starea Hartree-Fock folosind clasa InitializeModes. Pentru azot, procesul constă pur și simplu în aplicarea porților X pe primii num_elec_a qubiți și apoi pe num_elec_b qubiți, ambele fiind egale cu șapte pentru azot. Această stare reprezintă cei șapte electroni $\alpha$ și cei șapte $\beta$ ai azotului.

# SqDRIFT parameters
times = [1.0, 10.0]  # Total evolution times used for the subspace creation
num_exc = 10  # Number of excitation groups per circuit
num_circuits = 200  # Number of circuits to generate

init_circuits = []

hf_gate = InitializeModes.from_hartree_fock(norb, (num_elec_a, num_elec_b))

for time in times:
    evo_gate = Evolution(num_modes, canon, time)
    circ = FermionicCircuit(num_modes)
    circ.append(hf_gate, circ.modes)
    circ.append(evo_gate, circ.modes)
    init_circuits.append(circ)

Pasul 2: Optimizarea problemei pentru execuția pe hardware cuantic

Acum că avem circuitele noastre, vom folosi mai întâi pașii disponibili în qiskit-fermions pentru a efectua optimizări la nivel fermionic, urmate de transpilarea circuitului nostru pentru backend-ul ales. Deoarece acesta este un experiment cu simulator, vom face mai întâi acest lucru pentru AerSimulator.

Calculul ponderii pentru fiecare grup

În acest pas, efectuăm eșantionarea qDRIFT a termenilor stocastic, cu probabilități proporționale cu coeficienții lor din hamiltonian. Pasul de transpiler qDRIFT face acest lucru pentru noi. Putem crea acum circuite mai puțin adânci, care pot fi executate pe hardware mai eficient în ciuda conectivității limitate a qubiților, chiar și atunci când hamiltonianul conține cuplaje pe rază lungă și termeni de ordin mai mare decât pătratic. După gruparea termenilor, se eșantionează operatorii pe baza ponderilor lor. Pentru fiecare operator $h_i$, ponderea $W_{h_i}$ este definită după cum urmează:

$$W_{h_i} = |c_i| / \lambda$$

Deoarece termenii au fost grupați în Pasul 1, fiecare $h_i$ de aici este un grup întreg: $c_i$ este media coeficienților absoluți ai termenilor din grupul $i$, iar fiecare termen din grup este evoluat cu coeficientul său redus la semn.

Optimizări fermionice și native hardware-ului

Funcția generate_preset_jw_pass_manager() returnează un MultiStagePassManager care primește un FermionicCircuit și produce un circuit final optimizat pe care îl putem transpila pentru a rula pe hardware-ul nostru. Înlocuim etapa sa implicită de optimizare cu un FermionicPassManager care conține pasul nostru QDriftTrotterization:

• Pasul QDriftTrotterization folosește intern calculul ponderilor și eșantionarea pentru a genera circuitele pe care le vom folosi pentru eșantionare

• Pasul RelabelModes este un alt pas de optimizare care poate fi folosit pentru a permuta modurile fermionice pentru a optimiza conectivitatea între qubiți și a reduce adâncimea porților; citește mai multe în referința API

Etapele rămase ale MultiStagePassManager rulează automat și se ocupă de maparea completă fermion-qubit:

• F2QLayout: Managerul de pași presetat aplică pasul TrivialF2QLayout, care mapează trivial $n$ biți fermionici la $n$ qubiți.

• F2QSynth: Un pas de transpilare pentru a mapa instrucțiunile de circuit bazate pe fermioni la cele bazate pe qubiți.

qdrift = QDriftTrotterization(num_exc, rng=19)

pm = generate_preset_jw_pass_manager()
pm.optimization = FermionicPassManager([qdrift])

sqdrift_circuits = []
for circ in init_circuits:
    sqdrift_circuits += (pm.run(circ) for _ in range(num_circuits))

for circ in sqdrift_circuits:
    circ.measure_all()

print(len(sqdrift_circuits))

400

Acum că am terminat cu optimizările la nivel fermionic, putem transpila circuitele pentru execuție pe simulator.

simulator = AerSimulator()
shots = 100

transpiled_circuits = transpile(sqdrift_circuits, simulator)

Pasul 3: Execuția folosind primitivele Qiskit

Acum că avem circuitele noastre, le putem rula folosind primitivele Qiskit pe AerSimulator. Vom combina toate numărătorile (counts) din diferite circuite. Le convertim în vectori booleeni înainte de a le post-procesa în final cu SQD.

print(
    f"Executing {len(transpiled_circuits)} circuits with {shots} shots each..."
)

job = simulator.run(transpiled_circuits, shots=shots)
result = job.result()

all_counts = [result.get_counts(i) for i in range(len(transpiled_circuits))]

print(len(all_counts), "length before post processing")

Executing 400 circuits with 100 shots each...
400 length before post processing

Pasul 4: Post-procesarea și returnarea rezultatului în formatul clasic dorit

Folosirea șirurilor de biți pentru SQD

Putem rula acum schema de diagonalizare pe șirurile de biți selectate pentru a găsi cea mai mică valoare proprie, care va corespunde energiei stării fundamentale a moleculei. Creăm o funcție de callback, declarăm ocupările inițiale și setăm parametrii înainte de a rula, în final, schema de diagonalizare. Funcția de callback este folosită pentru a afișa iterația curentă și estimarea curentă a valorii proprii la fiecare iterație.

În final, pentru a obține estimarea stării fundamentale, adăugăm nuclear_repulsion_energy la energia rezultată.

Notă: Dimensiunea subspațiului nu este fixă între iterații, chiar și pe simulatorul fără zgomot — fiecare subeșantion extrage un set diferit de configurații, iar pasul de recuperare remodelează grupul între iterații, astfel încât dimensiunea raportată variază de la un subeșantion la altul. Eșantionarea fără zgomot nu fixează prin ea însăși dimensiunea subspațiului selectat. Rularea pe hardware, însă, tinde să dea sistematic subspații mai mari, deoarece shot-urile zgomotoase rup simetria numărului de particule, iar recuperarea configurațiilor le transformă în vectori de bază suplimentari. Din acest motiv, vom introduce și un alt pas pentru eliminarea (pruning) șirurilor de biți în secțiunea despre hardware.

combined_counts = {}
for counts in all_counts:
    for bitstring, count in counts.items():
        combined_counts[bitstring] = combined_counts.get(bitstring, 0) + count

bit_array = BitArray.from_counts(combined_counts)
print(bit_array.num_shots)

print(f"  Alpha electrons: {num_elec_a}")
print(f"  Beta electrons: {num_elec_b}")
print(f"  Number of orbitals: {norb}")
print(f"  Number of spin orbitals (qubits): {2*norb}")
print(f"Integral shapes: h1e={h1e.shape}, h2e={h2e.shape}")

# SQD parameters
samples_per_batch = 300
num_batches = 3
max_iterations = 5

initial_occupancies = (
    np.array([1] * num_elec_a + [0] * (norb - num_elec_a)),  # alpha
    np.array([1] * num_elec_b + [0] * (norb - num_elec_b)),  # beta
)

result_history = []

def callback(results: list[SCIResult]):
    result_history.append(results)
    iteration = len(result_history)
    print(f"Iteration {iteration}")
    for i, result in enumerate(results):
        print(f"\tSubsample {i}")
        print(f"\t\tEnergy: {result.energy + nuclear_repulsion_energy}")
        print(
            f"\t\tSubspace dimension: {np.prod(result.sci_state.amplitudes.shape)}"
        )

# Run SQD with configuration recovery
print("\nRunning SQD with configuration recovery...")
result = diagonalize_fermionic_hamiltonian(
    h1e,
    h2e,
    bit_array,
    samples_per_batch=samples_per_batch,
    norb=norb,
    nelec=(num_elec_a, num_elec_b),
    num_batches=num_batches,
    energy_tol=1e-3,
    occupancies_tol=1e-3,
    max_iterations=max_iterations,
    initial_occupancies=initial_occupancies,
    seed=42,
    callback=callback,
)

computed_energy = result.energy + nuclear_repulsion_energy

print("FINAL SQD RESULTS")
print(f"Orbital occupancies (alpha): {result.orbital_occupancies[0]}")
print(f"Orbital occupancies (beta): {result.orbital_occupancies[1]}")

energy_error = abs(computed_energy - reference_energy)
print(f"Reference Energy: {reference_energy:.10f} Ha")
print(f"Computed Energy:  {computed_energy:.10f} Ha")
print(f"Error:            {energy_error:.10e} Ha")

40000
  Alpha electrons: 7
  Beta electrons: 7
  Number of orbitals: 10
  Number of spin orbitals (qubits): 20
Integral shapes: h1e=(10, 10), h2e=(10, 10, 10, 10)

Running SQD with configuration recovery...
Iteration 1
	Subsample 0
		Energy: -107.64767025226178
		Subspace dimension: 5538
	Subsample 1
		Energy: -107.64772799119115
		Subspace dimension: 5670
	Subsample 2
		Energy: -107.64765512281548
		Subspace dimension: 5767
Iteration 2
	Subsample 0
		Energy: -107.64795948524682
		Subspace dimension: 6080
	Subsample 1
		Energy: -107.64806617355072
		Subspace dimension: 6300
	Subsample 2
		Energy: -107.64802260640258
		Subspace dimension: 6308
FINAL SQD RESULTS
Orbital occupancies (alpha): [0.99999464 0.99999643 0.99584631 0.99332984 0.96684652 0.96686712
 0.99301927 0.0373282  0.0373266  0.00944508]
Orbital occupancies (beta): [0.99999462 0.99999643 0.9958261  0.99332349 0.96684268 0.96686737
 0.99302145 0.03733536 0.03733399 0.0094585 ]
Reference Energy: -107.6481842917 Ha
Computed Energy:  -107.6480661736 Ha
Error:            1.1811817564e-04 Ha

Exemplu pe hardware

Acest exemplu folosește 20 de qubiți (10 orbitali spațiali). Această alegere este o comoditate pentru un tutorial care trebuie să ruleze rapid, nu o limită strictă a metodei.

Costul pasului clasic nu este stabilit direct de numărul de qubiți. SQD diagonalizează hamiltonianul proiectat pe subspațiul generat de configurațiile eșantionate, astfel încât ceea ce determină costul clasic este dimensiunea acelui subspațiu selectat — guvernată aici de samples_per_batch, num_batches și de câte configurații distincte produc de fapt circuitele — împreună cu algebra liniară rară necesară pentru a aplica hamiltonianul proiectat. Spațiul CI complet crește combinatorial cu numărul de orbitali și electroni, dar subspațiul selectat este o felie mică, reglabilă din acesta, iar noi îi controlăm direct dimensiunea. Prin urmare, numărul de qubiți și dificultatea clasică pot varia oarecum independent: un spațiu orbital mai larg eșantionat într-un subspațiu modest poate fi mai ieftin decât un sistem mai mic diagonalizat pe unul foarte mare.

Prin urmare, în practică, dimensiunea fezabilă a sistemului depinde de dimensiunea subspațiului de care ai nevoie pentru precizia dorită și de memoria și nucleele disponibile pentru rezolvatorul de valori proprii. Spațiile orbitale mai mari necesită de obicei un subspațiu mai mare pentru a atinge precizia chimică, iar acest lucru este ceea ce motivează în cele din urmă resursele distribuite — vezi qiskit-addon-sqd-hpc pentru extinderea acestui pas. În loc să presupui un prag fix, abordarea practică este să urmărești dimensiunea raportată a subspațiului și convergența energiei de-a lungul iterațiilor și să crești dimensiunea subspațiului până când energia încetează să se îmbunătățească sau epuizezi memoria disponibilă.

Notă: Din cauza erorii de eșantionare provocate de zgomotul din hardware, subspațiul creat pentru diagonalizare în rularea pe hardware va fi mai mare decât cel obținut atunci când folosim simulatorul. Deși acest lucru crește dimensiunea subspațiului pe care vrem să îl diagonalizăm, fluxul de lucru ne oferă totuși un răspuns precis, datorită robusteții SQD la zgomot.

Eliminarea șirurilor false

Aici putem alege să efectuăm un pas suplimentar. Când avem toate șirurile de biți din execuțiile circuitelor, putem fie să filtrăm șirurile de biți invalide înainte de a rula SQD, fie să mergem mai departe fără eliminare (pruning). Omiterea eliminării este în general preferabilă pentru rulările pe hardware, deoarece lasă măsurătorile (shots) cu simetrie ruptă disponibile pentru recuperarea configurațiilor, care le poate repara în configurații valide și astfel poate lărgi subspațiul în loc să elimine acele măsurători direct.

Deoarece azotul poate avea doar șapte electroni $\alpha$ și șapte $\beta$, orice șir de biți care are mai mult sau mai puțin de șapte de 1 în prima și în a doua jumătate a rezultatului poate fi eliminat. Definim o funcție care verifică dacă șirurile de biți sunt valide și, dacă nu sunt, le elimină. Odată ce filtrăm șirurile false, restul sunt trimise în schema de diagonalizare. Folosește steagul PRUNE de mai jos pentru a comuta între cele două comportamente.

Ține minte că eliminarea este doar una dintre mai multele alegeri care modelează subspațiul final, alături de numărul de circuite, setul de timpi de evoluție și filtrarea termenilor diagonali. Compararea unei rulări cu eliminare cu una fără eliminare este informativă doar dacă restul rămâne fix; fila C++ discută acest lucru mai în detaliu, deoarece acolo se postselectează în loc să se recupereze și, de asemenea, diferă în privința celorlalți parametri.

name = "fcidump_files/N2_sto_3g"

fcidump = tools.fcidump.read(name)

# Extract metadata from the FCIDump header
norb = fcidump["NORB"]  # number of spatial orbitals
nelec = fcidump["NELEC"]  # total number of electrons
e_nuc = fcidump["ECORE"]  # nuclear repulsion / core energy
ms2 = fcidump["MS2"]  # 2S (spin)

num_elec_a = (nelec + ms2) // 2  # alpha electrons
num_elec_b = (nelec - ms2) // 2  # beta  electrons

# Reconstruct full 4-index ERIs from the FCIDump (stored in 8-fold symmetry)
h1e = fcidump["H1"]  # shape (norb, norb)
h2e = ao2mo.restore(  # shape (norb, norb, norb, norb)
    1, fcidump["H2"], norb
)

cisolver = fci.direct_spin1.FCI()
cisolver.max_cycle = 200
cisolver.conv_tol = 1e-12

e_fci, _ = cisolver.kernel(
    h1e,
    h2e,
    norb,
    (num_elec_a, num_elec_b),
    ecore=e_nuc,  # adds nuclear repulsion to the final energy
)

reference_energy = e_fci

print(f"Reference FCI Energy  = {reference_energy:.10f} Ha")

nuclear_repulsion_energy = fcidump["ECORE"]
print(f"Nuclear Repulsion Energy = {nuclear_repulsion_energy:.10f} Ha")

fcidump = FCIDump.from_file(name)
hamiltonian = FermionOperator.from_fcidump(fcidump)
num_modes = 2 * fcidump.norb

# Apply automatic grouping
canon = canonical_order(hamiltonian.normal_ordered().simplify(atol=1e-16))
exit_code = group_terms_by_electronic_structure(
    canon, num_modes, two_body_physicist_order=False
)
filter_diagonal_terms(canon)

print(len(canon.groups))

# SqDRIFT parameters
times = [1.0, 10.0]  # Total evolution times used for the subspace creation
num_exc = 10  # Number of excitation groups per circuit
num_circuits = 200  # Number of circuits to generate

init_circuits = []
hf_gate = InitializeModes.from_hartree_fock(norb, (num_elec_a, num_elec_b))

for time in times:
    evo_gate = Evolution(num_modes, canon, time)
    circ = FermionicCircuit(num_modes)
    circ.append(hf_gate, circ.modes)
    circ.append(evo_gate, circ.modes)
    init_circuits.append(circ)

# Calculate weights for sampling (one per group)
qdrift = QDriftTrotterization(num_exc, rng=19)

pm = generate_preset_jw_pass_manager()
pm.optimization = FermionicPassManager([qdrift])

sqdrift_circuits = []
for circ in init_circuits:
    sqdrift_circuits += (pm.run(circ) for _ in range(num_circuits))

for circ in sqdrift_circuits:
    circ.measure_all()

print(len(sqdrift_circuits))

# This example assumes you have saved your IBM Quantum Platform account locally.
service = QiskitRuntimeService(channel="ibm_quantum_platform")

# Select backend (choose based on qubit requirements)
backend = service.least_busy(
    operational=True,
    simulator=False,
    min_num_qubits=2 * norb,
)

print(f"Selected backend: {backend.name} ({backend.num_qubits} qubits)")

# Transpile for hardware
transpiled_circuits = transpile(
    sqdrift_circuits,
    backend=backend,
    optimization_level=3,
    seed_transpiler=42,
)

shots = 100

sampler = Sampler(mode=backend)

sampler.options.environment.job_tags = ["TUT_SQDRIFT"]

job = sampler.run(transpiled_circuits, shots=shots)
result = job.result()

# Extract counts from SamplerV2 results
all_counts = [pub_result.data.meas.get_counts() for pub_result in result]

# Set to True to filter out bitstrings that violate electron-number conservation
PRUNE = False

def is_valid_bitstring(
    bitstring: str, norb: int, nelec: tuple[int, int]
) -> bool:
    n_alpha, n_beta = nelec
    return (
        len(bitstring) == 2 * norb
        and bitstring[norb:].count("1") == n_alpha
        and bitstring[:norb].count("1") == n_beta
    )

if PRUNE:
    all_counts_filtered = []
    for counts in all_counts:
        filtered_count = {}
        for key in counts:
            if not is_valid_bitstring(key, norb, (num_elec_a, num_elec_b)):
                continue
            elif key not in filtered_count.keys():
                filtered_count[key] = counts[key]
            else:
                filtered_count[key] += counts[key]
        all_counts_filtered.append(filtered_count)
    all_counts = all_counts_filtered

combined_counts = {}
for counts in all_counts:
    for bitstring, count in counts.items():
        combined_counts[bitstring] = combined_counts.get(bitstring, 0) + count

bit_array = BitArray.from_counts(combined_counts)
print(bit_array.num_shots)

print("Electron configuration:")
print(f"  Total electrons: {nelec}")
print(f"  Alpha electrons: {num_elec_a}")
print(f"  Beta electrons: {num_elec_b}")
print(f"  Number of orbitals: {norb}")
print(f"  Number of spin orbitals (qubits): {2*norb}")
print(f"Integral shapes: h1e={h1e.shape}, h2e={h2e.shape}")

# SQD parameters
samples_per_batch = 300
num_batches = 3
max_iterations = 5

initial_occupancies = (
    np.array([1] * num_elec_a + [0] * (norb - num_elec_a)),  # alpha
    np.array([1] * num_elec_b + [0] * (norb - num_elec_b)),  # beta
)

result_history = []

def callback(results: list[SCIResult]):
    result_history.append(results)
    iteration = len(result_history)
    print(f"Iteration {iteration}")
    for i, result in enumerate(results):
        print(f"\tSubsample {i}")
        print(f"\t\tEnergy: {result.energy + nuclear_repulsion_energy}")
        print(
            f"\t\tSubspace dimension: {np.prod(result.sci_state.amplitudes.shape)}"
        )

# Run SQD with configuration recovery
print("\nRunning SQD with configuration recovery...")
result = diagonalize_fermionic_hamiltonian(
    h1e,
    h2e,
    bit_array,
    samples_per_batch=samples_per_batch,
    norb=norb,
    nelec=(num_elec_a, num_elec_b),
    num_batches=num_batches,
    energy_tol=1e-3,
    occupancies_tol=1e-3,
    max_iterations=max_iterations,
    initial_occupancies=initial_occupancies,
    seed=42,
    callback=callback,
)

computed_energy = result.energy + nuclear_repulsion_energy

print("FINAL SQD RESULTS")
print(f"Orbital occupancies (alpha): {result.orbital_occupancies[0]}")
print(f"Orbital occupancies (beta): {result.orbital_occupancies[1]}")

energy_error = abs(computed_energy - reference_energy)
print(f"Reference Energy: {reference_energy:.10f} Ha")
print(f"Computed Energy:  {computed_energy:.10f} Ha")
print(f"Error:            {energy_error:.10e} Ha")

Reference FCI Energy  = -107.6481842917 Ha
Nuclear Repulsion Energy = 23.7887003074 Ha
5060
400

Selected backend: ibm_aachen (156 qubits)
40000
Electron configuration:
  Total electrons: 14
  Alpha electrons: 7
  Beta electrons: 7
  Number of orbitals: 10
  Number of spin orbitals (qubits): 20
Integral shapes: h1e=(10, 10), h2e=(10, 10, 10, 10)

Running SQD with configuration recovery...
Iteration 1
	Subsample 0
		Energy: -107.64593072647523
		Subspace dimension: 7221
	Subsample 1
		Energy: -107.6458270048177
		Subspace dimension: 7209
	Subsample 2
		Energy: -107.64007673117075
		Subspace dimension: 7138
Iteration 2
	Subsample 0
		Energy: -107.64757372124944
		Subspace dimension: 9009
	Subsample 1
		Energy: -107.64674060104392
		Subspace dimension: 8245
	Subsample 2
		Energy: -107.64731360491942
		Subspace dimension: 8178
Iteration 3
	Subsample 0
		Energy: -107.64765518770588
		Subspace dimension: 8835
	Subsample 1
		Energy: -107.64767975712016
		Subspace dimension: 8649
	Subsample 2
		Energy: -107.64761634415606
		Subspace dimension: 8648
FINAL SQD RESULTS
Orbital occupancies (alpha): [0.99999504 0.9999964  0.99590318 0.9932359  0.96697158 0.96696295
 0.99298797 0.03728154 0.03728186 0.00938359]
Orbital occupancies (beta): [0.9999946  0.99999641 0.99590413 0.99323077 0.96697361 0.96696174
 0.99298424 0.03728121 0.03728169 0.00939159]
Reference Energy: -107.6481842917 Ha
Computed Energy:  -107.6476797571 Ha
Error:            5.0453460619e-04 Ha

C++

Obține codul C++

Sursa C++, proiectul CMake și instrucțiunile de compilare se află în directorul docs/tutorials/assets/sqdrift/cpp al depozitului de documentație Qiskit. Datele de intrare FCIDump se află în directorul vecin fcidump_files. Clonează depozitul, sau descarcă ambele directoare, pentru a compila și rula singur acest tutorial.

Estimare de utilizare: 180 de secunde pe un procesor Heron r3 (NOTĂ: Aceasta este doar o estimare. Timpul tău de execuție poate varia.)

Rezultate ale învățării (C++)

• Învață cum să creezi circuite cu adâncime mai mică în comparație cu Trotterizarea

• Parcurge un flux de lucru complet pentru estimarea stării fundamentale folosind qDRIFT și SQD

• Învață cum să folosești qiskit-fermions împreună cu alte addon-uri Qiskit pentru a implementa un astfel de flux de lucru

Cerințe preliminare (C++)

• Citește prezentarea generală a Diagonalizării cuantice bazate pe eșantioane (SQD)

• Citește lecția Diagonalizare cuantică Krylov bazată pe eșantioane (SKQD)

Context (C++)

SqDRIFT este o variantă a SKQD care elimină nevoia de a alege un ansatz din care să eșantionezi șiruri de biți, folosind în schimb un ansamblu de circuite de evoluție temporală construite direct din Hamiltonianul țintă. Acest lucru se realizează prin subeșantionarea unor operatori de evoluție temporală mai mici din Hamiltonian, pe baza coeficienților săi, cunoscută sub numele de metoda de Trotterizare qDRIFT.

Ambele implementări folosesc Qiskit Fermions pentru a reprezenta Hamiltonianul molecular și a construi circuite pentru algoritmul qDRIFT. Implementarea Python folosește pass-uri fermionice de Transpiler, în timp ce implementarea C++ eșantionează explicit grupuri și le mapează la operatori de qubiți.

Fie Hamiltonianul de forma:

$$H = \sum_{i=1}^{N} c_i h_i$$

unde, fără a pierde din generalitate, cerem ca $c_i > 0$ și ca cea mai mare valoare proprie a lui $h_i$ să fie egală, în valoare absolută, cu $1$. Orice prefactor cu semn sau complex este absorbit în $h_i$, deci coeficienții $c_i$ sunt ponderi strict pozitive, în timp ce $h_i$ poartă direcția fiecărui termen. Aici $N$ este numărul de termeni (sau, după grupare, numărul de grupuri) din hamiltonian; este o proprietate a hamiltonianului și este distinct de numărul de operatori eșantionați într-un singur circuit, notat $n$ mai jos.

Algoritmul qDRIFT realizează apoi, pentru timpul țintă $t$, un anumit operator $V_k$, unde $k$ merge de la $1 \cdots K$ și semnifică al $k_{th}$ circuit SqDRIFT, definit ca:

$$V_k = \prod_{j=1}^{n} e^{-i h_{k_j} \lambda t / n }$$

Aici $n$ este numărul de operatori eșantionați per circuit, iar $K$ este numărul de circuite din ansamblu. Produsul se face peste cele $n$ trageri, nu peste toți cei $N$ termeni ai hamiltonianului, iar deoarece termenii sunt trași cu revenire, același $h_i$ poate apărea de mai multe ori într-un singur $V_k$.

Cantitatea:

$$\lambda = \sum_{i=1}^{N} c_i$$

este norma $L_1$ a coeficienților, deci fiecare dintre cei $n$ pași evoluează pe aceeași durată $\lambda t / n$, indiferent de termenul care a fost tras. Uniformitatea unghiului de pas este caracteristica esențială a qDRIFT: un coeficient influențează rezultatul prin cât de des este tras termenul său, nu prin cât de mult este rotit acel termen. Indicii sunt eșantionați din distribuția:

$$P[k_i] = \frac{c_i}{\lambda}$$

deci seria $(k_1, \ldots, k_n)$ este o secvență aleatoare de indici de termeni trași din această distribuție. Deoarece $c_i$ sunt pozitivi și însumează $\lambda$, aceasta este o distribuție de probabilitate normalizată, iar valoarea așteptată a canalului rezultat, peste tragerile aleatoare, aproximează evoluția sub $H$, cu o eroare care scade pe măsură ce $n$ crește. Observă că eroarea de aproximare depinde de $\lambda$, nu de numărul de termeni $N$.

(Articolul SqDRIFT scrie numărul de termeni ca $\mathcal{N}$ și lungimea secvenței ca $N$; noi folosim $N$ și $n$ aici pentru a-i menține clar distincți.)

Acest tutorial arată cum să generezi un ansamblu de astfel de circuite randomizate. După ce am creat aceste circuite, similar modului în care creăm un subspațiu Krylov pentru operatori diferiți, eșantionăm șiruri de biți din mai mulți astfel de operatori cu parametri de timp diferiți. Aceasta asigură o suprapunere mai mare între vectorii stării fundamentale și șirurile de biți eșantionate.

Cerințe (C++)

Această implementare rulează pe hardware IBM Quantum® și folosește un executabil SBD extern pentru diagonalizarea clasică. Nu include exemplul cu simulatorul Python.

Urmează instrucțiunile de compilare C++ pentru a instala compilatorul, CMake, Rust, Boost și dependențele Qiskit. Instalează Rust cu rustup: dependențele Qiskit își fixează propria versiune de Rust, pe care rustup o descarcă automat, dar un Rust instalat din Homebrew sau dintr-un manager de pachete Linux ar putea fi prea vechi. Pasul SBD final necesită în plus MPI, OpenMP și BLAS/LAPACK. Windows nu este suportat în prezent de această implementare.

Folosește sursa C++ completă și proiectul CMake din depozit. Păstrează directorul adiacent fcidump_files când descarci exemplul. Dintr-o copie locală a depozitului de documentație, compilează cu următoarele comenzi. Prima compilare descarcă și compilează toate dependențele și durează câteva minute.

cd docs/tutorials/assets/sqdrift/cpp
cmake -S . -B build -DCMAKE_BUILD_TYPE=Release
cmake --build build --config Release

Compilarea scrie executabilul SqDRIFT în docs/tutorials/assets/sqdrift/cpp (rădăcina proiectului), nu în build/. Rulează ./SqDRIFT din rădăcina proiectului, pentru ca acesta să poată găsi ../fcidump_files/N2_sto_3g. Toate celelalte comenzi din această filă se rulează tot din rădăcina proiectului, dacă nu se specifică altfel. Dacă ./SqDRIFT eșuează cu o eroare de bibliotecă partajată lipsă, urmează instrucțiunile specifice platformei despre calea bibliotecilor din ghidul de compilare; pe macOS, executabilul rulează de obicei fără ele.

Configurare (C++)

Fragmentele de mai jos sunt extrase din programul SqDRIFT.cpp complet. Compilează și rulează acel fișier ca un singur program.

Înainte de a rula SqDRIFT, asigură-te că tokenul tău API IBM Quantum este configurat pentru IBM Cloud Platform.

Dacă ai Python cu qiskit-ibm-runtime instalat, cea mai simplă cale de a crea un fișier de configurare formatat corect este:

from qiskit_ibm_runtime import QiskitRuntimeService
QiskitRuntimeService.save_account(channel="ibm_quantum_platform", token="YOUR_IBM_QUANTUM_TOKEN_HERE")

Aceasta scrie ~/.qiskit/qiskit-ibm.json (Linux/macOS) sau %USERPROFILE%\.qiskit\qiskit-ibm.json (Windows) cu structura corectă. Dacă preferi să creezi fișierul manual:

{
  "default-ibm-quantum-platform": {
    "channel": "ibm_quantum_platform",
    "token": "YOUR_IBM_QUANTUM_TOKEN_HERE"
  }
}

Obține-ți cheia API de la IBM Quantum Platform.

Notă: Un fișier de configurare lipsă sau malformat provoacă la pornire un panic Rust unwrap() greu de înțeles, în loc de o eroare clară. Verifică încă o dată numele cheii (default-ibm-quantum-platform) și valoarea canalului (ibm_quantum_platform) dacă întâmpini o problemă.

Implementare (C++)

Pasul 1: Mapează intrările clasice la o problemă cuantică (C++)

Încarcă Hamiltonianul molecular

char filename[] = "../fcidump_files/N2_sto_3g";
QfFCIDump* fcidump = qf_fcidump_from_file(filename);
if (fcidump == nullptr) {
    std::cerr << "Failed to load FCIDump!" << std::endl;
    return 1;
}

uint32_t norb = qf_fcidump_norb(fcidump);
uint32_t nelec = qf_fcidump_nelec(fcidump);
uint32_t n_alpha = nelec / 2;
uint32_t n_beta  = nelec - n_alpha;
uint32_t num_modes = 2 * norb;
uint32_t num_qubits = num_modes;
std::cout << "✓ Loaded N2 molecule (" << norb << " orbitals, "
          << num_qubits << " qubits, " << nelec << " electrons)" << std::endl;

Output:

✓ Loaded N2 molecule (10 orbitals, 20 qubits, 14 electrons)

Încarcă Hamiltonianul moleculei $N_2$ dintr-un fișier FCIDump. Numărul de qubiți este egal cu dublul numărului de orbitali (spin sus și spin jos).

Creează Hamiltonianul și aplică ordonarea normală

// Create and normal order Hamiltonian
QfFermionOperator* hamiltonian = qf_ferm_op_from_fcidump(fcidump);
QfFermionOperator* normal_ordered = qf_ferm_op_normal_ordered(hamiltonian, nullptr);

Convertește datele FCIDump într-un operator fermionic și aplică ordonarea normală pentru a simplifica structura operatorului.

Grupează termenii după structura electronică

// Group terms by electronic structure
QfExitCode exit_code = qf_ferm_op_group_terms_by_electronic_structure(
    normal_ordered, num_modes, false
);
if (exit_code != QfExitCode_Success) {
    std::cerr << "Failed to group terms!" << std::endl;
    return 1;
}

uint32_t num_groups = qf_ferm_op_num_groups(normal_ordered);
std::cout << "Grouped into " << num_groups << " groups" << std::endl;

// Split into group operators
QfFermionOperator** group_ops = new QfFermionOperator*[num_groups];
qf_ferm_op_split_out_groups(normal_ordered, nullptr, 0, group_ops);

Output:

Grouped into 1902 groups

Grupează termenii înrudiți ai Hamiltonianului după structura electronică înainte de a construi circuitele lor de evoluție.

Calculează ponderile de eșantionare, apoi normalizează și mapează la operatori de qubiți

Ordinea contează aici. Normalizarea rescalează fiecare coeficient la magnitudine unitară pe loc, deci ponderile de eșantionare trebuie calculate mai întâi. Dacă normalizarea rulează prima, fiecare $|c_j|$ se citește ca exact 1.0, iar ponderea fiecărui grup se reduce la simplul număr al termenilor săi — distribuția de eșantionare ar fi atunci determinată de mărimea grupului, nu de coeficienții fizici.

Ponderea fiecărui grup este media coeficienților absoluți ai termenilor săi. Pasul QDriftTrotterization din fila Python folosește aceeași pondere pentru un Hamiltonian grupat. Folosirea sumei în schimb ar face ca $\lambda$ să crească odată cu mărimea grupului, ceea ce ar denatura timpul de evoluție $t_k = \lambda\tau/n$.

// Calculate sampling weights BEFORE normalization.
std::vector<double> weights(num_groups);
double total_weight = 0.0;

for (uint32_t i = 0; i < num_groups; i++) {
    QkComplex64* coeffs;
    uint64_t num_terms;
    qf_ferm_op_get_coeffs(group_ops[i], &coeffs, &num_terms);

    double abs_sum = 0.0;
    for (uint64_t j = 0; j < num_terms; j++) {
        abs_sum += std::sqrt(coeffs[j].re * coeffs[j].re +
                             coeffs[j].im * coeffs[j].im);
    }
    weights[i] = (num_terms > 0) ? (abs_sum / static_cast<double>(num_terms)) : 0.0;
    total_weight += weights[i];
}

std::cout << "Total weight (λ): " << total_weight << std::endl;

// Normalize the evolution operators term-by-term, now that the sampling
// weights above have already captured the original magnitudes.
for (uint32_t i = 0; i < num_groups; i++) {
    QkComplex64* coeffs;
    uint64_t num_terms;
    qf_ferm_op_get_coeffs(group_ops[i], &coeffs, &num_terms);

    for (uint64_t j = 0; j < num_terms; j++) {
        const double magnitude = std::sqrt(coeffs[j].re * coeffs[j].re +
                                           coeffs[j].im * coeffs[j].im);
        if (magnitude > 0.0) {
            coeffs[j].re /= magnitude;
            coeffs[j].im /= magnitude;
        }
    }
}

QkObs** qubit_ops = new QkObs*[num_groups];
for (uint32_t i = 0; i < num_groups; i++) {
    QfExitCode jw_exit = qf_ferm_op_jordan_wigner(group_ops[i], num_qubits, &qubit_ops[i]);
    if (jw_exit != QfExitCode_Success) {
        std::cerr << "Failed to map group " << i << std::endl;
        return 1;
    }
}
std::cout << "Mapped all " << num_groups << " normalized groups to qubit operators" << std::endl;

Output:

Total weight (λ): 337.261
Mapped all 1902 normalized groups to qubit operators

Ponderea $\lambda_i = \frac{1}{|G_i|}\sum_{j \in G_i} |c_{ij}|$ este media coeficienților absoluți ai termenilor grupului $G_i$, citiți din coeficienții originali. Normalizarea rescalează apoi fiecare termen fermionic la magnitudine unitară înainte de transformarea Jordan-Wigner, păstrând în operatorul evoluat doar faza/semnul său. Aceasta corespunde evoluției cu termeni normalizați din qDRIFT: magnitudinile coeficienților determină probabilitățile de eșantionare prin $\lambda$, în timp ce evoluția circuitului folosește termenul normalizat, astfel încât coeficienții mari să nu fie numărați de două ori.

Eșantionează seturile de operatori

// SqDRIFT Sampling: Create operator sets
    const int num_circuits = 100;    // Number of circuits to create
    const int ops_per_circuit = 10;  // Operators per circuit
    const double time_step = 1;    // Time step for evolution
    
    std::cout << "\n Generating " << num_circuits << " circuits with " 
              << ops_per_circuit << " operators each..." << std::endl;
    
    std::mt19937 gen(42);
    std::discrete_distribution<> dist(weights.begin(), weights.end());
    
    // Store sampled operator indices
    std::vector<std::vector<int>> operator_sets(num_circuits);
    for (int i = 0; i < num_circuits; i++) {
        operator_sets[i].resize(ops_per_circuit);
        for (int j = 0; j < ops_per_circuit; j++) {
            operator_sets[i][j] = dist(gen);
        }
    }

Output:

 Generating 100 circuits with 10 operators each...

Extrage seturile de operatori pentru 100 de circuite aleatoare, fiecare conținând 10 grupuri eșantionate din distribuția ponderată. Pasul următor transformă aceste seturi în ansamblul de circuite qDRIFT.

Pasul 2: Construiește circuitele qDRIFT (C++)

În versiunea C++, acest pas construiește circuitele din seturile de operatori eșantionate. Transpilarea pentru Backend-ul țintă are loc în Pasul 3, imediat înainte ca fiecare circuit să fie trimis.

    // Create Suzuki-Trotter circuits by composing all operators
    std::vector<Qiskit::circuit::QuantumCircuit> circuits;
    circuits.reserve(num_circuits);

    const uint32_t trotter_order = 1;
    const uint32_t trotter_reps = 1;
    const bool preserve_order = false;
    const bool insert_barriers = false;

    // Gate name mapping is static — compute once before all loops
    auto name_map = Qiskit::circuit::get_standard_gate_name_mapping();

    for (int circuit_idx = 0; circuit_idx < num_circuits; circuit_idx++) {
        // Create initial circuit with state preparation using C++ API
        Qiskit::circuit::QuantumCircuit qc(num_qubits, num_qubits);

        // State preparation: apply X gates to spin-up qubits 0..n_alpha-1
        // and spin-down qubits norb..norb+n_beta-1 (derived from FCIDump NELEC)
        for (uint32_t q = 0; q < n_alpha && q < num_qubits; q++) {
            qc.x(q);
        }
        for (uint32_t q = norb; q < norb + n_beta && q < num_qubits; q++) {
            qc.x(q);
        }

        // For each operator in this circuit's set
        for (int op_idx = 0; op_idx < ops_per_circuit; op_idx++) {
            int group_idx = operator_sets[circuit_idx][op_idx];
            double evolution_time = (total_weight * time_step) / ops_per_circuit;
            
            // Create Suzuki-Trotter evolution circuit
            QkCircuit* evolution_circuit = qk_circuit_library_suzuki_trotter(
                qubit_ops[group_idx], trotter_order, trotter_reps,
                evolution_time, preserve_order, insert_barriers
            );
            
            if (evolution_circuit != nullptr) {
                // Get number of instructions
                uint32_t num_ops = qk_circuit_num_instructions(evolution_circuit);
                
                // Manually append each instruction
                for (uint32_t i = 0; i < num_ops; i++) {
                    QkCircuitInstruction inst;
                    qk_circuit_get_instruction(evolution_circuit, i, &inst);

                    // Prepare qubit vector
                    std::vector<uint32_t> qubits(inst.num_qubits);
                    for (uint32_t j = 0; j < inst.num_qubits; j++) {
                        qubits[j] = inst.qubits[j];
                    }

                    // Get operation kind
                    QkOperationKind kind = qk_circuit_instruction_kind(evolution_circuit, i);

                    // Get mutable pointer to circuit
                    QkCircuit* mutable_circuit = qc.get_rust_circuit().get();

                    // Append based on operation type
                    if (kind == QkOperationKind_Gate) {
                        // Convert gate name string to QkGate enum
                        std::string gate_name(inst.name);
                        auto gate_it = name_map.find(gate_name);
                        if (gate_it != name_map.end()) {
                            QkGate gate_enum = gate_it->second.gate_map();
                            qk_circuit_parameterized_gate(
                                mutable_circuit,
                                gate_enum,           // Use QkGate enum, not string
                                qubits.data(),
                                inst.params
                            );
                        } else {
                            std::cerr << "Unknown gate in evolution circuit: " << gate_name << std::endl;
                            qk_circuit_instruction_clear(&inst);
                            qk_circuit_free(evolution_circuit);
                            return 1;
                        }
                    } else if (kind == QkOperationKind_Barrier) {
                        qk_circuit_barrier(
                            mutable_circuit,
                            qubits.data(), 
                            inst.num_qubits
                        );
                    } else if (kind == QkOperationKind_Reset) {
                        qk_circuit_reset(
                            mutable_circuit,
                            qubits[0]
                        );
                    }

                    qk_circuit_instruction_clear(&inst);
                }

                qk_circuit_free(evolution_circuit);
            }
        }

        circuits.push_back(qc);

        // Display progress
        if ((circuit_idx + 1) % 10 == 0) {
            std::cout << "  Created " << (circuit_idx + 1) << "/" 
                    << num_circuits << " circuits" << std::endl;
        }
    }

    std::cout << "✓ Created all " << num_circuits << " Suzuki-Trotter circuits" << std::endl;

Output:

  Created 10/100 circuits
  Created 20/100 circuits
  Created 30/100 circuits
  Created 40/100 circuits
  Created 50/100 circuits
  Created 60/100 circuits
  Created 70/100 circuits
  Created 80/100 circuits
  Created 90/100 circuits
  Created 100/100 circuits
✓ Created all 100 Suzuki-Trotter circuits

Programul complet afișează apoi statistici ale circuitelor și indicii operatorilor circuitului 0, care sunt omiși aici.

Pentru fiecare operator eșantionat, creează un circuit de evoluție Suzuki-Trotter cu timpul $t_k = \frac{\lambda \tau}{n}$, unde $n$ este numărul de operatori eșantionați per circuit și $\tau$ este pasul de timp. Deoarece fiecare grup eșantionat a fost normalizat termen cu termen înainte de mapare, această evoluție aplică în timpul rotației doar faza/semnul fiecărui termen; magnitudinile contribuie doar prin distribuția de eșantionare. Starea inițială Hartree-Fock este pregătită aplicând porți X pe qubiții 0-6 și 10-16 (șapte electroni în fiecare sector de spin). Instrucțiunile fiecărui circuit de evoluție sunt adăugate manual la circuitul principal folosind API-ul C++.

Pasul 3: Execuția folosind primitivele Qiskit (C++)

std::cout << "\n Connecting to IBM Quantum Cloud..." << std::endl;

// Initialize IBM Quantum Compute Service using qiskit-cpp
Qiskit::service::QiskitRuntimeService service;
std::cout << "✓ Connected to IBM Quantum" << std::endl;

// Get backend
const std::string backend_name = "ibm_fez";
auto backend = service.backend(backend_name);

if (backend.name().empty()) {
    std::cerr << "Backend " << backend_name
              << " is not available to this account." << std::endl;
    return 1;
}

std::cout << "✓ Selected backend: " << backend.name() << std::endl;

// Create sampler primitive
const int32_t shots = 100;  // Number of shots per circuit
Qiskit::primitives::BackendSamplerV2 sampler(backend, shots);

std::cout << "\n Submitting " << num_circuits << " circuits ("
            << shots << " shots each)..." << std::endl;

// Submit one job per circuit.
std::vector<std::shared_ptr<Qiskit::primitives::BasePrimitiveJob>> jobs;

for (int i = 0; i < num_circuits; i++) {
    for (uint32_t q = 0; q < num_qubits; q++) {
        circuits[i].measure(q, q);
    }
    Qiskit::circuit::QuantumCircuit transpiled_qc =
        Qiskit::compiler::transpile(circuits[i], backend, 2, 1.0, 42);
    std::vector<Qiskit::primitives::SamplerPub> pubs;
    pubs.push_back(Qiskit::primitives::SamplerPub(transpiled_qc, shots));
    jobs.push_back(sampler.run(pubs));
}

std::cout << "\n⏳ Waiting for " << num_circuits << " jobs to complete..." << std::endl;

// Collect results — poll each job and extract bitstrings
std::vector<boost::dynamic_bitset<>> all_bitstrings;

for (int i = 0; i < num_circuits; i++) {
    while (!jobs[i]->in_final_state()) {
        sleep(5);
    }
    if (jobs[i]->done()) {
        auto result = jobs[i]->result();
        if (result.size() > 0) {
            auto& pub_result = result[0];
            auto& bit_array  = pub_result.data();
            for (size_t s = 0; s < bit_array.num_shots(); s++) {
                auto sample = bit_array[s].to_string();
                boost::dynamic_bitset<> bs(num_qubits);
                for (uint32_t bit = 0; bit < num_qubits; bit++) {
                    bs[bit] = (sample[num_qubits - 1 - bit] == '1');
                }
                all_bitstrings.push_back(bs);
            }
        }
    } else {
        std::cout << "  Circuit " << (i+1) << " did not complete successfully" << std::endl;
    }
}

std::cout << "\n✓ Collected " << all_bitstrings.size()
            << " total bitstrings from all circuits" << std::endl;

Se conectează la IBM Quantum, selectează Backend-ul după nume, transpilează fiecare circuit pentru acel Backend și trimite un job per circuit: 100 de joburi Sampler separate, de câte 100 de măsurători (shots) fiecare, pe un registru de 20 de qubiți. Apoi interoghează periodic fiecare job până se finalizează și colectează șirurile de biți din result[0]. Fiecare job consumă timp QPU și așteaptă separat în coadă, deci te poți aștepta ca acest pas să dureze mult mai mult decât un singur job în lot.

Rezultatele pe hardware citate în această pagină provin dintr-o rulare pe ibm_pittsburgh (100 de circuite × 100 de măsurători), cu backend_name schimbat față de valoarea implicită ibm_fez. Eșantioanele salvate au fost postprocesate cu ordinea corectată a biților determinanților, arătată în Pasul 4. Zgomotul hardware variază de la o rulare la alta, deci numărătorile și energia finală ale tale vor diferi.

Output:

 Connecting to IBM Quantum Cloud...
✓ Connected to IBM Quantum
✓ Selected backend: ibm_pittsburgh

 Submitting 100 circuits (100 shots each)...

⏳ Waiting for 100 jobs to complete...

✓ Collected 10000 total bitstrings from all circuits

Schimbarea Backend-ului: Schimbă backend_name în SqDRIFT.cpp pentru a viza un alt dispozitiv; service.backends() listează ce poate accesa contul tău.

Pasul 4: Post-procesarea și returnarea rezultatului în formatul clasic dorit (C++)

Postselectează șirurile de biți și scrie șirurile CI

// Postselect bitstrings by Hamming weight (derived from FCIDump NELEC)
std::vector<double> bitstring_weights(all_bitstrings.size(), 1.0);
auto [filtered_bitstrings, filtered_weights] = Qiskit::addon::sqd::postselect_bitstrings(
    all_bitstrings,
    bitstring_weights,
    Qiskit::addon::sqd::MatchesRightLeftHamming<uint32_t>(n_alpha, n_beta)
);
std::ignore = filtered_weights;

std::cout << "✓ Postselected " << filtered_bitstrings.size()
          << " bitstrings with Hamming weight ("
          << n_alpha << "," << n_beta << ")" << std::endl;

// Convert bitstrings to CI strings using SQD addon
auto ci_strings = Qiskit::addon::sqd::bitstrings_to_ci_strings_symmetrize_spin(
    filtered_bitstrings,
    std::nullopt  // No dimension limit
);

std::cout << "✓ Generated " << ci_strings.size() << " CI strings" << std::endl;

// Write CI strings to file for SBD
std::ofstream alpha_file("alphadets_from_sqd.txt");
for (const auto& ci_string : ci_strings) {
    std::string bitstr;
    // SBD expects the highest orbital index first (bit zero on the right).
    boost::to_string(ci_string, bitstr);
    alpha_file << bitstr << "\n";
}
alpha_file.close();

std::cout << "✓ Wrote " << ci_strings.size()
          << " CI strings to alphadets_from_sqd.txt" << std::endl;
std::cout << " Ready for SBD diagonalization!" << std::endl;

Output:

✓ Postselected 1137 bitstrings with Hamming weight (7,7)
✓ Generated 49 CI strings
✓ Wrote 49 CI strings to alphadets_from_sqd.txt
 Ready for SBD diagonalization!

Postselectează șirurile de biți cu ponderea Hamming corectă (n_alpha electroni cu spin sus, n_beta electroni cu spin jos, derivate din câmpul NELEC al FCIDump), le convertește în șiruri de interacțiune de configurații (CI) folosind simetrizarea de spin și le scrie în alphadets_from_sqd.txt pentru Diagonalizarea ulterioară a Bazei Selectate.

Eliberează resursele

// Cleanup
// circuits vector will be automatically cleaned up (RAII)
// Free qubit operators
for (uint32_t i = 0; i < num_groups; i++) {
    qk_obs_free(qubit_ops[i]);
}
delete[] qubit_ops;
// Free fermionic operators
for (uint32_t i = 0; i < num_groups; i++) {
    qf_ferm_op_free(group_ops[i]);
}
delete[] group_ops;
qf_ferm_op_free(normal_ordered);
qf_ferm_op_free(hamiltonian);
qf_fcidump_free(fcidump);
return 0;

Eliberează corect toată memoria și toate resursele alocate. std::vector<QuantumCircuit> din C++ este curățat automat prin RAII. Operatorii de qubiți, operatorii fermionici și datele FCIDump sunt eliberate manual.

Acum că am creat baza pe care să proiectăm Hamiltonianul, putem continua cu procesul de diagonalizare pentru a obține estimarea stării fundamentale. Rulează următoarele comenzi din rădăcina proiectului, docs/tutorials/assets/sqdrift/cpp.

Compilează binarul SBD diag

CMake descarcă SBD doar ca sursă în deps/sbd/ (compilarea principală îl omite intenționat). Aplicația SBD folosește un Makefile simplu, iar SBD însuși necesită MPI, OpenMP și BLAS/LAPACK — deci un simplu make nu este suficient până când acestea sunt instalate și fișierul Configuration al aplicației corespunde lanțului tău de instrumente.

Makefile citește opțiunile de compilare și de link din fișierul Configuration din același director. Versiunea livrată în amonte este ajustată pentru o structură specifică macOS/Homebrew și nu se va compila ca atare pe o mașină standard — mai ales că -I/opt/homebrew/opt/llvm/include nu conține omp.h (acesta se află într-un director de resurse clang specific versiunii), iar clang++ de la Apple respinge categoric -fopenmp:

clang++: error: unsupported option '-fopenmp'

Rețeta următoare este verificată ca funcțională pe macOS (Apple Silicon, Darwin 25.6, Open MPI 5.x + Homebrew LLVM).

macOS (Apple Silicon) — testat

# 1. Install MPI and an OpenMP-capable compiler.
#    Apple clang cannot do OpenMP, so Homebrew's LLVM provides both
#    clang++ and libomp. BLAS/LAPACK come from Apple's Accelerate
#    framework, which is already part of macOS — no openblas needed.
brew install open-mpi llvm

Apoi înlocuiește deps/sbd/apps/chemistry_tpb_selected_basis_diagonalization/Configuration cu:

# Path to the SBD library
SBD_PATH=../..

# MPI C++ compiler. OMPI_CXX points Open MPI's wrapper at Homebrew's clang++,
# which (unlike Apple clang) supports -fopenmp.
CCCOM=OMPI_CXX=/opt/homebrew/opt/llvm/bin/clang++ mpicxx

# Build flags. omp.h lives in a version-specific clang resource dir, so let
# Homebrew's clang++ find it itself rather than hardcoding an -I path.
CCFLAGS= -std=c++17 -fopenmp -O3

# Link flags: OpenMP runtime from Homebrew LLVM, BLAS/LAPACK from the
# Accelerate framework that ships with macOS (no openblas install needed).
SYSLIB= -L/opt/homebrew/opt/llvm/lib -lomp -framework Accelerate

# 2. Build
cd deps/sbd/apps/chemistry_tpb_selected_basis_diagonalization
make
cd ../../../..

O compilare reușită afișează cele două linii de compilare/link și lasă un binar diag lângă Makefile. Un singur ld: warning: ignoring duplicate libraries: '-lomp' este inofensiv.

macOS Intel: Prefixul Homebrew este /usr/local, nu /opt/homebrew, deci înlocuiește-l în ambele căi anterioare.

Linux (Ubuntu/Debian) — netestat

Următoarele cerințe prealabile acoperă MPI, OpenMP și BLAS/LAPACK. GCC suportă nativ -fopenmp, deci nu este nevoie de nicio suprascriere a compilatorului. Acest lucru nu a fost verificat pe Linux:

sudo apt install -y libopenmpi-dev libomp-dev libblas-dev liblapack-dev

SBD_PATH=../..
CCCOM=mpicxx
CCFLAGS= -std=c++17 -fopenmp -O3
SYSLIB= -llapack -lblas

Rulează SqDRIFT pentru a genera șirurile CI

Dacă ai rulat deja ./SqDRIFT prin Pașii 1–4, alphadets_from_sqd.txt există deja în rădăcina proiectului și poți sări peste acest pas. Altfel, rulează-l din rădăcina proiectului:

./SqDRIFT
# Writes alphadets_from_sqd.txt in the project root

Rulează diagonalizarea

cd deps/sbd/apps/chemistry_tpb_selected_basis_diagonalization

# The FCIDump lives in the shared sqdrift/fcidump_files/ directory, one level
# above this tutorial's project root, hence five levels up from here.
ln -sf ../../../../../fcidump_files/N2_sto_3g fcidump.txt

# alphadets_from_sqd.txt is written by ./SqDRIFT into the project root itself,
# which is four levels up.
ln -sf ../../../../alphadets_from_sqd.txt alphadets.txt

./diag \
--fcidump fcidump.txt \
--adetfile alphadets.txt \
--method 0 \
--iteration 100 \
--block 10 \
--tolerance 1e-8 \
| tee ../../../../sbd_output.txt

cd ../../../..

Parametri:

• --fcidump: Calea către fișierul FCIDump care conține Hamiltonianul
• --adetfile: Calea către fișierul de determinanți alfa (șirurile CI de la SqDRIFT)
• --method 0: Metoda de diagonalizare Davidson
• --iteration 100: Numărul maxim de iterații Davidson
• --block 10: Dimensiunea blocului pentru algoritmul Davidson
• --tolerance 1e-8: Toleranța de convergență pentru energie

 Elapsed time for helper construction 0.001512 (sec)
 Elapsed time for init 1e-06 (sec)
 Davidson iteration 0.0 (tol=0.1857493937033264): -107.493531425221
 Davidson iteration 0.1 (tol=0.0455699976495391): -107.5098150493379 -105.93368397386
 Davidson iteration 0.2 (tol=0.009291923190816917): -107.5108356028703 -106.1245977311231 -105.6823366883346
 Davidson iteration 0.3 (tol=0.001488674426419356): -107.5108711084406 -106.3219246016818 -105.7376912083761 -105.4392974334491
 Davidson iteration 0.4 (tol=0.0003188697654737439): -107.5108722237947 -106.5935152393201 -106.1108017566682 -105.4932621297472
 Davidson iteration 0.5 (tol=3.080273116687244e-05): -107.5108722587463 -106.7245930484005 -106.2313070085243 -105.6423237790271
 Davidson iteration 0.6 (tol=2.981303155550973e-06): -107.5108722590151 -106.8008116251279 -106.330411722324 -105.69489630447
 Davidson iteration 0.7 (tol=4.565075388900609e-07): -107.5108722590186 -106.8090569201597 -106.3345525819528 -106.0578325555555
 Davidson iteration 0.8 (tol=4.870689101320012e-08): -107.5108722590186 -106.8098939189322 -106.3973756614146 -106.1617813840531
 Davidson iteration 0.9 (tol=6.224330513865847e-09): -107.5108722590187 -106.8159869231416 -106.4315891597535 -106.2222021964013
 Elapsed time for davidson 0.031276 (sec)
 Elapsed time for diagonalization 0.031286 (sec)
 Elapsed time for mult 0.002511 (sec)
 Energy = -107.5108722590187
 Elapsed time for measurement 4.6e-05 (sec)
 Sample-based diagonalization: Energy = -107.5108722590187
 Sample-based diagonalization: density = [1.999989340109303,1.99999662423008,1.999199121034164,1.981735389885314,1.999825046204946,1.999847567448985,1.995347419300126,0.01136512984757364,0.009844780833945199,0.002849581105572152
 Sample-based diagonalization: carryover bitstrings = [], size = 0

Linia Energy = -107.5108722590187 este estimarea stării fundamentale pe subspațiul eșantionat.

Cum diferă de exemplul Python

Cele două exemple diferă în mai multe privințe simultan, nu doar prin modul în care tratează șirurile de biți cu simetrie greșită, deci nu ne așteptăm să producă energii identice:

	Fila C++	Fila Python
Șiruri de biți cu simetrie greșită	Eliminate de postselect_bitstrings după ponderea Hamming	Reparate prin recuperare iterativă a configurațiilor
Structura buclei	O singură trecere: eșantionare → postselecție → o singură diagonalizare	Buclă exterioară: diagonalizare, citirea ocupărilor orbitale, recuperarea configurațiilor, rediagonalizare
Numărul de circuite	100 de circuite (num_circuits = 100)	400 de circuite (200 per timp de evoluție)
Timpi de evoluție	Un singur $\tau$ (time_step = 1)	Doi timpi, times = [1.0, 10.0]
Termeni diagonali	Păstrați în operatorul eșantionat	Eliminați de filter_diagonal_terms, ceea ce schimbă distribuția de eșantionare și circuitele
Ordonarea termenilor	Doar grupați	canonical_order aplicat înainte de grupare pentru reproductibilitate
Solver de diagonalizare	SBD (diag, binar extern MPI/OpenMP)	Solver qiskit-addon-sqd în același proces
Backend	Fixat după nume (backend_name = "ibm_fez"; rularea de pe această pagină a folosit ibm_pittsburgh)	Backend operațional least_busy cu suficienți qubiți
Transpilare	Fiecare circuit separat, nivel de optimizare 2	Toate circuitele împreună, nivel de optimizare 3
Trimiterea joburilor	100 de joburi, câte un circuit fiecare	Un job Sampler care conține toate cele 400 de circuite
Rulare pe simulator	Niciuna; doar hardware	Exemplu cu simulator fără zgomot înainte de rularea pe hardware

Ambele versiuni folosesc același Hamiltonian, își grupează termenii după structura electronică, ponderează fiecare grup după media coeficienților absoluți ai termenilor săi, pornesc de la starea Hartree-Fock și eșantionează 10 grupuri și 100 de măsurători (shots) per circuit.

Fiecare diferență afectează subspațiul eșantionat, și deci energia, în mod independent. Ansamblul mai mare și al doilea timp de evoluție, mai lung, lărgesc amândouă configurațiile explorate de exemplul Python; eșantionarea la mai mult de un timp oferă SqDRIFT acoperirea sa de tip Krylov a spațiului de energie joasă. Filtrarea termenilor diagonali schimbă operatorii pe care eșantionatorul îi poate extrage și ponderile lor, deci cele două rulări eșantionează distribuții diferite.

Postselecția este un filtru strict, iar o măsurătoare cu numărul greșit de particule este eliminată pur și simplu. În rularea de mai sus, aceasta a aruncat cea mai mare parte a datelor: 1137 din 10000 de măsurători au trecut verificarea ponderii Hamming (7,7) (aproximativ 11%), reducându-se la doar 49 de $\alpha$-determinanți distincți. Recuperarea configurațiilor folosește în schimb ocupările orbitale medii dintr-o diagonalizare anterioară pentru a inversa biți și a repara astfel de măsurători, reciclând ceea ce postselecția ar elimina. Această reacție o face intrinsec iterativă: recover_configurations din qiskit-addon-sqd-hpc necesită un argument avg_occupancies care există doar după ce a rulat o diagonalizare.

Această rulare C++ produce un subspațiu de 49 de $\alpha$-determinanți și o energie de -107.5109 Ha. Spațiul complet al $\alpha$-determinanților pentru șapte electroni în 10 orbitali are doar $\binom{10}{7} = 120$ de șiruri, deci furnizarea tuturor celor 120 aceluiași binar diag dă energia FCI exactă -107.6482 Ha. Rularea se află așadar cu aproximativ 0.1373 Ha peste valoarea exactă, acoperind 49 din cei 120 de $\alpha$-determinanți. Deoarece SQD este variațional în subspațiul selectat, un subspațiu mai mic poate doar să ridice energia, deci subspațiul incomplet este cea mai directă explicație a diferenței.

Această verificare cu FCI exact este posibilă doar pentru că exemplul este mic; nu face parte din fluxul de lucru SqDRIFT, iar pentru dimensiunile de sistem la care SQD este cu adevărat necesar, nu există o astfel de referință. Aici servește doar pentru a cuantifica cât din spațiu a recuperat această rulare.

O comparație cu exemplul Python nu este controlată, deci o diferență nu ar trebui atribuită doar absenței recuperării configurațiilor; reține că numărul de circuite, timpul de evoluție unic, filtrarea termenilor diagonali și zgomotul hardware diferă toate simultan și toate modelează subspațiul.

Pentru a reduce diferența, iată câțiva pași următori posibili de încercat după finalizarea acestui tutorial: mărește shots, num_circuits sau ops_per_circuit; eșantionează la mai mult de un time_step; și adaugă o buclă de recuperare peste Qiskit::addon::sqd::recover_configurations.

Pașii următori

Recomandări

• Diagonalizarea cuantică Krylov bazată pe eșantioane a unui model fermionic pe rețea — un tutorial înrudit care folosește circuite de evoluție în timp în loc de un ansatz variațional.

• Diagonalizarea cuantică bazată pe eșantioane a unui Hamiltonian chimic — un tutorial despre cum să construiești un circuit local unitary cluster Jastrow (LUCJ) pentru simularea chimiei cuantice.

• Lucrarea SqDRIFT — literatura pe care se bazează acest tutorial. (Reține că unele dintre optimizările discutate în această lucrare sunt în prezent în curs de dezvoltare, iar acest tutorial poate fi modificat în viitor în funcție de evoluția bibliotecilor folosite.)