Simulasi Kicked Ising Hamiltonian dengan dynamic circuits

Estimasi penggunaan: 7,5 menit pada prosesor Heron r3. (CATATAN: Ini hanya estimasi. Waktu aktual kamu bisa berbeda.) Dynamic circuits adalah circuit dengan classical feedforward — dengan kata lain, pengukuran di tengah circuit yang diikuti oleh operasi logika klasik yang menentukan operasi kuantum berdasarkan output klasik tersebut. Dalam tutorial ini, kita mensimulasikan model Ising kicked pada kisi heksagonal spin dan menggunakan dynamic circuits untuk merealisasikan interaksi melampaui konektivitas fisik perangkat keras.

Model Ising telah dipelajari secara ekstensif di berbagai bidang fisika. Model ini menggambarkan spin yang mengalami interaksi Ising antar situs kisi, serta tendangan dari medan magnet lokal di setiap situs. Evolusi waktu Trotterisasi spin yang dipertimbangkan dalam tutorial ini, diambil dari [1], diberikan oleh unitary berikut:

$$U(\theta)=\left(\prod_{\langle j, k\rangle} \exp \left(i \frac{\pi}{8} Z_j Z_k\right)\right)\left(\prod_j \exp \left(-i \frac{\theta}{2} X_j\right)\right)$$

Untuk menelaah dinamika spin, kita mempelajari rata-rata magnetisasi spin di setiap situs sebagai fungsi dari langkah Trotter. Oleh karena itu, kita membangun observable berikut:

$$\langle O\rangle = \frac{1}{N} \sum_i \langle Z_i \rangle$$

Untuk merealisasikan interaksi ZZ antar situs kisi, kita menyajikan solusi menggunakan fitur dynamic circuit, yang menghasilkan kedalaman dua-qubit yang jauh lebih pendek dibandingkan metode routing standar dengan gate SWAP. Di sisi lain, operasi classical feedforward dalam dynamic circuits biasanya memiliki waktu eksekusi yang lebih lama daripada gate kuantum; sehingga dynamic circuits memiliki keterbatasan dan trade-off. Kita juga menyajikan cara menambahkan urutan dynamical decoupling pada qubit yang menganggur selama operasi classical feedforward menggunakan durasi stretch.

Persyaratan

Sebelum memulai tutorial ini, pastikan kamu telah menginstal:

• Qiskit SDK v2.0 atau lebih baru dengan dukungan visualisasi
• Qiskit Runtime v0.37 atau lebih baru dengan dukungan visualisasi (pip install 'qiskit-ibm-runtime[visualization]')
• Library graf Rustworkx (pip install rustworkx)
• Qiskit Aer (pip install qiskit-aer)

Pengaturan

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

import numpy as np
from typing import List
import rustworkx as rx
import matplotlib.pyplot as plt
from rustworkx.visualization import mpl_draw
from qiskit.circuit import (
    Parameter,
    QuantumCircuit,
    QuantumRegister,
    ClassicalRegister,
)
from qiskit.transpiler import CouplingMap
from qiskit.quantum_info import SparsePauliOp
from qiskit.circuit.classical import expr
from qiskit.transpiler.preset_passmanagers import (
    generate_preset_pass_manager,
)
from qiskit.transpiler import PassManager
from qiskit.circuit.library import RZGate, XGate
from qiskit.transpiler.passes import (
    ALAPScheduleAnalysis,
    PadDynamicalDecoupling,
)

from qiskit.transpiler.basepasses import TransformationPass
from qiskit.circuit.measure import Measure
from qiskit.transpiler.passes.utils.remove_final_measurements import (
    calc_final_ops,
)
from qiskit.circuit import Instruction

from qiskit.visualization import plot_circuit_layout
from qiskit.circuit.tools import pi_check

from qiskit_aer import AerSimulator
from qiskit_aer.primitives import SamplerV2 as Aer_Sampler

from qiskit_ibm_runtime import (
    QiskitRuntimeService,
    Batch,
    SamplerV2 as Sampler,
)
from qiskit_ibm_runtime.exceptions import QiskitBackendNotFoundError
from qiskit_ibm_runtime.visualization import (
    draw_circuit_schedule_timing,
)

Langkah 1: Petakan input klasik ke circuit kuantum

Kita mulai dengan mendefinisikan kisi yang akan disimulasikan. Kita memilih untuk bekerja dengan kisi honeycomb (disebut juga heksagonal), yaitu graf planar dengan simpul berderajat 3. Di sini, kita menentukan ukuran kisi, parameter circuit yang relevan dalam dinamika Trotterisasi. Kita mensimulasikan evolusi waktu Trotterisasi di bawah model Ising dengan tiga nilai $\theta$ berbeda dari medan magnet lokal.

hex_rows = 3  # specify lattice size
hex_cols = 5
depths = range(9)  # specify Trotter steps
zz_angle = np.pi / 8  # parameter for ZZ interaction
max_angle = np.pi / 2  # max theta angle
points = 3  # number of theta parameters

θ = Parameter("θ")
params = np.linspace(0, max_angle, points)

def make_hex_lattice(hex_rows=1, hex_cols=1):
    """Define hexagon lattice."""
    hex_cmap = CouplingMap.from_hexagonal_lattice(
        hex_rows, hex_cols, bidirectional=False
    )
    data = list(hex_cmap.physical_qubits)
    graph = hex_cmap.graph.to_undirected(multigraph=False)
    edge_colors = rx.graph_misra_gries_edge_color(graph)
    layer_edges = {color: [] for color in edge_colors.values()}
    for edge_index, color in edge_colors.items():
        layer_edges[color].append(graph.edge_list()[edge_index])
    return data, layer_edges, hex_cmap, graph

Mari mulai dengan contoh uji kecil:

hex_rows_test = 1
hex_cols_test = 2

data_test, layer_edges_test, hex_cmap_test, graph_test = make_hex_lattice(
    hex_rows=hex_rows_test, hex_cols=hex_cols_test
)

# display a small example for illustration
node_colors_test = ["lightblue"] * len(graph_test.node_indices())
pos = rx.graph_spring_layout(
    graph_test,
    k=5 / np.sqrt(len(graph_test.nodes())),
    repulsive_exponent=1,
    num_iter=150,
)
mpl_draw(graph_test, node_color=node_colors_test, pos=pos)

Kita akan menggunakan contoh kecil ini untuk ilustrasi dan simulasi. Di bawah ini kita juga membangun contoh besar untuk menunjukkan bahwa alur kerja ini bisa diperluas ke ukuran yang lebih besar.

data, layer_edges, hex_cmap, graph = make_hex_lattice(
    hex_rows=hex_rows, hex_cols=hex_cols
)
num_qubits = len(data)
print(f"num_qubits = {num_qubits}")

# display the honeycomb lattice to simulate
node_colors = ["lightblue"] * len(graph.node_indices())
pos = rx.graph_spring_layout(
    graph,
    k=5 / np.sqrt(num_qubits),
    repulsive_exponent=1,
    num_iter=150,
)
mpl_draw(graph, node_color=node_colors, pos=pos)
plt.show()

num_qubits = 46

Membangun unitary circuits

Dengan ukuran masalah dan parameter yang sudah ditentukan, kita siap membangun circuit berparameter yang mensimulasikan evolusi waktu Trotterisasi $U(\theta)$ dengan berbagai langkah Trotter, yang ditentukan oleh argumen depth. Circuit yang kita bangun memiliki lapisan bergantian antara gate Rx($\theta$) dan gate Rzz. Gate Rzz merealisasikan interaksi ZZ antar spin yang terhubung, yang akan ditempatkan di antara setiap situs kisi yang ditentukan oleh argumen layer_edges.

def gen_hex_unitary(
    num_qubits=6,
    zz_angle=np.pi / 8,
    layer_edges=[
        [(0, 1), (2, 3), (4, 5)],
        [(1, 2), (3, 4), (5, 0)],
    ],
    θ=Parameter("θ"),
    depth=1,
    measure=False,
    final_rot=True,
):
    """Build unitary circuit."""
    circuit = QuantumCircuit(num_qubits)
    # Build trotter layers
    for _ in range(depth):
        for i in range(num_qubits):
            circuit.rx(θ, i)
        circuit.barrier()
        for coloring in layer_edges.keys():
            for e in layer_edges[coloring]:
                circuit.rzz(zz_angle, e[0], e[1])
        circuit.barrier()
    # Optional final rotation, set True to be consistent with Ref. [1]
    if final_rot:
        for i in range(num_qubits):
            circuit.rx(θ, i)
    if measure:
        circuit.measure_all()

    return circuit

Visualisasikan circuit uji kecil:

circ_unitary_test = gen_hex_unitary(
    num_qubits=len(data_test),
    layer_edges=layer_edges_test,
    θ=Parameter("θ"),
    depth=1,
    measure=True,
)
circ_unitary_test.draw(output="mpl", fold=-1)

Demikian juga, bangun unitary circuits contoh besar pada berbagai langkah Trotter beserta observable untuk memperkirakan nilai ekspektasi.

circuits_unitary = []
for depth in depths:
    circ = gen_hex_unitary(
        num_qubits=num_qubits,
        layer_edges=layer_edges,
        θ=Parameter("θ"),
        depth=depth,
        measure=True,
    )
    circuits_unitary.append(circ)

observables_unitary = SparsePauliOp.from_sparse_list(
    [("Z", [i], 1 / num_qubits) for i in range(num_qubits)],
    num_qubits=num_qubits,
)

Membangun implementasi dynamic circuit

Bagian ini mendemonstrasikan implementasi dynamic circuit utama untuk mensimulasikan evolusi waktu Trotterisasi yang sama. Perlu diperhatikan bahwa kisi honeycomb yang ingin kita simulasikan tidak cocok dengan kisi heavy dari qubit perangkat keras. Salah satu cara langsung untuk memetakan circuit ke perangkat keras adalah dengan memperkenalkan serangkaian operasi SWAP untuk membawa qubit yang berinteraksi agar berdekatan satu sama lain, guna merealisasikan interaksi ZZ. Di sini kita menyoroti pendekatan alternatif menggunakan dynamic circuits sebagai solusi, yang mengilustrasikan bahwa kita bisa menggunakan kombinasi komputasi kuantum dan komputasi klasik real-time dalam sebuah circuit di Qiskit untuk merealisasikan interaksi melampaui nearest-neighbor.

Dalam implementasi dynamic circuit, interaksi ZZ secara efektif diimplementasikan menggunakan qubit ancilla, pengukuran mid-circuit, dan feedforward. Untuk memahami ini, perlu diingat bahwa rotasi ZZ menerapkan faktor fase $e^{i\theta}$ pada state berdasarkan paritas-nya. Untuk dua qubit, state basis komputasi adalah $|00\rangle$, $|01\rangle$, $|10\rangle$, dan $|11\rangle$. Gate rotasi ZZ menerapkan faktor fase pada state $|01\rangle$ dan $|10\rangle$ yang paritasnya (jumlah satu dalam state) ganjil dan membiarkan state paritas genap tidak berubah. Berikut ini menjelaskan bagaimana kita bisa mengimplementasikan interaksi ZZ secara efektif pada dua qubit menggunakan dynamic circuits.

1. Hitung paritas ke dalam qubit ancilla: alih-alih menerapkan ZZ langsung ke dua qubit, kita memperkenalkan qubit ketiga, yaitu qubit ancilla, untuk menyimpan informasi paritas dari dua qubit data. Kita mengaitkan ancilla dengan setiap qubit data menggunakan gate CX dari qubit data ke qubit ancilla.

2. Terapkan rotasi Z satu-qubit ke qubit ancilla: ini karena ancilla memiliki informasi paritas dari dua qubit data, yang secara efektif mengimplementasikan rotasi ZZ pada qubit data.

3. Ukur qubit ancilla dalam basis X: ini adalah langkah kunci yang meruntuhkan state qubit ancilla, dan hasil pengukuran memberi tahu kita apa yang terjadi:

   • Ukur 0: ketika hasil 0 diamati, kita sebenarnya telah menerapkan rotasi $ZZ(\theta)$ dengan benar pada qubit data kita.

   • Ukur 1: ketika hasil 1 diamati, kita telah menerapkan $ZZ(\theta + \pi)$ sebagai gantinya.

4. Terapkan gate koreksi ketika mengukur 1: jika kita mengukur 1, kita menerapkan gate Z ke qubit data untuk "memperbaiki" fase $\pi$ ekstra.

Circuit yang dihasilkan adalah sebagai berikut:

Ketika kita mengadopsi pendekatan ini untuk mensimulasikan kisi honeycomb, circuit yang dihasilkan tertanam sempurna ke dalam perangkat keras dengan kisi heavy-hex: semua qubit data berada di situs berderajat-3 dari kisi, yang membentuk kisi heksagonal. Setiap pasang qubit data berbagi satu qubit ancilla yang berada di situs berderajat-2. Di bawah ini, kita membangun kisi qubit untuk implementasi dynamic circuit, dengan memperkenalkan qubit ancilla (ditunjukkan dalam lingkaran ungu yang lebih gelap).

def make_lattice(hex_rows=1, hex_cols=1):
    """Define heavy-hex lattice and corresponding lists of data and ancilla nodes."""
    hex_cmap = CouplingMap.from_hexagonal_lattice(
        hex_rows, hex_cols, bidirectional=False
    )
    data = list(hex_cmap.physical_qubits)

    heavyhex_cmap = CouplingMap()
    for d in data:
        heavyhex_cmap.add_physical_qubit(d)

    # make coupling map
    a = len(data)
    for edge in hex_cmap.get_edges():
        heavyhex_cmap.add_physical_qubit(a)
        heavyhex_cmap.add_edge(edge[0], a)
        heavyhex_cmap.add_edge(edge[1], a)
        a += 1
    ancilla = list(range(len(data), a))
    qubits = data + ancilla

    # color edges
    graph = heavyhex_cmap.graph.to_undirected(multigraph=False)
    edge_colors = rx.graph_misra_gries_edge_color(graph)
    layer_edges = {color: [] for color in edge_colors.values()}
    for edge_index, color in edge_colors.items():
        layer_edges[color].append(graph.edge_list()[edge_index])

    # construct observable
    obs_hex = SparsePauliOp.from_sparse_list(
        [("Z", [i], 1 / len(data)) for i in data],
        num_qubits=len(qubits),
    )

    return (data, qubits, ancilla, layer_edges, heavyhex_cmap, graph, obs_hex)

Visualisasikan kisi heavy-hex untuk qubit data dan qubit ancilla pada skala kecil:

(data, qubits, ancilla, layer_edges, heavyhex_cmap, graph, obs_hex) = (
    make_lattice(hex_rows=hex_rows, hex_cols=hex_cols)
)

print(f"number of data qubits = {len(data)}")
print(f"number of ancilla qubits = {len(ancilla)}")

node_colors = []
for node in graph.node_indices():
    if node in ancilla:
        node_colors.append("purple")
    else:
        node_colors.append("lightblue")

pos = rx.graph_spring_layout(
    graph,
    k=1 / np.sqrt(len(qubits)),
    repulsive_exponent=2,
    num_iter=200,
)

# Visualize the graph, blue circles are data qubits and purple circles are ancillas
mpl_draw(graph, node_color=node_colors, pos=pos)
plt.show()

number of data qubits = 46
number of ancilla qubits = 60

Di bawah ini, kita membangun dynamic circuit untuk evolusi waktu Trotterisasi. Gate RZZ digantikan dengan implementasi dynamic circuit menggunakan langkah-langkah yang dijelaskan di atas.

def gen_hex_dynamic(
    depth=1,
    zz_angle=np.pi / 8,
    θ=Parameter("θ"),
    hex_rows=1,
    hex_cols=1,
    measure=False,
    add_dd=True,
):
    """Build dynamic circuits."""
    (data, qubits, ancilla, layer_edges, heavyhex_cmap, graph, obs_hex) = (
        make_lattice(hex_rows=hex_rows, hex_cols=hex_cols)
    )
    # Initialize circuit
    qr = QuantumRegister(len(qubits), "qr")
    cr = ClassicalRegister(len(ancilla), "cr")
    circuit = QuantumCircuit(qr, cr)

    for k in range(depth):
        # Single-qubit Rx layer
        for d in data:
            circuit.rx(θ, d)
        circuit.barrier()

        # CX gates from data qubits to ancilla qubits
        for same_color_edges in layer_edges.values():
            for e in same_color_edges:
                circuit.cx(e[0], e[1])
        circuit.barrier()

        # Apply Rz rotation on ancilla qubits and rotate into X basis
        for a in ancilla:
            circuit.rz(zz_angle, a)
            circuit.h(a)
        # Add barrier to align terminal measurement
        circuit.barrier()

        # Measure ancilla qubits
        for i, a in enumerate(ancilla):
            circuit.measure(a, i)
        d2ros = {}
        a2ro = {}
        # Retrieve ancilla measurement outcomes
        for a in ancilla:
            a2ro[a] = cr[ancilla.index(a)]

        # For each data qubit, retrieve measurement outcomes of neighboring ancilla qubits
        for d in data:
            ros = [a2ro[a] for a in heavyhex_cmap.neighbors(d)]
            d2ros[d] = ros

        # Build classical feedforward operations (optionally add DD on idling data qubits)
        for d in data:
            if add_dd:
                circuit = add_stretch_dd(circuit, d, f"data_{d}_depth_{k}")

            # # XOR the neighboring readouts of the data qubit; if True, apply Z to it
            ros = d2ros[d]
            parity = ros[0]
            for ro in ros[1:]:
                parity = expr.bit_xor(parity, ro)
            with circuit.if_test(expr.equal(parity, True)):
                circuit.z(d)

        # Reset the ancilla if its readout is 1
        for a in ancilla:
            with circuit.if_test(expr.equal(a2ro[a], True)):
                circuit.x(a)
        circuit.barrier()

    # Final single-qubit Rx layer to match the unitary circuits
    for d in data:
        circuit.rx(θ, d)

    if measure:
        circuit.measure_all()
    return circuit, obs_hex

def add_stretch_dd(qc, q, name):
    """Add XpXm DD sequence."""
    s = qc.add_stretch(name)
    qc.delay(s, q)
    qc.x(q)
    qc.delay(s, q)
    qc.delay(s, q)
    qc.rz(np.pi, q)
    qc.x(q)
    qc.rz(-np.pi, q)
    qc.delay(s, q)
    return qc

Dynamical decoupling (DD) dan dukungan untuk durasi stretch

Salah satu kelemahan penggunaan implementasi dynamic circuit untuk merealisasikan interaksi ZZ adalah bahwa pengukuran mid-circuit dan operasi classical feedforward biasanya membutuhkan waktu eksekusi yang lebih lama daripada gate kuantum. Untuk menekan dekoherensi qubit selama waktu menganggur saat operasi klasik berlangsung, kita menambahkan urutan dynamical decoupling (DD) setelah operasi pengukuran pada qubit ancilla, dan sebelum operasi Z kondisional pada qubit data, sebelum pernyataan if_test.

Urutan DD ditambahkan oleh fungsi add_stretch_dd(), yang menggunakan durasi stretch untuk menentukan interval waktu antara gate DD. Durasi stretch adalah cara untuk menentukan durasi waktu yang bisa diregangkan untuk operasi delay sehingga durasi penundaan bisa bertambah untuk mengisi waktu idle qubit. Variabel durasi yang ditentukan oleh stretch diselesaikan pada waktu kompilasi menjadi durasi yang diinginkan yang memenuhi batasan tertentu. Ini sangat berguna ketika timing urutan DD sangat penting untuk mencapai performa penekanan error yang baik. Untuk detail lebih lanjut tentang tipe stretch, lihat dokumentasi OpenQASM. Saat ini, dukungan untuk tipe stretch di Qiskit Runtime bersifat eksperimental. Untuk detail tentang batasan penggunaannya, silakan lihat bagian keterbatasan dari dokumentasi stretch.

Menggunakan fungsi-fungsi yang didefinisikan di atas, kita membangun circuit evolusi waktu Trotterisasi, dengan dan tanpa DD, beserta observable yang sesuai. Kita mulai dengan memvisualisasikan dynamic circuit dari contoh kecil:

hex_rows_test = 1
hex_cols_test = 1

(
    data_test,
    qubits_test,
    ancilla_test,
    layer_edges_test,
    heavyhex_cmap_test,
    graph_test,
    obs_hex_test,
) = make_lattice(hex_rows=hex_rows_test, hex_cols=hex_cols_test)

node_colors = []
for node in graph_test.node_indices():
    if node in ancilla_test:
        node_colors.append("purple")
    else:
        node_colors.append("lightblue")
pos = rx.graph_spring_layout(
    graph_test,
    k=5 / np.sqrt(len(qubits_test)),
    repulsive_exponent=2,
    num_iter=150,
)

# display a small example for illustration
node_colors_test = ["lightblue"] * len(graph_test.node_indices())
mpl_draw(graph_test, node_color=node_colors, pos=pos)

circuit_dynamic_test, obs_dynamic_test = gen_hex_dynamic(
    depth=1,
    θ=Parameter("θ"),
    hex_rows=hex_rows_test,
    hex_cols=hex_cols_test,
    measure=False,
    add_dd=False,
)
circuit_dynamic_test.draw("mpl", fold=-1)

circuit_dynamic_dd_test, _ = gen_hex_dynamic(
    depth=1,
    θ=Parameter("θ"),
    hex_rows=hex_rows_test,
    hex_cols=hex_cols_test,
    measure=False,
    add_dd=True,
)
circuit_dynamic_dd_test.draw("mpl", fold=-1)

Demikian juga, bangun dynamic circuits untuk contoh besar:

circuits_dynamic = []
circuits_dynamic_dd = []
observables_dynamic = []
for depth in depths:
    circuit, obs = gen_hex_dynamic(
        depth=depth,
        θ=Parameter("θ"),
        hex_rows=hex_rows,
        hex_cols=hex_cols,
        measure=True,
        add_dd=False,
    )
    circuits_dynamic.append(circuit)

    circuit_dd, _ = gen_hex_dynamic(
        depth=depth,
        θ=Parameter("θ"),
        hex_rows=hex_rows,
        hex_cols=hex_cols,
        measure=True,
        add_dd=True,
    )
    circuits_dynamic_dd.append(circuit_dd)
    observables_dynamic.append(obs)

Langkah 2: Optimalkan masalah untuk eksekusi perangkat keras

Kita sekarang siap untuk mentranspile circuit ke perangkat keras. Kita akan mentranspile implementasi standar unitary maupun implementasi dynamic circuit ke perangkat keras.

Untuk mentranspile ke perangkat keras, kita pertama-tama membuat instans Backend. Jika tersedia, kita akan memilih Backend yang mendukung instruksi MidCircuitMeasure (measure_2).

service = QiskitRuntimeService()
try:
    backend = service.least_busy(
        operational=True,
        simulator=False,
        use_fractional_gates=True,
        filters=lambda b: "measure_2" in b.supported_instructions,
    )
except QiskitBackendNotFoundError:
    backend = service.least_busy(
        operational=True,
        simulator=False,
        use_fractional_gates=True,
    )

Transpilasi untuk dynamic circuits

Pertama, kita mentranspile dynamic circuits, dengan dan tanpa menambahkan urutan DD. Untuk memastikan kita menggunakan set qubit fisik yang sama di semua circuit demi hasil yang lebih konsisten, kita pertama-tama mentranspile circuit satu kali, lalu menggunakan layout-nya untuk semua circuit berikutnya, yang ditentukan oleh initial_layout dalam pass manager. Kemudian kita membangun primitive unified blocs (PUB) sebagai input primitif Sampler.

pm_temp = generate_preset_pass_manager(
    optimization_level=3,
    backend=backend,
)
isa_temp = pm_temp.run(circuits_dynamic[-1])
dynamic_layout = isa_temp.layout.initial_index_layout(filter_ancillas=True)

pm = generate_preset_pass_manager(
    optimization_level=3, backend=backend, initial_layout=dynamic_layout
)

dynamic_isa_circuits = [pm.run(circ) for circ in circuits_dynamic]
dynamic_pubs = [(circ, params) for circ in dynamic_isa_circuits]

dynamic_isa_circuits_dd = [pm.run(circ) for circ in circuits_dynamic_dd]
dynamic_pubs_dd = [(circ, params) for circ in dynamic_isa_circuits_dd]

Kita bisa memvisualisasikan layout qubit dari circuit yang telah ditranspile di bawah ini. Lingkaran hitam menunjukkan qubit data dan qubit ancilla yang digunakan dalam implementasi dynamic circuit.

def _heron_coords_r2():
    cord_map = np.array(
        [
            [
                0,
                1,
                2,
                3,
                4,
                5,
                6,
                7,
                8,
                9,
                10,
                11,
                12,
                13,
                14,
                15,
                3,
                7,
                11,
                15,
                0,
                1,
                2,
                3,
                4,
                5,
                6,
                7,
                8,
                9,
                10,
                11,
                12,
                13,
                14,
                15,
                1,
                5,
                9,
                13,
                0,
                1,
                2,
                3,
                4,
                5,
                6,
                7,
                8,
                9,
                10,
                11,
                12,
                13,
                14,
                15,
                3,
                7,
                11,
                15,
                0,
                1,
                2,
                3,
                4,
                5,
                6,
                7,
                8,
                9,
                10,
                11,
                12,
                13,
                14,
                15,
                1,
                5,
                9,
                13,
                0,
                1,
                2,
                3,
                4,
                5,
                6,
                7,
                8,
                9,
                10,
                11,
                12,
                13,
                14,
                15,
                3,
                7,
                11,
                15,
                0,
                1,
                2,
                3,
                4,
                5,
                6,
                7,
                8,
                9,
                10,
                11,
                12,
                13,
                14,
                15,
                1,
                5,
                9,
                13,
                0,
                1,
                2,
                3,
                4,
                5,
                6,
                7,
                8,
                9,
                10,
                11,
                12,
                13,
                14,
                15,
                3,
                7,
                11,
                15,
                0,
                1,
                2,
                3,
                4,
                5,
                6,
                7,
                8,
                9,
                10,
                11,
                12,
                13,
                14,
                15,
            ],
            -1
            * np.array([j for i in range(15) for j in [i] * [16, 4][i % 2]]),
        ],
        dtype=int,
    )

    hcords = []
    ycords = cord_map[0]
    xcords = cord_map[1]
    for i in range(156):
        hcords.append([xcords[i] + 1, np.abs(ycords[i]) + 1])

    return hcords

plot_circuit_layout(
    dynamic_isa_circuits_dd[8],
    backend,
    qubit_coordinates=_heron_coords_r2(),
    view="virtual",
)

catatan

Jika kamu mendapatkan error tentang neato tidak ditemukan dari plot_circuit_layout(), pastikan kamu telah menginstal paket graphviz dan tersedia di PATH kamu. Jika diinstal di lokasi non-default (misalnya, menggunakan homebrew di MacOS), kamu mungkin perlu memperbarui variabel environment PATH. Ini bisa dilakukan di dalam notebook ini menggunakan perintah berikut:

import os
os.environ['PATH'] = f"path/to/neato{os.pathsep}{os.environ['PATH']}"

dynamic_isa_circuits[1].draw(fold=-1, output="mpl", idle_wires=False)

dynamic_isa_circuits_dd[1].draw(fold=-1, output="mpl", idle_wires=False)

Transpilasi menggunakan MidCircuitMeasure

MidCircuitMeasure adalah tambahan untuk operasi pengukuran yang tersedia, dikalibrasi khusus untuk melakukan pengukuran mid-circuit. Instruksi MidCircuitMeasure memetakan ke instruksi measure_2 yang didukung oleh Backend. Perlu diperhatikan bahwa measure_2 tidak didukung pada semua Backend. Kamu bisa menggunakan service.backends(filters=lambda b: "measure_2" in b.supported_instructions) untuk menemukan Backend yang mendukungnya. Di sini, kita menunjukkan cara mentranspile circuit sehingga pengukuran mid-circuit yang didefinisikan dalam circuit dieksekusi menggunakan operasi MidCircuitMeasure, jika Backend mendukungnya.

Di bawah ini, kita mencetak durasi untuk instruksi measure_2 dan instruksi measure standar.

print(
    f'Mid-circuit measurement `measure_2` duration: {backend.instruction_durations.get('measure_2',0) * backend.dt * 1e9/1e3} μs'
)
print(
    f'Terminal measurement `measure` duration: {backend.instruction_durations.get('measure',0) * backend.dt *1e9/1e3} μs'
)

Mid-circuit measurement `measure_2` duration:  1.624 μs
Terminal measurement `measure` duration: 2.2 μs

"""Pass that replaces terminal measures in the middle of the circuit with
MidCircuitMeasure instructions."""

class ConvertToMidCircuitMeasure(TransformationPass):
    """This pass replaces terminal measures in the middle of the circuit with
    MidCircuitMeasure instructions.
    """

    def __init__(self, target):
        super().__init__()
        self.target = target

    def run(self, dag):
        """Run the pass on a dag."""
        mid_circ_measure = None
        for inst in self.target.instructions:
            if isinstance(inst[0], Instruction) and inst[0].name.startswith(
                "measure_"
            ):
                mid_circ_measure = inst[0]
                break
        if not mid_circ_measure:
            return dag

        final_measure_nodes = calc_final_ops(dag, {"measure"})
        for node in dag.op_nodes(Measure):
            if node not in final_measure_nodes:
                dag.substitute_node(node, mid_circ_measure, inplace=True)

        return dag

pm = PassManager(ConvertToMidCircuitMeasure(backend.target))

dynamic_isa_circuits_meas2 = [pm.run(circ) for circ in dynamic_isa_circuits]
dynamic_pubs_meas2 = [(circ, params) for circ in dynamic_isa_circuits_meas2]

dynamic_isa_circuits_dd_meas2 = [
    pm.run(circ) for circ in dynamic_isa_circuits_dd
]
dynamic_pubs_dd_meas2 = [
    (circ, params) for circ in dynamic_isa_circuits_dd_meas2
]

Transpilasi untuk circuit unitary

Untuk membuat perbandingan yang adil antara dynamic circuits dan padanan unitarynya, kita menggunakan set qubit fisik yang sama yang digunakan pada dynamic circuits untuk data qubit sebagai layout saat mentranspilasi unitary circuits.

init_layout = [
    dynamic_layout[ind] for ind in range(circuits_unitary[0].num_qubits)
]

pm = generate_preset_pass_manager(
    target=backend.target,
    initial_layout=init_layout,
    optimization_level=3,
)

def transpile_minimize(circ: QuantumCircuit, pm: PassManager, iterations=10):
    """Transpile circuits for specified number of iterations and return the one with smallest two-qubit gate depth"""
    circs = [pm.run(circ) for i in range(iterations)]
    circs_sorted = sorted(
        circs,
        key=lambda x: x.depth(lambda x: x.operation.num_qubits == 2),
    )
    return circs_sorted[0]

unitary_isa_circuits = []
for circ in circuits_unitary:
    circ_t = transpile_minimize(circ, pm, iterations=100)
    unitary_isa_circuits.append(circ_t)

unitary_pubs = [(circ, params) for circ in unitary_isa_circuits]

Kita memvisualisasikan layout qubit dari unitary circuits yang sudah ditranspilasi. Lingkaran hitam menunjukkan qubit fisik yang digunakan untuk mentranspilasi unitary circuits, dan indeksnya sesuai dengan indeks virtual qubit. Dengan membandingkan ini dengan layout yang diplot untuk dynamic circuits, kita bisa memastikan bahwa unitary circuits menggunakan set qubit fisik yang sama dengan data qubit pada dynamic circuits.

plot_circuit_layout(
    unitary_isa_circuits[-1],
    backend,
    qubit_coordinates=_heron_coords_r2(),
    view="virtual",
)

Sekarang kita menambahkan urutan DD ke circuit yang sudah ditranspilasi dan membuat PUB yang sesuai untuk pengiriman job.

pm_dd = PassManager(
    [
        ALAPScheduleAnalysis(target=backend.target),
        PadDynamicalDecoupling(
            dd_sequence=[
                XGate(),
                RZGate(np.pi),
                XGate(),
                RZGate(-np.pi),
            ],
            spacing=[1 / 4, 1 / 2, 0, 0, 1 / 4],
            target=backend.target,
        ),
    ]
)

unitary_isa_circuits_dd = pm_dd.run(unitary_isa_circuits)
unitary_pubs_dd = [(circ, params) for circ in unitary_isa_circuits_dd]

Bandingkan kedalaman gate dua-qubit dari unitary dan dynamic circuits

# compare circuit depth of unitary and dynamic circuit implementations
unitary_depth = [
    unitary_isa_circuits[i].depth(lambda x: x.operation.num_qubits == 2)
    for i in range(len(unitary_isa_circuits))
]

dynamic_depth = [
    dynamic_isa_circuits[i].depth(lambda x: x.operation.num_qubits == 2)
    for i in range(len(dynamic_isa_circuits))
]

plt.plot(
    list(range(len(unitary_depth))),
    unitary_depth,
    label="unitary circuits",
    color="#be95ff",
)
plt.plot(
    list(range(len(dynamic_depth))),
    dynamic_depth,
    label="dynamic circuits",
    color="#ff7eb6",
)
plt.xlabel("Trotter steps")
plt.ylabel("Two-qubit depth")
plt.legend()

<matplotlib.legend.Legend at 0x374225760>

Keunggulan utama dari circuit berbasis pengukuran adalah bahwa, saat mengimplementasikan beberapa interaksi ZZ, layer CX bisa diparalelkan dan pengukuran bisa terjadi secara bersamaan. Ini karena semua interaksi ZZ bersifat komutatif, sehingga komputasi bisa dilakukan dengan kedalaman pengukuran 1. Setelah mentranspilasi circuit, kita mengamati bahwa pendekatan dynamic circuit menghasilkan kedalaman dua-qubit yang jauh lebih pendek dibandingkan pendekatan unitary standar, dengan catatan bahwa pengukuran mid-circuit tambahan dan classical feedforward sendiri membutuhkan waktu dan memperkenalkan sumber error mereka sendiri.

Langkah 3: Eksekusi menggunakan Qiskit primitives

Mode pengujian lokal

Sebelum mengirim job ke hardware, kita bisa menjalankan simulasi kecil dari dynamic circuit menggunakan mode pengujian lokal.

aer_sim = AerSimulator()
pm = generate_preset_pass_manager(backend=aer_sim, optimization_level=1)
circuit_dynamic_test.measure_all()
isa_qc = pm.run(circuit_dynamic_test)
with Batch(backend=aer_sim) as batch:
    sampler = Sampler(mode=batch)
    result = sampler.run([(isa_qc, params)]).result()

print(
    "Simulated average magnetization at trotter step = 1 at three theta values"
)
result[0].data["meas"].expectation_values(obs_dynamic_test[0])

Simulated average magnetization at trotter step = 1 at three theta values

array([ 0.16666667,  0.01855469, -0.13476562])

Simulasi MPS

Untuk circuit yang besar, kita bisa menggunakan simulator matrix_product_state (MPS), yang memberikan hasil perkiraan dari nilai ekspektasi sesuai dengan dimensi bond yang dipilih. Nanti kita menggunakan hasil simulasi MPS sebagai baseline untuk membandingkan hasil dari hardware.

# The MPS simulation below took approximately 7 minutes to run on a laptop with Apple M1 chip

mps_backend = AerSimulator(
    method="matrix_product_state",
    matrix_product_state_truncation_threshold=1e-5,
    matrix_product_state_max_bond_dimension=100,
)
mps_sampler = Aer_Sampler.from_backend(mps_backend)

shots = 4096

data_sim = []
for j in range(points):
    circ_list = [
        circ.assign_parameters([params[j]]) for circ in circuits_unitary
    ]

    mps_job = mps_sampler.run(circ_list, shots=shots)
    result = mps_job.result()

    point_data = [
        result[d].data["meas"].expectation_values(observables_unitary)
        for d in depths
    ]

    data_sim.append(point_data)  # data at one theta value

data_sim = np.array(data_sim)

Dengan circuit dan observable yang sudah disiapkan, sekarang kita mengeksekusinya di hardware menggunakan Sampler primitive.

Di sini kita mengirimkan tiga job untuk unitary_pubs, dynamic_pubs, dan dynamic_pubs_dd. Masing-masing merupakan daftar circuit berparameter yang sesuai dengan sembilan Trotter step berbeda dengan tiga parameter $\theta$ yang berbeda.

shots = 10000

with Batch(backend=backend) as batch:
    sampler = Sampler(mode=batch)

    sampler.options.experimental = {
        "execution": {
            "scheduler_timing": True
        },  # set to True to retrieve circuit timing info
    }

    job_unitary = sampler.run(unitary_pubs, shots=shots)
    print(f"unitary: {job_unitary.job_id()}")

    job_unitary_dd = sampler.run(unitary_pubs_dd, shots=shots)
    print(f"unitary_dd: {job_unitary_dd.job_id()}")

    job_dynamic = sampler.run(dynamic_pubs, shots=shots)
    print(f"dynamic: {job_dynamic.job_id()}")

    job_dynamic_dd = sampler.run(dynamic_pubs_dd, shots=shots)
    print(f"dynamic_dd: {job_dynamic_dd.job_id()}")

    job_dynamic_meas2 = sampler.run(dynamic_pubs_meas2, shots=shots)
    print(f"dynamic_meas2: {job_dynamic_meas2.job_id()}")

    job_dynamic_dd_meas2 = sampler.run(dynamic_pubs_dd_meas2, shots=shots)
    print(f"dynamic_dd_meas2: {job_dynamic_dd_meas2.job_id()}")

unitary: d5dtt0ldq8ts73fvbhj0
unitary: d5dtt11smlfc739onuag
dynamic: d5dtt1hsmlfc739onuc0
dynamic_dd: d5dtt25jngic73avdne0
dynamic_meas2: d5dtt2ldq8ts73fvbhm0
dynamic_dd_meas2: d5dtt2tjngic73avdnf0

Langkah 4: Pasca-proses dan kembalikan hasil dalam format klasik yang diinginkan

Setelah job selesai, kita bisa mengambil durasi circuit dari metadata hasil job dan memvisualisasikan informasi jadwal circuit. Untuk membaca lebih lanjut tentang memvisualisasikan informasi penjadwalan circuit, lihat halaman ini.

# Circuit durations is reported in the unit of `dt` which can be retrieved from `Backend` object
unitary_durations = [
    job_unitary.result()[i].metadata["compilation"]["scheduler_timing"][
        "circuit_duration"
    ]
    for i in depths
]

dynamic_durations = [
    job_dynamic.result()[i].metadata["compilation"]["scheduler_timing"][
        "circuit_duration"
    ]
    for i in depths
]

dynamic_durations_meas2 = [
    job_dynamic_meas2.result()[i].metadata["compilation"]["scheduler_timing"][
        "circuit_duration"
    ]
    for i in depths
]

result_dd = job_dynamic_dd.result()[1]
circuit_schedule_dd = result_dd.metadata["compilation"]["scheduler_timing"][
    "timing"
]

# to visualize the circuit schedule, one can show the figure below
fig_dd = draw_circuit_schedule_timing(
    circuit_schedule=circuit_schedule_dd,
    included_channels=None,
    filter_readout_channels=False,
    filter_barriers=False,
    width=1000,
)

# Save to a file since the figure is large
fig_dd.write_html("scheduler_timing_dd.html")

Kita memplot durasi circuit untuk unitary circuits dan dynamic circuits. Dari plot di bawah, kita bisa melihat bahwa meskipun membutuhkan waktu untuk pengukuran mid-circuit dan operasi klasik, implementasi dynamic circuit dengan measure_2 menghasilkan durasi circuit yang sebanding dengan implementasi unitary.

# visualize circuit durations

def convert_dt_to_microseconds(circ_duration: List, backend_dt: float):
    dt = backend_dt * 1e6  # dt in microseconds
    return list(map(lambda x: x * dt, circ_duration))

dt = backend.target.dt
plt.plot(
    depths,
    convert_dt_to_microseconds(unitary_durations, dt),
    color="#be95ff",
    linestyle=":",
    label="unitary",
)
plt.plot(
    depths,
    convert_dt_to_microseconds(dynamic_durations, dt),
    color="#ff7eb6",
    linestyle="-.",
    label="dynamic",
)
plt.plot(
    depths,
    convert_dt_to_microseconds(dynamic_durations_meas2, dt),
    color="#ff7eb6",
    linestyle="-.",
    marker="s",
    mfc="none",
    label="dynamic w/ meas2",
)

plt.xlabel("Trotter steps")
plt.ylabel(r"Circuit durations in $\mu$s")
plt.legend()

<matplotlib.legend.Legend at 0x17f73c6e0>

Setelah job selesai, kita mengambil data di bawah dan menghitung magnetisasi rata-rata yang diperkirakan oleh observable observables_unitary atau observables_dynamic yang sudah kita buat sebelumnya.

runs = {
    "unitary": (
        job_unitary,
        [observables_unitary] * len(circuits_unitary),
    ),
    "unitary_dd": (
        job_unitary_dd,
        [observables_unitary] * len(circuits_unitary),
    ),
    # Omitting Dyn w/o DD and Dynamic w/ DD plots for better readability
    # "dynamic": (job_dynamic, observables_dynamic),
    # "dynamic_dd": (job_dynamic_dd, observables_dynamic),
    "dynamic_meas2": (job_dynamic_meas2, observables_dynamic),
    "dynamic_dd_meas2": (
        job_dynamic_dd_meas2,
        observables_dynamic,
    ),
}

data_dict = {}
for key, (job, obs) in runs.items():
    data = []
    for i in range(points):
        data.append(
            [
                job.result()[ind].data["meas"].expectation_values(obs[ind])[i]
                for ind in depths
            ]
        )
    data_dict[key] = data

Di bawah ini kita memplot magnetisasi spin sebagai fungsi dari Trotter steps pada nilai $\theta$ yang berbeda, yang sesuai dengan kekuatan medan magnet lokal yang berbeda. Kita memplot hasil simulasi MPS yang sudah dihitung sebelumnya untuk unitary ideal circuits, bersama dengan hasil eksperimental dari:

1. menjalankan unitary circuits dengan DD
2. menjalankan dynamic circuits dengan DD dan MidCircuitMeasure

plt.figure(figsize=(10, 6))

colors = ["#0f62fe", "#be95ff", "#ff7eb6"]
for i in range(points):
    plt.plot(
        depths,
        data_sim[i],
        color=colors[i],
        linestyle="solid",
        label=f"θ={pi_check(i*max_angle/(points-1))} (MPS)",
    )
    # plt.plot(
    #     depths,
    #     data_dict["unitary"][i],
    #     color=colors[i],
    #     linestyle=":",
    #     label=f"θ={pi_check(i*max_angle/(points-1))} (Unitary)",
    # )

    plt.plot(
        depths,
        data_dict["unitary_dd"][i],
        color=colors[i],
        marker="o",
        mfc="none",
        linestyle=":",
        label=f"θ={pi_check(i*max_angle/(points-1))} (Unitary w/DD)",
    )

    # Omitting Dyn w/o DD and Dynamic w/ DD plots for better readability
    # plt.plot(
    #     depths,
    #     data_dict["dynamic"][i],
    #     color=colors[i],
    #     linestyle="-.",
    #     label=f"θ={pi_check(i*max_angle/(points-1))} (Dyn w/o DD)",
    # )
    # plt.plot(
    #     depths,
    #     data_dict["dynamic_dd"][i],
    #     marker="D",
    #     mfc="none",
    #     color=colors[i],
    #     linestyle="-.",
    #     label=f"θ={pi_check(i*max_angle/(points-1))} (Dynamic w/ DD)",
    # )

    # plt.plot(
    #     depths,
    #     data_dict["dynamic_meas2"][i],
    #     color=colors[i],
    #     marker="s",
    #     mfc="none",
    #     linestyle=':',
    #     label=f"θ={pi_check(i*max_angle/(points-1))} (Dynamic w/ MidCircuitMeas)",
    # )

    plt.plot(
        depths,
        data_dict["dynamic_dd_meas2"][i],
        color=colors[i],
        marker="*",
        markersize=8,
        linestyle=":",
        label=f"θ={pi_check(i*max_angle/(points-1))} (Dynamic w/ DD & MidCircuitMeas)",
    )

plt.xlabel("Trotter steps", fontsize=16)
plt.ylabel("Average magnetization", fontsize=16)
plt.xticks(rotation=45)
handles, labels = plt.gca().get_legend_handles_labels()
plt.legend(
    handles,
    labels,
    loc="upper right",
    bbox_to_anchor=(1.46, 1.0),
    shadow=True,
    ncol=1,
)
plt.title(
    f"{hex_rows}x{hex_cols} hex ring, {num_qubits} data qubits, {len(ancilla)} ancilla qubits \n{backend.name}: Sampler"
)
plt.show()

Saat kita membandingkan hasil eksperimental dengan simulasi, kita melihat bahwa implementasi dynamic circuit (garis putus-putus dengan bintang) secara keseluruhan memiliki performa yang lebih baik dibandingkan implementasi unitary standar (garis putus-putus dengan lingkaran). Singkatnya, kita mempresentasikan dynamic circuits sebagai solusi untuk mensimulasikan model spin Ising pada kisi honeycomb, topologi yang tidak native untuk hardware. Solusi dynamic circuit memungkinkan interaksi ZZ antara qubit yang bukan tetangga terdekat, dengan kedalaman gate dua-qubit yang lebih pendek dibandingkan menggunakan gate SWAP, dengan konsekuensi memperkenalkan ancilla qubit tambahan dan operasi classical feedforward.

Referensi

[1] Quantum computing with Qiskit, oleh Javadi-Abhari, A., Treinish, M., Krsulich, K., Wood, C.J., Lishman, J., Gacon, J., Martiel, S., Nation, P.D., Bishop, L.S., Cross, A.W. dan Johnson, B.R., 2024. arXiv preprint arXiv:2405.08810 (2024)