動的回路によるキック付きアイジングハミルトニアンのシミュレーション

使用時間推定値: Heron r3 プロセッサ上で 7.5 分。（注：これはあくまで概算です。）実行時間は異なる場合があります。

動的回路とは、古典的なフィードフォワードを備えた回路である。言い換えれば、回路の中間測定に続いて、その古典的出力に基づいて量子操作を決定する古典論理演算が行われる回路である。このチュートリアルでは、六角格子上のスピン系におけるキックド・アイジングモデルをシミュレートし、動的回路を用いてハードウェアの物理的接続性を超えた相互作用を実現する。

アイジングモデルは物理学の様々な分野で広く研究されてきた。格子点間でアイジング相互作用を受けるスピンと、各点における局所磁場からのキックをモデル化する。本チュートリアルで扱うスピンのトロッター化時間発展は、 [1] より引用した以下のユニタリ演算子によって与えられる：

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)

スピンダイナミクスを探るため、我々は各サイトにおけるスピンの平均磁化をトロッターステップの関数として研究する。したがって、我々は以下の観測量を構築する：

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

格子サイト間のZZ相互作用を実現するため、動的回路機能を用いた解法を提案する。これにより、SWAPゲートを用いた標準的なルーティング手法と比較して、2量子ビット深度が大幅に短縮される。一方、動的回路における古典的なフィードフォワード演算は、量子ゲートよりも実行時間が長い傾向にある。したがって、動的回路には限界とトレードオフが存在する。また、ストレッチ持続時間を利用して、古典的なフィードフォワード操作中にアイドル状態の量子ビットに対して動的デカップリングシーケンスを追加する方法も提示する。

要件

このチュートリアルを始める前に、以下のものがインストールされていることを確認してください：

Qiskit SDK v2.0 または、可視化サポート付きの後続バージョン
Qiskit Runtime v0.37 またはそれ以降のバージョンで可視化サポート付き (pip install 'qiskit-ibm-runtime[visualization]')
Rustworkx グラフライブラリ (pip install rustworkx)
Qiskit Aer (pip install qiskit-aer)

のセットアップ

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.providers.exceptions import QiskitBackendNotFoundError
from qiskit_ibm_runtime.visualization import (
    draw_circuit_schedule_timing,
)

ステップ1：古典的な入力を量子回路にマッピングする

シミュレーション対象の格子を定義することから始めます。我々はハニカム（六角形とも呼ばれる）格子を用いて研究を行うことを選択した。これは次数3の頂点を持つ平面グラフである。ここでは、格子のサイズと、トロッター化されたダイナミクスにおいて関心のある関連回路パラメータを指定する。我々は、局所磁場に対して3つの異なる $\theta$ 値のもとで、アイジングモデルにおけるトロッター化時間発展をシミュレートする。

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

まずは小さなテスト例から始めましょう：

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)

Output:

この小さな例を用いて説明とシミュレーションを行います。以下では、ワークフローを大規模なサイズに拡張できることを示すため、大規模な例も構築します。

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()

Output:

num_qubits = 46

単一回路を構築する

問題の規模とパラメータが指定されたので、 $U(\theta)$ のトロッター化時間発展をシミュレートするパラメータ化回路を構築する準備が整った。この depth 回路は引数によって指定される異なるトロッターステップを用いる。構築する回路は、 $\theta$ ゲートと Rzz ゲートが交互に積層 Rxされた層構造を持つ。ゲート Rzz は結合したスピン間のZZ相互作用を実現し、これらは引数 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

小さなテスト回路を可視化する：

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)

Output:

同様に、大規模な例題のユニタリー回路を異なるトロッター化ステップで構築し、期待値を推定するための観測量を構築する。

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,
)

動的回路実装を構築する

このセクションでは、同じトロッター化時間発展をシミュレートするための主要な動的回路実装を示す。注意：シミュレートしたいハニカム格子は、ハードウェア量子ビットの重格子とは一致しません。回路をハードウェアにマッピングする単純な方法の一つは、相互作用する量子ビットを隣り合わせにする一連のSWAP操作を導入し、ZZ相互作用を実現することである。ここでは、動的回路を用いた代替アプローチを解決策として提示する。これは、Qiskit内の回路内で量子計算とリアルタイムの古典計算を組み合わせることで、最近接相互作用を超えた相互作用を実現できることを示している。

動的回路実装において、ZZ相互作用は補助量子ビット、回路途中測定、およびフィードフォワードを用いて効果的に実装される。これを理解するには、ZZ回転が状態のパリティに基づいて位相因子 $e^{i\theta}$ を適用することに留意されたい。 2量子ビットの場合、計算基底状態は $|00\rangle$ 、 $|01\rangle$ 、 $|10\rangle$ 、および $|11\rangle$ である。ZZ回転ゲートは、状態 $|01\rangle$ および $|10\rangle$ （状態内の1の数が奇数である状態）に位相因子を適用し、偶数パリティの状態は変化させない。以下では、動的回路を用いて2つの量子ビット間でZZ相互作用を効果的に実装する方法について説明する。

アンシラ量子ビットにパリティを計算する：2つの量子ビットに直接ZZ演算を適用する代わりに、3つ目の量子ビットであるアンシラ量子ビットを導入し、2つのデータ量子ビットのパリティ情報を格納する。各データ量子ビットから補助量子ビットへCXゲートを用いて、データ量子ビットと補助量子ビットを絡み合わせる。
補助量子ビットに単一量子ビットZ回転を適用する：これは補助量子ビットが2つのデータ量子ビットのパリティ情報を保持しているためであり、これによりデータ量子ビットに対してZZ回転が効果的に実現される。
補助量子ビットをX基底で測定する：これが補助量子ビットの状態を収縮させる重要なステップであり、測定結果は起こったことを示す：
- 測定0：結果が0となる場合、我々は実際にデータ量子ビットに対して $ZZ(\theta)$ 回転を正しく適用したことになる。
- 対策1：結果1が観測された場合、代わりに $ZZ(\theta + \pi)$ を適用した。
補正ゲートを適用するタイミング（測定時）：1. 測定値が1の場合、データ量子ビットにZゲートを適用し、余分な $\pi$ 位相を「修正」する。

結果として得られた回路は以下の通りです：

この手法を用いてハニカム格子をシミュレートすると、結果として得られる回路はヘビーヘックス格子を持つハードウェアに完全に埋め込まれる：全てのデータ量子ビットは格子の degree-3 サイト上に配置され、これが六角格子を形成する。各データ量子ビットのペアは、 degree-2 サイト上に存在する補助量子ビットを共有する。以下に、動的回路実装のための量子ビット格子を構築し、補助量子ビット（濃い紫色の円で示される）を導入する。

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)

データ量子ビットと補助量子ビットの重六角格子を小規模で可視化する：

(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()

Output:

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

以下に、トロッター化時間発展のための動的回路を構築する。ゲート RZZ は、上記の手順を用いて動的回路実装に置き換えられる。

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

動的デカップリング（DD）と `stretch` 持続時間のサポート

動的回路実装を用いてZZ相互作用を実現する際の注意点として、回路中間測定と古典的フィードフォワード操作は、量子ゲートよりも実行に通常より長い時間を要する点が挙げられる。量子ビットのデコヒーレンスを抑制するため、古典演算が行われる待機時間中に、補助量子ビットに対する測定演算の後、データ量子ビットに対する条件付きZ演算の前、すなわち if_test ``文の前に動的デカップリング（DD）シーケンスを追加した。

DDシーケンスは、この関数によって追加されます add_stretch_dd()。この関数は、各DDゲートの間の時間間隔を決定するために、各ゲートの持続 stretch 時間を使用します。「持続 stretch 時間」とは、遅延時間が量子ビットのアイドル時間を埋めるまで長くなるように、その delay 操作に対して伸縮可能な時間間隔を指定する方法である。で指定された持続時間変数は、コンパイル時に、特定の制約を満たす所望の持続時間へと解決 stretch されます。これは、優れた誤差抑制性能を実現するためにDDシーケンスのタイミングが極めて重要となる場合に、非常に有用です。この型 stretch に関する詳細については、 OpenQASM のドキュメントを参照してください。現在、この型の stretch サポートは実験的な段階にあります。使用上の制約に関する詳細については、ドキュメント stretch の「制限事項」のセクションをご参照ください。

上記で定義した関数を用いて、DDの有無に応じたトロッター化時間発展回路と、それに対応する観測量を構築する。

まず、小さな例題の動的回路を可視化することから始めます：

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)

Output:

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)

Output:

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)

Output:

同様に、大規模な例に対する動的回路を構築する：

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)

ステップ2: ハードウェア実行に向けた問題の最適化

回路をハードウェアへトランスパイルする準備が整いました。単一標準実装と動的回路実装の両方をハードウェアへトランスパイルします。

ハードウェアへトランスパイルするには、まずバックエンドをインスタンス化します。利用可能な場合、(measure_2) MidCircuitMeasure 命令をサポートするバックエンドを選択します。

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,
    )

動的回路のためのトランスパイル

まず、動的回路をトランスパイルします。DDシーケンスを追加する場合と追加しない場合の両方で。一貫した結果を得るため、全ての回路で同一の物理量子ビットセットを使用することを保証するには、まず回路を一度トランスパイルし、その後パ initial_layout スマネージャーで指定される後続の全ての回路に対してそのレイアウトを使用します。次に、プリミティブ統一ブロック（PUB）をサンプリングプリミティブの入力として構築する。

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]

以下のトランスパイルされた回路の量子ビット配置を可視化できます。黒い円は、動的回路実装で使用されるデータ量子ビットと補助量子ビットを示している。

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",
)

Output:

Note

neato/usr/lib/lib64/lib64-lib64-lib64-lib64-lib64-lib64-lib64-lib64-lib64-lib64-lib64-lib64-lib64-lib64-lib64-lib64-lib64-lib64-lib64-lib64-lib64-libgraphviz``plot_circuit_layout()デフォルト以外の場所（例： MacOS を使用homebrew）にインストールした場合、環境変数PATH` を更新する必要がある場合があります。このノートブック内で以下の方法で行うことができます:

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)

Output:

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

Output:

トランスパイルを使用して `MidCircuitMeasure`

MidCircuitMeasure これは、既存の測定機能に追加されたもので、回路中間部の測定を行うために特別に校正されています。その MidCircuitMeasure 命令は、バックエンドでサポートされている命令 measure_2 に対応しています。なお、 measure_2 この機能はすべてのバックエンドでサポートされているわけではありません。を使用して service.backends(filters=lambda b: "measure_2" in b.supported_instructions) 、それをサポートしているバックエンドを検索できます。ここでは、バックエンドが対応している場合、回路内で定義された回路中間測定が `` MidCircuitMeasure 演算子を使用して実行されるように、回路をトランスパイルする方法を示す。

以下に、命令 measure_2 の所要時間と標準 measure 命令の所要時間を印刷します。

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

Output:

Mid-circuit measurement `measure_2` duration: 1.3800000000000003 μs
Terminal measurement `measure` duration: 2.1800000000000006 μ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
]

単一回路のためのトランスパイレーション

動的回路とそのユニタリー対応物との公正な比較を確立するため、データ量子ビットとして動的回路で使用された物理量子ビットの同一セットを、ユニタリー回路のトランスパイル用レイアウトとして使用する。

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]

トランスパイルされたユニタリ回路の量子ビット配置を可視化する。黒い円はユニタリー回路をトランスパイルするために使用される物理量子ビットを示し、それらのインデックスは仮想量子ビットのインデックスに対応する。これを動的回路用に描かれたレイアウトと比較することで、ユニタリー回路が動的回路のデータ量子ビットと同じ物理量子ビットのセットを使用していることを確認できる。

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

Output:

次に、トランスパイルされた回路にDDシーケンスを追加し、ジョブ提出用の対応するPUBを構築します。

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]

ユニタリ回路と動的回路における2量子ビットゲートの深さの比較

# 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()

Output:

<matplotlib.legend.Legend at 0x12628b0e0>

測定ベース回路の主な利点は、複数のZZ相互作用を実装する際に、CX層を並列化でき、測定を同時に行える点である。これは、すべてのZZ相互作用が可換であるため、測定深度1で計算を実行できるからである。回路をトランスパイルした後、動的回路アプローチは標準的なユニタリアプローチに比べて2量子ビット深度が大幅に短いことが確認された。ただし、追加の中間回路測定と古典的なフィードフォワード自体が時間を要し、独自の誤差源をもたらす点に留意が必要である。

ステップ3: `Qiskit primitives`を使用して実行する

ローカルテストモード

ハードウェアにジョブを送信する前に、ローカルテストモードを使用して動的回路の簡易テストシミュレーションを実行できます。

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])

Output:

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

array([ 0.16666667,  0.01529948, -0.14290365])

MPSシミュレーション

大規模回路の場合、選択した結合次元に応じて期待値に対する近似結果を提供するMPS matrix_product_state シミュレータを使用できます。その後、MPSシミュレーション結果をベースラインとして使用し、ハードウェアからの結果と比較する。

# 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)

回路と観測量が準備されたので、サンプラープリミティブを用いてハードウェア上でそれらを実行する。

ここでは、 dynamic_pubs、、 unitary_pubs に対して3つの dynamic_pubs_ddジョブを提出します。それぞれは、9種類の異なるトロッターステップと3種類の異なる $\theta$ パラメータに対応する、パラメータ化された回路のリストである。

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()}")

Output:

unitary: d96s4b52su3c739hakrg
unitary_dd: d96s4bt2su3c739haksg
dynamic: d96s4c0tcv6s73dk55mg
dynamic_dd: d96s4ckqp3as739qvid0
dynamic_meas2: d96s4csqp3as739qvie0
dynamic_dd_meas2: d96s4daf47jc73a5v8s0

ステップ4：後処理を行い、結果を希望の古典的な形式で返す

ジョブが完了した後、ジョブ結果のメタデータから実行時間を取得し、実行スケジュールの情報を可視化することができます。回路のスケジューリング情報の可視化について詳しくは、こちらのページをご覧ください。

# 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")

ユニタリー回路と動的回路の回路持続時間をプロットする。下図から、中間回路測定と古典的操作に要する時間にもかかわらず、動的回路実装はユニタリー実装と同等の measure_2 回路実行時間をもたらすことがわかる。

# 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()

Output:

<matplotlib.legend.Legend at 0x12bfde270>

ジョブが完了した後、以下のデータを取得し、先に構築した観測量 observables_unitary``observables_dynamic またはによって推定される平均磁化を計算する。

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

以下に、局所磁場の強度に対応する異なる $\theta$ 値におけるトロッターステップ数に対するスピン磁化をプロットする。ユニタリー理想回路に対する事前計算されたMPSシミュレーション結果と、以下の実験結果を共にプロットする：

DDを用いたユニタリー回路の実行
DDを用いた動的回路の実行 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))} "
        f"(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, "
    f"{len(ancilla)} ancilla qubits \n{backend.name}: Sampler"
)
plt.show()

Output:

実験結果とシミュレーション結果を比較すると、動的回路実装（星印の付いた点線）が標準ユニタリー実装（丸印の付いた点線）よりも全体的に優れた性能を示していることがわかる。要約すると、我々はハニカム格子上でのアイジングスピンモデルのシミュレーション手法として動的回路を提案する。このトポロジーはハードウェアに固有のものではない。動的回路ソリューションは、最隣接でない量子ビット間のZZ相互作用を可能とし、SWAPゲートを使用する場合よりも短い2量子ビットゲート深度を実現する。ただし、追加の補助量子ビットと古典的なフィードフォワード操作を導入する代償を伴う。

参照

[1] Qiskitを用いた量子コンピューティング、Javadi-Abhari, A. 著トレイニッシュ, M., クルシリッチ, K., ウッド、 C.J リシュマン, J., ガコン, J., マルティエル, S., ネイション、 P.D ビショップ、 L.S クロス、 A.W。およびジョンソン、 B.R 2024. arXiv プレプリント arXiv:2405.08810 (2024)

このページは役に立ちましたか？

バグや誤字の報告、またはコンテンツの要求はGitHubで行ってください。

動的回路によるキック付きアイジングハミルトニアンのシミュレーション

要件

のセットアップ

ステップ1：古典的な入力を量子回路にマッピングする

単一回路を構築する

動的回路実装を構築する

動的デカップリング（DD）と stretch 持続時間のサポート

ステップ2: ハードウェア実行に向けた問題の最適化

動的回路のためのトランスパイル

トランスパイルを使用して MidCircuitMeasure

単一回路のためのトランスパイレーション

ユニタリ回路と動的回路における2量子ビットゲートの深さの比較

ステップ3: Qiskit primitivesを使用して実行する

ローカルテストモード

MPSシミュレーション

ステップ4：後処理を行い、結果を希望の古典的な形式で返す

参照

動的デカップリング（DD）と `stretch` 持続時間のサポート

トランスパイルを使用して `MidCircuitMeasure`

ステップ3: `Qiskit primitives`を使用して実行する