Lewati ke konten utama
Proyek ini masih dalam tahap beta. Menemukan bug? Buat isu GitHub. Punya pertanyaan atau ide? Tanya di Discussions — kami sangat menghargai masukan kamu!
Bagaimana kualitas terjemahan di halaman ini?

Meningkatkan estimasi energi model kisi Fermionik dengan SQD

Dalam tutorial ini kita mengimplementasikan sebuah pola Qiskit yang menunjukkan cara memproses sampel kuantum yang bising untuk menemukan pendekatan ground state dari Hamiltonian kisi Fermionik yang dikenal sebagai model Anderson impuritas tunggal. Kita akan mengikuti pendekatan diagonalisasi kuantum berbasis sampel untuk memproses sampel yang diambil dari sekumpulan 16-Qubit keadaan basis Krylov pada interval waktu yang semakin meningkat. Keadaan-keadaan ini direalisasikan pada perangkat kuantum menggunakan Trotterisasi evolusi waktu. Untuk memperhitungkan efek noise kuantum, teknik pemulihan konfigurasi digunakan. Dengan asumsi keadaan awal yang baik dan sparsitas dari ground state, pendekatan ini terbukti konvergen secara efisien.

Polanya dapat digambarkan dalam empat langkah:

  1. Langkah 1: Petakan ke masalah kuantum

    • Hasilkan sekumpulan keadaan basis Krylov (yaitu, Circuit evolusi waktu Trotterisasi) pada interval waktu yang semakin meningkat untuk mengestimasi ground state

  2. Langkah 2: Optimalkan masalah

    • Transpile Circuit untuk Backend

  3. Langkah 3: Jalankan eksperimen

    • Ambil sampel dari Circuit menggunakan primitif Sampler

  4. Langkah 4: Pasca-proses hasil

    • Loop pemulihan konfigurasi swakonsisten

      • Pasca-proses seluruh kumpulan sampel bitstring, menggunakan pengetahuan awal tentang jumlah partikel dan rata-rata hunian orbital yang dihitung pada iterasi terbaru

      • Buat batch subsampel secara probabilistik dari bitstring yang dipulihkan

      • Proyeksikan dan diagonalisasikan Hamiltonian kisi Fermionik pada setiap subruang yang diambil sampelnya

      • Simpan energi ground state minimum yang ditemukan di semua batch dan perbarui rata-rata hunian orbital

Langkah 1: Petakan masalah ke Circuit kuantum​

Pertama, kita akan membuat Hamiltonian satu-benda dan dua-benda yang menggambarkan model Anderson impuritas tunggal (SIAM) satu dimensi dengan 7 situs bath (8 elektron dalam 8 orbital). Model ini digunakan untuk menggambarkan impuritas magnetik yang tertanam dalam logam.

Kemudian kita akan membuat Circuit Trotter 16-Qubit yang digunakan untuk menghasilkan subruang Krylov kuantum.

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

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

n_bath = 7  # number of bath sites

V = 1  # hybridization amplitude
t = 1  # bath hopping amplitude
U = 10  # Impurity onsite repulsion
eps = -U / 2  # Chemical potential for the impurity

# Place the impurity on the first qubit
impurity_index = 0

# One body matrix elements in the "position" basis
h1e = -t * np.diag(np.ones(n_bath), k=1) - t * np.diag(np.ones(n_bath), k=-1)
h1e[impurity_index, impurity_index + 1] = -V
h1e[impurity_index + 1, impurity_index] = -V
h1e[impurity_index, impurity_index] = eps

# Two body matrix elements in the "position" basis
h2e = np.zeros((n_bath + 1, n_bath + 1, n_bath + 1, n_bath + 1))
h2e[impurity_index, impurity_index, impurity_index, impurity_index] = U

Selanjutnya, kita akan menghasilkan subruang Krylov kuantum dengan sekumpulan Circuit kuantum yang di-Trotterisasi. Di sini kita membuat fungsi pembantu untuk menghasilkan keadaan awal (referensi) serta evolusi waktu untuk bagian satu-benda dan dua-benda dari Hamiltonian. Untuk deskripsi yang lebih rinci tentang model ini dan bagaimana Circuit dirancang, silakan merujuk ke makalah tersebut.

import ffsim
import scipy
from qiskit import QuantumCircuit, QuantumRegister
from qiskit.circuit.library import CPhaseGate, XGate, XXPlusYYGate

n_modes = n_bath + 1
nelec = (n_modes // 2, n_modes // 2)

dt = 0.2
Utar = scipy.linalg.expm(-1j * dt * h1e)

# The reference state
def initial_state(q_circuit, norb, nocc):
    """Prepare an initial state."""
    for i in range(nocc):
        q_circuit.append(XGate(), [i])
        q_circuit.append(XGate(), [norb + i])
    rot = XXPlusYYGate(np.pi / 2, -np.pi / 2)

    for i in range(3):
        for j in range(nocc - i - 1, nocc + i, 2):
            q_circuit.append(rot, [j, j + 1])
            q_circuit.append(rot, [norb + j, norb + j + 1])
    q_circuit.append(rot, [j + 1, j + 2])
    q_circuit.append(rot, [norb + j + 1, norb + j + 2])

# The one-body time evolution
free_fermion_evolution = ffsim.qiskit.OrbitalRotationJW(n_modes, Utar)

# The two-body time evolution
def append_diagonal_evolution(dt, U, impurity_qubit, num_orb, q_circuit):
    """Append two-body time evolution to a quantum circuit."""
    if U != 0:
        q_circuit.append(
            CPhaseGate(-dt / 2 * U),
            [impurity_qubit, impurity_qubit + num_orb],
        )

Hasilkan d keadaan yang berevolusi terhadap waktu yang menentukan subruang Krylov kuantum. Di sini, kita telah memilih d=8. Error dari pengambilan sampel keadaan basis Krylov konvergen seiring meningkatnya d. Perlu dicatat bahwa kekhususan instans masalah ini memungkinkan kompilasi yang efisien dari evolusi satu-benda menggunakan OrbitalRotationJW; namun, secara umum, kita dapat menggunakan PauliEvolutionGate dari Qiskit untuk mengimplementasikan evolusi waktu Trotterisasi dari Hamiltonian lengkap.

# Generate the initial state
qubits = QuantumRegister(2 * n_modes, name="q")
init_state = QuantumCircuit(qubits)
initial_state(init_state, n_modes, n_modes // 2)
init_state.draw("mpl", scale=0.4, fold=-1)

d = 8  # Number of Krylov basis states
circuits = []
for i in range(d):
    circ = init_state.copy()
    circuits.append(circ)
    for _ in range(i):
        append_diagonal_evolution(dt, U, impurity_index, n_modes, circ)
        circ.append(free_fermion_evolution, qubits)
        append_diagonal_evolution(dt, U, impurity_index, n_modes, circ)
    circ.measure_all()
circuits[0].draw("mpl", scale=0.4, fold=-1)

Quantum circuit diagram

circuits[-1].draw("mpl", scale=0.4, fold=-1)

Quantum circuit diagram

Langkah 2: Optimalkan masalah​

Setelah kita membuat Circuit Trotterisasi, kita akan mengoptimalkannya untuk perangkat keras target. Kita perlu memilih perangkat keras yang akan digunakan sebelum optimasi. Kita akan menggunakan Backend 127-Qubit palsu dari qiskit_ibm_runtime untuk mengemulasi perangkat nyata. Untuk menjalankan pada QPU nyata, cukup ganti Backend palsu dengan Backend nyata. Lihat dokumentasi Qiskit IBM Runtime untuk info lebih lanjut.

from qiskit_ibm_runtime.fake_provider.backends import FakeSherbrooke

backend = FakeSherbrooke()

Selanjutnya, kita akan men-transpile Circuit ke Backend target menggunakan Qiskit.

from qiskit.transpiler import generate_preset_pass_manager

# The circuit needs to be transpiled to the AerSimulator target
pass_manager = generate_preset_pass_manager(optimization_level=3, backend=backend)
isa_circuits = pass_manager.run(circuits)

Langkah 3: Jalankan eksperimen​

Setelah mengoptimalkan Circuit untuk eksekusi perangkat keras, kita siap menjalankannya pada perangkat keras target dan mengumpulkan sampel untuk estimasi energi ground state. Di sini kita menggunakan SamplerV2 dari qiskit-ibm-runtime untuk mensimulasikan sampel bising yang diambil dari Backend ibm_sherbrooke. Kemudian kita menggabungkan counts dari setiap keadaan basis Krylov menjadi satu kamus counts dan memplot 20 bitstring yang paling sering diambil sampelnya.

Catatan: Qiskit Aer diperlukan untuk mensimulasikan sampel dari Circuit yang telah di-transpile.

from qiskit.primitives import BitArray
from qiskit.visualization import plot_histogram
from qiskit_ibm_runtime import SamplerV2 as Sampler

# Sample from the circuits
noisy_sampler = Sampler(backend, options={"simulator": {"seed_simulator": 24}})
job = noisy_sampler.run(isa_circuits, shots=500)

# Combine the counts from the individual Trotter circuits
bit_array = BitArray.concatenate_shots([result.data.meas for result in job.result()])

plot_histogram(bit_array.get_counts(), number_to_keep=20)

Quantum circuit diagram

Langkah 4: Pasca-proses hasil​

Sekarang, kita menjalankan algoritma SQD menggunakan fungsi diagonalize_fermionic_hamiltonian. Lihat dokumentasi API untuk penjelasan argumen fungsi ini.

from qiskit_addon_sqd.fermion import SCIResult, diagonalize_fermionic_hamiltonian

# List to capture intermediate results
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}")
        print(f"\t\tSubspace dimension: {np.prod(result.sci_state.amplitudes.shape)}")

rng = np.random.default_rng(24)
result = diagonalize_fermionic_hamiltonian(
    h1e,
    h2e,
    bit_array,
    samples_per_batch=300,
    norb=n_modes,
    nelec=nelec,
    num_batches=3,
    max_iterations=10,
    symmetrize_spin=True,
    callback=callback,
    seed=rng,
)
Iteration 1
	Subsample 0
		Energy: -13.257128325607997
		Subspace dimension: 3969
	Subsample 1
		Energy: -13.257128325607997
		Subspace dimension: 3969
	Subsample 2
		Energy: -13.257128325607997
		Subspace dimension: 3969
Iteration 2
	Subsample 0
		Energy: -13.319666127542039
		Subspace dimension: 4096
	Subsample 1
		Energy: -13.420534292304595
		Subspace dimension: 4624
	Subsample 2
		Energy: -9.136171594591085
		Subspace dimension: 4624
Iteration 3
	Subsample 0
		Energy: -13.422491814612833
		Subspace dimension: 4900
	Subsample 1
		Energy: -13.422491814612833
		Subspace dimension: 4900
	Subsample 2
		Energy: -13.422491814612833
		Subspace dimension: 4900
Iteration 4
	Subsample 0
		Energy: -13.422491814612833
		Subspace dimension: 4900
	Subsample 1
		Energy: -13.422491814612833
		Subspace dimension: 4900
	Subsample 2
		Energy: -13.422491814612833
		Subspace dimension: 4900

Sekarang, kita memplot hasilnya. Plot pertama menunjukkan bahwa setelah beberapa iterasi, kita mendapatkan energi ground state yang tepat.

Contoh ini cukup kecil sehingga kita bisa menjelajahi seluruh ruang Hilbert, seperti yang terlihat pada pernyataan print di atas. Ingat, seluruh ruang Hilbert memiliki dimensi (num_orbitals choose nelec_a) * (num_orbitals choose nelec_b). Jadi untuk masalah ini: (8 choose 4)**2 = 4900. Dimensi subruang meningkat karena adanya enhanced configuration recovery, dan juga karena algoritma SQD membawa konfigurasi penting dari satu iterasi ke iterasi berikutnya. Pada iterasi terakhir, kita melakukan diagonalisasi di seluruh ruang Hilbert.

Plot kedua menunjukkan rata-rata hunian setiap orbital spasial di semua solusi batch. Kita melihat bahwa dengan probabilitas tinggi setiap orbital mengandung satu elektron.

import matplotlib.pyplot as plt

exact_energy = -13.422491814605827
min_es = [min(result, key=lambda res: res.energy).energy for result in result_history]
min_id, min_e = min(enumerate(min_es), key=lambda x: x[1])

# Data for energies plot
x1 = range(len(result_history))
yt1 = list(np.arange(-13.5, -13.1, 0.1))
ytl = [f"{i:.1f}" for i in yt1]

# Data for avg spatial orbital occupancy
y2 = np.sum(result.orbital_occupancies, axis=0)
x2 = range(len(y2))

fig, axs = plt.subplots(1, 2, figsize=(12, 6))

# Plot energies
axs[0].plot(x1, min_es, label="SQD energy", marker="o")
axs[0].set_xticks(x1)
axs[0].set_xticklabels(x1)
axs[0].set_yticks(yt1)
axs[0].set_yticklabels(ytl)
axs[0].axhline(y=exact_energy, color="#BF5700", linestyle="--", label="Exact energy")
axs[0].set_title("Approximated Ground State Energy vs SQD Iterations")
axs[0].set_xlabel("Iteration Index", fontdict={"fontsize": 12})
axs[0].set_ylabel("Energy", fontdict={"fontsize": 12})
axs[0].legend()

# Plot orbital occupancy
axs[1].bar(x2, y2, width=0.8)
axs[1].set_xticks(x2)
axs[1].set_xticklabels(x2)
axs[1].set_title("Avg Occupancy per Spatial Orbital")
axs[1].set_xlabel("Orbital Index", fontdict={"fontsize": 12})
axs[1].set_ylabel("Avg Occupancy", fontdict={"fontsize": 12})

print(f"Exact energy: {exact_energy:.5f}")
print(f"SQD energy: {min_e:.5f}")
print(f"Absolute error: {abs(min_e - exact_energy):.5f}")
plt.tight_layout()
plt.show()
Exact energy: -13.42249
SQD energy: -13.42249
Absolute error: 0.00000

Plot output

Masalah situs atau terjemahan?Edit halaman ini di doQumentation