Penyepelasan kuantum berasaskan sampel bagi Hamiltonian kimia

Anggaran penggunaan: kurang dari satu minit pada pemproses Heron r2 (NOTA: Ini adalah anggaran sahaja. Masa jalan anda mungkin berbeza.)

Latar belakang

Dalam tutorial ini, kami menunjukkan cara memproses sampel kuantum yang bermasalah untuk mencari penghampiran kepada keadaan dasar molekul nitrogen $\text{N}_2$ pada panjang ikatan keseimbangan, menggunakan algoritma penyepelasan kuantum berasaskan sampel (SQD), untuk sampel yang diambil daripada Circuit kuantum 59-Qubit (52 Qubit sistem + 7 Qubit ancilla). Pelaksanaan berasaskan Qiskit tersedia dalam SQD Qiskit addons, maklumat lanjut boleh didapati dalam dokumen yang berkaitan dengan contoh mudah untuk memulakan. SQD ialah teknik untuk mencari nilai eigen dan vektor eigen bagi operator kuantum, seperti Hamiltonian sistem kuantum, menggunakan pengkomputeran kuantum dan pengkomputeran klasik teragih bersama-sama. Pengkomputeran klasik teragih digunakan untuk memproses sampel yang diperoleh daripada pemproses kuantum, dan untuk mengunjur dan menyepelasi Hamiltonian sasaran dalam subruang yang dibentang oleh sampel tersebut. Ini membolehkan SQD bersifat tahan terhadap sampel yang rosak akibat bunyi kuantum dan menangani Hamiltonian yang besar, seperti Hamiltonian kimia dengan jutaan sebutan interaksi, yang melampaui kemampuan sebarang kaedah penyepelasan tepat. Aliran kerja berasaskan SQD mempunyai langkah-langkah berikut:

1. Pilih ansatz Circuit dan gunakannya pada komputer kuantum kepada keadaan rujukan (dalam kes ini, keadaan Hartree-Fock).
2. Ambil sampel rentetan bit daripada keadaan kuantum yang terhasil.
3. Jalankan prosedur pemulihan konfigurasi kendiri-konsisten pada rentetan bit untuk mendapatkan penghampiran kepada keadaan dasar.

SQD diketahui berfungsi dengan baik apabila eigenstate sasaran adalah jarang: fungsi gelombang disokong dalam set keadaan asas $\mathcal{S} = \{|x\rangle \}$ yang saiznya tidak meningkat secara eksponen dengan saiz masalah.

Kimia kuantum

Sifat-sifat molekul sebahagian besarnya ditentukan oleh struktur elektron di dalamnya. Sebagai zarah fermionic, elektron boleh diterangkan menggunakan formalisme matematik yang dipanggil kuantisasi kedua. Ideanya ialah terdapat sejumlah orbital, yang masing-masing boleh kosong atau dihuni oleh fermion. Sistem dengan $N$ orbital diterangkan oleh set operator pemusnahan fermionic $\{\hat{a}_p\}_{p=1}^N$ yang memenuhi hubungan antikomutasi fermionic,

$$\begin{align*} \hat{a}_p \hat{a}_q + \hat{a}_q \hat{a}_p &= 0, \\ \hat{a}_p \hat{a}_q^\dagger + \hat{a}_q^\dagger \hat{a}_p &= \delta_{pq}. \end{align*}$$

Konjugat adjoin $\hat{a}_p^\dagger$ dipanggil operator penciptaan.

Setakat ini, penerangan kami belum mengambil kira spin, yang merupakan sifat asas fermion. Apabila mengambil kira spin, orbital datang berpasangan yang dipanggil orbital ruang. Setiap orbital ruang terdiri daripada dua orbital spin, satu berlabel "spin-$\alpha$" dan satu berlabel "spin-$\beta$". Kami kemudian menulis $\hat{a}_{p\sigma}$ untuk operator pemusnahan yang berkaitan dengan orbital spin dengan spin $\sigma$ ($\sigma \in \{\alpha, \beta\}$) dalam orbital ruang $p$. Jika kami mengambil $N$ sebagai bilangan orbital ruang, maka terdapat sejumlah $2N$ orbital spin. Ruang Hilbert sistem ini dibentang oleh $2^{2N}$ vektor asas ortonormal yang dilabelkan dengan rentetan bit dua bahagian $\lvert z \rangle = \lvert z_\beta z_\alpha \rangle = \lvert z_{\beta, N} \cdots z_{\beta, 1} z_{\alpha, N} \cdots z_{\alpha, 1} \rangle$.

Hamiltonian sistem molekul boleh ditulis sebagai

$$\hat{H} = \sum_{ \substack{pr\\\sigma} } h_{pr} \, \hat{a}^\dagger_{p\sigma} \hat{a}_{r\sigma} + \frac12 \sum_{ \substack{prqs\\\sigma\tau} } h_{prqs} \, \hat{a}^\dagger_{p\sigma} \hat{a}^\dagger_{q\tau} \hat{a}_{s\tau} \hat{a}_{r\sigma},$$

di mana $h_{pr}$ dan $h_{prqs}$ ialah nombor kompleks yang dipanggil integral molekul yang boleh dikira daripada spesifikasi molekul menggunakan program komputer. Dalam tutorial ini, kami mengira integral menggunakan pakej perisian PySCF.

Untuk butiran tentang bagaimana Hamiltonian molekul diterbitkan, rujuk buku teks kimia kuantum (contohnya, Modern Quantum Chemistry oleh Szabo dan Ostlund). Untuk penjelasan peringkat tinggi tentang cara masalah kimia kuantum dipetakan ke komputer kuantum, semak kuliah Mapping Problems to Qubits dari Qiskit Global Summer School 2024.

Ansatz kluster unitari tempatan Jastrow (LUCJ)

SQD memerlukan ansatz Circuit kuantum untuk mengambil sampel. Dalam tutorial ini, kami akan menggunakan ansatz kluster unitari tempatan Jastrow (LUCJ) kerana gabungan motivasi fizikal dan keramahan perkakasannya.

Ansatz LUCJ ialah bentuk khusus daripada ansatz kluster unitari umum Jastrow (UCJ), yang mempunyai bentuk

$$\lvert \Psi \rangle = \prod_{\mu=1}^L e^{\hat{K}_\mu} e^{i \hat{J}_\mu} e^{-\hat{K}_\mu} | \Phi_0 \rangle$$

di mana $\lvert \Phi_0 \rangle$ ialah keadaan rujukan, sering diambil sebagai keadaan Hartree-Fock, dan $\hat{K}_\mu$ serta $\hat{J}_\mu$ mempunyai bentuk

$$\hat{K}_\mu = \sum_{pq, \sigma} K_{pq}^\mu \, \hat{a}^\dagger_{p \sigma} \hat{a}^{\phantom{\dagger}}_{q \sigma} \;,\; \hat{J}_\mu = \sum_{pq, \sigma\tau} J_{pq,\sigma\tau}^\mu \, \hat{n}_{p \sigma} \hat{n}_{q \tau} \;,$$

di mana kami telah mentakrifkan operator nombor

$$\hat{n}_{p \sigma} = \hat{a}^\dagger_{p \sigma} \hat{a}^{\phantom{\dagger}}_{p \sigma}.$$

Operator $e^{\hat{K}_\mu}$ ialah putaran orbital, yang boleh dilaksanakan menggunakan algoritma yang diketahui dalam kedalaman linear dan menggunakan kesambungan linear. Melaksanakan sebutan $e^{i \mathcal{J}_k}$ bagi ansatz UCJ memerlukan sama ada kesambungan semua-ke-semua atau penggunaan rangkaian tukar fermionic, menjadikannya mencabar untuk pemproses kuantum pra-toleran-kesalahan yang bermasalah dengan kesambungan terhad. Idea ansatz UCJ tempatan adalah untuk mengenakan kekangan kejarangan pada matriks $\mathbf{J}^{\alpha\alpha}$ dan $\mathbf{J}^{\alpha\beta}$ yang membolehkan ia dilaksanakan dalam kedalaman malar pada topologi Qubit dengan kesambungan terhad. Kekangan ditentukan oleh senarai indeks yang menunjukkan entri matriks dalam segitiga atas yang dibenarkan bukan sifar (memandangkan matriks adalah simetri, hanya segitiga atas perlu ditentukan). Indeks ini boleh ditafsirkan sebagai pasangan orbital yang dibenarkan berinteraksi.

Sebagai contoh, pertimbangkan topologi Qubit kekisi persegi. Kami boleh meletakkan orbital $\alpha$ dan $\beta$ dalam baris selari pada kekisi, dengan sambungan antara baris ini membentuk "anak tangga" seperti ini:

Dengan persediaan ini, orbital berspin sama disambungkan dengan topologi garisan, manakala orbital berspin berbeza disambungkan apabila mereka berkongsi orbital ruang yang sama. Ini menghasilkan kekangan indeks berikut pada matriks $\mathbf{J}$:

$$\begin{align*} \mathbf{J}^{\alpha\alpha} &: \set{(p, p+1) \; , \; p = 0, \ldots, N-2} \\ \mathbf{J}^{\alpha\beta} &: \set{(p, p) \;, \; p = 0, \ldots, N-1} \end{align*}$$

Dengan kata lain, jika matriks $\mathbf{J}$ bukan sifar hanya pada indeks yang ditetapkan dalam segitiga atas, maka sebutan $e^{i \mathcal{J}_k}$ boleh dilaksanakan pada topologi persegi tanpa menggunakan sebarang gate tukar, dalam kedalaman malar. Sudah tentu, mengenakan kekangan sedemikian pada ansatz menjadikannya kurang ekspresif, jadi lebih banyak ulangan ansatz mungkin diperlukan.

Perkakasan IBM® mempunyai topologi Qubit kekisi segi enam berat, di mana kami boleh menggunakan corak "zig-zag", seperti yang digambarkan di bawah. Dalam corak ini, orbital berspin sama dipetakan ke Qubit dengan topologi garisan (bulatan merah dan biru), dan sambungan antara orbital berspin berbeza wujud pada setiap orbital ruang ke-4, dengan sambungan difasilitasi oleh Qubit ancilla (bulatan ungu). Dalam kes ini, kekangan indeks adalah

$$\begin{align*} \mathbf{J}^{\alpha\alpha} &: \set{(p, p+1) \; , \; p = 0, \ldots, N-2} \\ \mathbf{J}^{\alpha\beta} &: \set{(p, p) \;, \; p = 0, 4, 8, \ldots (p \leq N-1)} \end{align*}$$

Pemulihan konfigurasi kendiri-konsisten

Prosedur pemulihan konfigurasi kendiri-konsisten direka untuk mengekstrak sebanyak mungkin isyarat daripada sampel kuantum yang bermasalah. Kerana Hamiltonian molekul memelihara nombor zarah dan spin Z, masuk akal untuk memilih ansatz Circuit yang juga memelihara simetri ini. Apabila digunakan pada keadaan Hartree-Fock, keadaan yang terhasil mempunyai nombor zarah dan spin Z yang tetap dalam tetapan tanpa bunyi. Oleh itu, bahagian spin-$\alpha$ dan spin-$\beta$ bagi mana-mana rentetan bit yang disampling daripada keadaan ini harus mempunyai berat Hamming yang sama seperti dalam keadaan Hartree-Fock. Disebabkan kehadiran bunyi dalam pemproses kuantum semasa, sesetengah rentetan bit yang diukur akan melanggar sifat ini. Bentuk pemilihan pasca yang mudah adalah membuang rentetan bit ini, tetapi ini membazir kerana rentetan bit tersebut mungkin masih mengandungi beberapa isyarat. Prosedur pemulihan kendiri-konsisten cuba memulihkan sebahagian isyarat itu dalam pasca-pemprosesan. Prosedur ini adalah berulang dan memerlukan sebagai input anggaran purata penghunian setiap orbital dalam keadaan dasar, yang pertama kali dikira daripada sampel mentah. Prosedur dijalankan dalam gelung, dan setiap iterasi mempunyai langkah-langkah berikut:

1. Untuk setiap rentetan bit yang melanggar simetri yang ditetapkan, balikkan bitnya dengan prosedur kebarangkalian yang direka untuk membawa rentetan bit lebih dekat kepada anggaran semasa purata penghunian orbital, untuk mendapatkan rentetan bit baharu.
2. Kumpulkan semua rentetan bit lama dan baharu yang memenuhi simetri, dan subsampel subset bersaiz tetap, dipilih terlebih dahulu.
3. Untuk setiap subset rentetan bit, unjurkan Hamiltonian ke dalam subruang yang dibentang oleh vektor asas yang sepadan (lihat bahagian sebelumnya untuk penerangan vektor asas ini), dan kira anggaran keadaan dasar bagi Hamiltonian yang diunjurkan pada komputer klasik.
4. Kemas kini anggaran purata penghunian orbital dengan anggaran keadaan dasar yang mempunyai tenaga paling rendah.

Rajah aliran kerja SQD

Aliran kerja SQD digambarkan dalam rajah berikut:

Keperluan

Sebelum memulakan tutorial ini, pastikan anda telah memasang yang berikut:

• Qiskit SDK v1.0 atau lebih baru, dengan sokongan visualisasi
• Qiskit Runtime v0.22 atau lebih baru (pip install qiskit-ibm-runtime)
• SQD Qiskit addon v0.11 atau lebih baru (pip install qiskit-addon-sqd)
• ffsim v0.0.58 atau lebih baru (pip install ffsim)

Persediaan

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

import pyscf
import pyscf.cc
import pyscf.mcscf
import ffsim
import numpy as np
import matplotlib.pyplot as plt

from qiskit import QuantumCircuit, QuantumRegister
from qiskit.transpiler.preset_passmanagers import generate_preset_pass_manager

from qiskit_ibm_runtime import QiskitRuntimeService
from qiskit_ibm_runtime import SamplerV2 as Sampler

Langkah 1: Petakan input klasik ke masalah kuantum

Dalam tutorial ini, kami akan mencari penghampiran kepada keadaan dasar molekul pada keseimbangan dalam set asas cc-pVDZ. Pertama, kami menentukan molekul dan sifat-sifatnya.

# Specify molecule properties
open_shell = False
spin_sq = 0

# Build N2 molecule
mol = pyscf.gto.Mole()
mol.build(
    atom=[["N", (0, 0, 0)], ["N", (1.0, 0, 0)]],
    basis="cc-pvdz",
    symmetry="Dooh",
)

# Define active space
n_frozen = 2
active_space = range(n_frozen, mol.nao_nr())

# Get molecular integrals
scf = pyscf.scf.RHF(mol).run()
num_orbitals = len(active_space)
n_electrons = int(sum(scf.mo_occ[active_space]))
num_elec_a = (n_electrons + mol.spin) // 2
num_elec_b = (n_electrons - mol.spin) // 2
cas = pyscf.mcscf.CASCI(scf, num_orbitals, (num_elec_a, num_elec_b))
mo = cas.sort_mo(active_space, base=0)
hcore, nuclear_repulsion_energy = cas.get_h1cas(mo)
eri = pyscf.ao2mo.restore(1, cas.get_h2cas(mo), num_orbitals)

# Store reference energy from SCI calculation performed separately
exact_energy = -109.22690201485733

converged SCF energy = -108.929838385609

Sebelum membina Circuit ansatz LUCJ, kami terlebih dahulu menjalankan pengiraan CCSD dalam sel kod berikut. Amplitud $t_1$ dan $t_2$ daripada pengiraan ini akan digunakan untuk memulakan parameter ansatz.

# Get CCSD t2 amplitudes for initializing the ansatz
ccsd = pyscf.cc.CCSD(
    scf, frozen=[i for i in range(mol.nao_nr()) if i not in active_space]
).run()
t1 = ccsd.t1
t2 = ccsd.t2

E(CCSD) = -109.2177884185543  E_corr = -0.2879500329450045

Kini, kami menggunakan ffsim untuk mencipta Circuit ansatz. Memandangkan molekul kami mempunyai keadaan Hartree-Fock cangkerang tertutup, kami menggunakan varian ansatz UCJ yang seimbang spin, UCJOpSpinBalanced. Kami menghantar pasangan interaksi yang sesuai untuk topologi Qubit kekisi segi enam berat (lihat bahagian latar belakang tentang ansatz LUCJ). Kami menetapkan optimize=True dalam kaedah from_t_amplitudes untuk mendayakan pemfaktoran berganda "dimampatkan" bagi amplitud $t_2$ (lihat The local unitary cluster Jastrow (LUCJ) ansatz dalam dokumentasi ffsim untuk butiran).

n_reps = 1
alpha_alpha_indices = [(p, p + 1) for p in range(num_orbitals - 1)]
alpha_beta_indices = [(p, p) for p in range(0, num_orbitals, 4)]

ucj_op = ffsim.UCJOpSpinBalanced.from_t_amplitudes(
    t2=t2,
    t1=t1,
    n_reps=n_reps,
    interaction_pairs=(alpha_alpha_indices, alpha_beta_indices),
    # Setting optimize=True enables the "compressed" factorization
    optimize=True,
    # Limit the number of optimization iterations to prevent the code cell from running
    # too long. Removing this line may improve results.
    options=dict(maxiter=1000),
)

nelec = (num_elec_a, num_elec_b)

# create an empty quantum circuit
qubits = QuantumRegister(2 * num_orbitals, name="q")
circuit = QuantumCircuit(qubits)

# prepare Hartree-Fock state as the reference state and append it to the quantum circuit
circuit.append(ffsim.qiskit.PrepareHartreeFockJW(num_orbitals, nelec), qubits)

# apply the UCJ operator to the reference state
circuit.append(ffsim.qiskit.UCJOpSpinBalancedJW(ucj_op), qubits)
circuit.measure_all()

Langkah 2: Optimumkan masalah untuk pelaksanaan perkakasan kuantum

Seterusnya, kami mengoptimumkan Circuit untuk perkakasan sasaran.

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

print(f"Using backend {backend.name}")

Using backend ibm_fez

Kami mengesyorkan langkah-langkah berikut untuk mengoptimumkan ansatz dan menjadikannya serasi perkakasan.

• Pilih Qubit fizikal (initial_layout) daripada perkakasan sasaran yang mematuhi corak "zig-zag" yang diterangkan sebelumnya. Menyusun Qubit dalam corak ini menghasilkan Circuit yang cekap dan serasi perkakasan dengan gate yang lebih sedikit. Di sini, kami menyertakan beberapa kod untuk mengautomasikan pemilihan corak "zig-zag" yang menggunakan heuristik pemarkahan untuk meminimumkan ralat yang berkaitan dengan susun atur yang dipilih.
• Jana pengurus laluan berstingkat menggunakan fungsi generate_preset_pass_manager daripada qiskit dengan pilihan backend dan initial_layout anda.
• Tetapkan peringkat pre_init bagi pengurus laluan berstingkat anda kepada ffsim.qiskit.PRE_INIT. ffsim.qiskit.PRE_INIT menyertakan laluan Transpiler qiskit yang mengurai gate kepada putaran orbital dan kemudian menggabungkan putaran orbital, menghasilkan lebih sedikit gate dalam Circuit akhir.
• Jalankan pengurus laluan pada Circuit anda.

Kod untuk pemilihan automatik susun atur "zig-zag"

from typing import Sequence

import rustworkx
from qiskit.providers import BackendV2
from rustworkx import NoEdgeBetweenNodes, PyGraph

IBM_TWO_Q_GATES = {"cx", "ecr", "cz"}

def create_linear_chains(num_orbitals: int) -> PyGraph:
    """In zig-zag layout, there are two linear chains (with connecting qubits between
    the chains). This function creates those two linear chains: a rustworkx PyGraph
    with two disconnected linear chains. Each chain contains `num_orbitals` number
    of nodes, that is, in the final graph there are `2 * num_orbitals` number of nodes.

    Args:
        num_orbitals (int): Number orbitals or nodes in each linear chain. They are
            also known as alpha-alpha interaction qubits.

    Returns:
        A rustworkx.PyGraph with two disconnected linear chains each with `num_orbitals`
            number of nodes.
    """
    G = rustworkx.PyGraph()

    for n in range(num_orbitals):
        G.add_node(n)

    for n in range(num_orbitals - 1):
        G.add_edge(n, n + 1, None)

    for n in range(num_orbitals, 2 * num_orbitals):
        G.add_node(n)

    for n in range(num_orbitals, 2 * num_orbitals - 1):
        G.add_edge(n, n + 1, None)

    return G

def create_lucj_zigzag_layout(
    num_orbitals: int, backend_coupling_graph: PyGraph
) -> tuple[PyGraph, int]:
    """This function creates the complete zigzag graph that 'can be mapped' to an IBM QPU with
    heavy-hex connectivity (the zigzag must be an isomorphic sub-graph to the QPU/backend
    coupling graph for it to be mapped).
    The zigzag pattern includes both linear chains (alpha-alpha interactions) and connecting
    qubits between the linear chains (alpha-beta interactions).

    Args:
        num_orbitals (int): Number of orbitals, that is, number of nodes in each alpha-alpha linear chain.
        backend_coupling_graph (PyGraph): The coupling graph of the backend on which the LUCJ ansatz
            will be mapped and run. This function takes the coupling graph as a undirected
            `rustworkx.PyGraph` where there is only one 'undirected' edge between two nodes,
            that is, qubits. Usually, the coupling graph of a IBM backend is directed (for example, Eagle devices
            such as ibm_brisbane) or may have two edges between two nodes (for example, Heron `ibm_torino`).
            A user needs to be make such graphs undirected and/or remove duplicate edges to make them
            compatible with this function.

    Returns:
        G_new (PyGraph): The graph with IBM backend compliant zigzag pattern.
        num_alpha_beta_qubits (int): Number of connecting qubits between the linear chains
            in the zigzag pattern. While we want as many connecting (alpha-beta) qubits between
            the linear (alpha-alpha) chains, we cannot accommodate all due to qubit and connectivity
            constraints of backends. This is the maximum number of connecting qubits the zigzag pattern
            can have while being backend compliant (that is, isomorphic to backend coupling graph).
    """
    isomorphic = False
    G = create_linear_chains(num_orbitals=num_orbitals)

    num_iters = num_orbitals
    while not isomorphic:
        G_new = G.copy()
        num_alpha_beta_qubits = 0
        for n in range(num_iters):
            if n % 4 == 0:
                new_node = 2 * num_orbitals + num_alpha_beta_qubits
                G_new.add_node(new_node)
                G_new.add_edge(n, new_node, None)
                G_new.add_edge(new_node, n + num_orbitals, None)
                num_alpha_beta_qubits = num_alpha_beta_qubits + 1
        isomorphic = rustworkx.is_subgraph_isomorphic(
            backend_coupling_graph, G_new
        )
        num_iters -= 1

    return G_new, num_alpha_beta_qubits

def lightweight_layout_error_scoring(
    backend: BackendV2,
    virtual_edges: Sequence[Sequence[int]],
    physical_layouts: Sequence[int],
    two_q_gate_name: str,
) -> list[list[list[int], float]]:
    """Lightweight and heuristic function to score isomorphic layouts. There can be many zigzag patterns,
    each with different set of physical qubits, that can be mapped to a backend. Some of them may
    include less noise qubits and couplings than others. This function computes a simple error score
    for each such layout. It sums up 2Q gate error for all couplings in the zigzag pattern (layout) and
    measurement of errors of physical qubits in the layout to compute the error score.

    Note:
        This lightweight scoring can be refined using concepts such as mapomatic.

    Args:
        backend (BackendV2): A backend.
        virtual_edges (Sequence[Sequence[int]]): Edges in the device compliant zigzag pattern where
            nodes are numbered from 0 to (2 * num_orbitals + num_alpha_beta_qubits).
        physical_layouts (Sequence[int]): All physical layouts of the zigzag pattern that are isomorphic
            to each other and to the larger backend coupling map.
        two_q_gate_name (str): The name of the two-qubit gate of the backend. The name is used for fetching
            two-qubit gate error from backend properties.

    Returns:
        scores (list): A list of lists where each sublist contains two items. First item is the layout, and
            second item is a float representing error score of the layout. The layouts in the `scores` are
            sorted in the ascending order of error score.
    """
    props = backend.properties()
    scores = []
    for layout in physical_layouts:
        total_2q_error = 0
        for edge in virtual_edges:
            physical_edge = (layout[edge[0]], layout[edge[1]])
            try:
                ge = props.gate_error(two_q_gate_name, physical_edge)
            except Exception:
                ge = props.gate_error(two_q_gate_name, physical_edge[::-1])
            total_2q_error += ge
        total_measurement_error = 0
        for qubit in layout:
            meas_error = props.readout_error(qubit)
            total_measurement_error += meas_error
        scores.append([layout, total_2q_error + total_measurement_error])

    return sorted(scores, key=lambda x: x[1])

def _make_backend_cmap_pygraph(backend: BackendV2) -> PyGraph:
    graph = backend.coupling_map.graph
    if not graph.is_symmetric():
        graph.make_symmetric()
    backend_coupling_graph = graph.to_undirected()

    edge_list = backend_coupling_graph.edge_list()
    removed_edge = []
    for edge in edge_list:
        if set(edge) in removed_edge:
            continue
        try:
            backend_coupling_graph.remove_edge(edge[0], edge[1])
            removed_edge.append(set(edge))
        except NoEdgeBetweenNodes:
            pass

    return backend_coupling_graph

def get_zigzag_physical_layout(
    num_orbitals: int, backend: BackendV2, score_layouts: bool = True
) -> tuple[list[int], int]:
    """The main function that generates the zigzag pattern with physical qubits that can be used
    as an `intial_layout` in a preset passmanager/transpiler.

    Args:
        num_orbitals (int): Number of orbitals.
        backend (BackendV2): A backend.
        score_layouts (bool): Optional. If `True`, it uses the `lightweight_layout_error_scoring`
            function to score the isomorphic layouts and returns the layout with less erroneous qubits.
            If `False`, returns the first isomorphic subgraph.

    Returns:
        A tuple of device compliant layout (list[int]) with zigzag pattern and an int representing
            number of alpha-beta-interactions.
    """
    backend_coupling_graph = _make_backend_cmap_pygraph(backend=backend)

    G, num_alpha_beta_qubits = create_lucj_zigzag_layout(
        num_orbitals=num_orbitals,
        backend_coupling_graph=backend_coupling_graph,
    )

    isomorphic_mappings = rustworkx.vf2_mapping(
        backend_coupling_graph, G, subgraph=True
    )
    isomorphic_mappings = list(isomorphic_mappings)

    edges = list(G.edge_list())

    layouts = []
    for mapping in isomorphic_mappings:
        initial_layout = [None] * (2 * num_orbitals + num_alpha_beta_qubits)
        for key, value in mapping.items():
            initial_layout[value] = key
        layouts.append(initial_layout)

    two_q_gate_name = IBM_TWO_Q_GATES.intersection(
        backend.configuration().basis_gates
    ).pop()

    if score_layouts:
        scores = lightweight_layout_error_scoring(
            backend=backend,
            virtual_edges=edges,
            physical_layouts=layouts,
            two_q_gate_name=two_q_gate_name,
        )

        return scores[0][0][:-num_alpha_beta_qubits], num_alpha_beta_qubits

    return layouts[0][:-num_alpha_beta_qubits], num_alpha_beta_qubits

initial_layout, _ = get_zigzag_physical_layout(num_orbitals, backend=backend)

pass_manager = generate_preset_pass_manager(
    optimization_level=3, backend=backend, initial_layout=initial_layout
)

# without PRE_INIT passes
isa_circuit = pass_manager.run(circuit)
print(f"Gate counts (w/o pre-init passes): {isa_circuit.count_ops()}")

# with PRE_INIT passes
# We will use the circuit generated by this pass manager for hardware execution
pass_manager.pre_init = ffsim.qiskit.PRE_INIT
isa_circuit = pass_manager.run(circuit)
print(f"Gate counts (w/ pre-init passes): {isa_circuit.count_ops()}")

Gate counts (w/o pre-init passes): OrderedDict({'rz': 12438, 'sx': 12169, 'cz': 3986, 'x': 573, 'measure': 52, 'barrier': 1})
Gate counts (w/ pre-init passes): OrderedDict({'sx': 7059, 'rz': 6962, 'cz': 1876, 'measure': 52, 'x': 35, 'barrier': 1})

Langkah 3: Jalankan menggunakan primitif Qiskit

Selepas mengoptimumkan Circuit untuk pelaksanaan perkakasan, kami sedia untuk menjalankannya pada perkakasan sasaran dan mengumpul sampel untuk anggaran tenaga keadaan dasar. Memandangkan kami hanya mempunyai satu Circuit, kami akan menggunakan mod pelaksanaan Job Qiskit Runtime dan melaksanakan Circuit kami.

sampler = Sampler(mode=backend)
job = sampler.run([isa_circuit], shots=100_000)

primitive_result = job.result()
pub_result = primitive_result[0]

Langkah 4: Pasca-proses dan kembalikan keputusan dalam format klasik yang dikehendaki

Kini, kami menganggarkan tenaga keadaan dasar Hamiltonian menggunakan fungsi diagonalize_fermionic_hamiltonian. Fungsi ini menjalankan prosedur pemulihan konfigurasi kendiri-konsisten untuk secara berulang memperhalusi sampel kuantum yang bermasalah bagi meningkatkan anggaran tenaga. Kami menghantar fungsi panggil balik supaya kami dapat menyimpan keputusan pertengahan untuk analisis kemudian. Lihat dokumentasi API untuk penjelasan argumen bagi diagonalize_fermionic_hamiltonian.

from functools import partial

from qiskit_addon_sqd.fermion import (
    SCIResult,
    diagonalize_fermionic_hamiltonian,
    solve_sci_batch,
)

# SQD options
energy_tol = 1e-3
occupancies_tol = 1e-3
max_iterations = 5

# Eigenstate solver options
num_batches = 3
samples_per_batch = 300
symmetrize_spin = True
carryover_threshold = 1e-4
max_cycle = 200

# Pass options to the built-in eigensolver. If you just want to use the defaults,
# you can omit this step, in which case you would not specify the sci_solver argument
# in the call to diagonalize_fermionic_hamiltonian below.
sci_solver = partial(solve_sci_batch, spin_sq=0.0, max_cycle=max_cycle)

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

result = diagonalize_fermionic_hamiltonian(
    hcore,
    eri,
    pub_result.data.meas,
    samples_per_batch=samples_per_batch,
    norb=num_orbitals,
    nelec=nelec,
    num_batches=num_batches,
    energy_tol=energy_tol,
    occupancies_tol=occupancies_tol,
    max_iterations=max_iterations,
    sci_solver=sci_solver,
    symmetrize_spin=symmetrize_spin,
    carryover_threshold=carryover_threshold,
    callback=callback,
    seed=12345,
)

Iteration 1
	Subsample 0
		Energy: -108.97781410104506
		Subspace dimension: 28561
	Subsample 1
		Energy: -108.97781410104506
		Subspace dimension: 28561
	Subsample 2
		Energy: -108.97781410104506
		Subspace dimension: 28561
Iteration 2
	Subsample 0
		Energy: -109.05150860754739
		Subspace dimension: 287296
	Subsample 1
		Energy: -109.08152283263908
		Subspace dimension: 264196
	Subsample 2
		Energy: -109.11636893067873
		Subspace dimension: 284089
Iteration 3
	Subsample 0
		Energy: -109.15878555367885
		Subspace dimension: 429025
	Subsample 1
		Energy: -109.16462831275786
		Subspace dimension: 473344
	Subsample 2
		Energy: -109.15895026995382
		Subspace dimension: 435600
Iteration 4
	Subsample 0
		Energy: -109.17784051223317
		Subspace dimension: 622521
	Subsample 1
		Energy: -109.1774651326829
		Subspace dimension: 657721
	Subsample 2
		Energy: -109.18085212360191
		Subspace dimension: 617796
Iteration 5
	Subsample 0
		Energy: -109.18636242518915
		Subspace dimension: 815409
	Subsample 1
		Energy: -109.18451014767518
		Subspace dimension: 837225
	Subsample 2
		Energy: -109.18333728638888
		Subspace dimension: 857476

Visualisasi keputusan

Plot pertama menunjukkan bahawa selepas beberapa iterasi, kami menganggarkan tenaga keadaan dasar dalam ~41 mH (ketepatan kimia biasanya diterima sebagai 1 kcal/mol $\approx$ 1.6 mH). Tenaga boleh ditingkatkan dengan membenarkan lebih banyak iterasi pemulihan konfigurasi atau meningkatkan bilangan sampel setiap batch.

Plot kedua menunjukkan purata penghunian setiap orbital ruang selepas iterasi akhir. Kami dapat melihat bahawa kedua-dua elektron spin-atas dan spin-bawah menghuni lima orbital pertama dengan kebarangkalian yang tinggi dalam penyelesaian kami.

# Data for energies plot
x1 = range(len(result_history))
min_e = [
    min(result, key=lambda res: res.energy).energy + nuclear_repulsion_energy
    for result in result_history
]
e_diff = [abs(e - exact_energy) for e in min_e]
yt1 = [1.0, 1e-1, 1e-2, 1e-3, 1e-4]

# Chemical accuracy (+/- 1 milli-Hartree)
chem_accuracy = 0.001

# 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, e_diff, label="energy error", marker="o")
axs[0].set_xticks(x1)
axs[0].set_xticklabels(x1)
axs[0].set_yticks(yt1)
axs[0].set_yticklabels(yt1)
axs[0].set_yscale("log")
axs[0].set_ylim(1e-4)
axs[0].axhline(
    y=chem_accuracy,
    color="#BF5700",
    linestyle="--",
    label="chemical accuracy",
)
axs[0].set_title("Approximated Ground State Energy Error vs SQD Iterations")
axs[0].set_xlabel("Iteration Index", fontdict={"fontsize": 12})
axs[0].set_ylabel("Energy Error (Ha)", 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})

plt.tight_layout()
plt.show()

Tinjauan tutorial

Sila ambil tinjauan ringkas ini untuk memberikan maklum balas tentang tutorial ini. Pandangan anda akan membantu kami meningkatkan kandungan dan pengalaman pengguna kami.

Pautan ke tinjauan

Note: This survey is provided by IBM Quantum and relates to the original English content. To give feedback on doQumentation's website, translations, or code execution, please open a GitHub issue.