Skip to main content
IBM Quantum Platform

QESEM関数による傾斜場イジングモデル( 2D )のシミュレーション

Note

Qiskit Functions は、 IBM Quantum® Premium Plan、Flexプラン、 On-Prem ( IBM Quantum Platform API経由)プランのユーザーだけが利用できる実験的な機能です。 これらはプレビューリリースであり、変更される可能性がある。

使用時間の目安:Heron r2 プロセッサーで20分。 (注:これはあくまでも目安です。 実行時間は異なるかもしれない)。


背景

このチュートリアルでは、Qedma の Qiskit 関数である QESEM を使用して、非クリフォード角を持つ 2D 傾斜場イジング (TFI) モデルという標準的な量子スピン モデルのダイナミクスをシミュレートする方法を説明します。

H=Ji,jZiZj+gxiXi+gziZi,H = J \sum_{\langle i,j \rangle} Z_i Z_j + g_x \sum_i X_i + g_z \sum_i Z_i ,

ここで、 i,j\langle i,j \rangle は格子上の最近傍を表す。 多体量子系の時間発展をシミュレートすることは、古典的なコンピュータにとって計算上困難な課題である。 これに対して量子コンピューターは、当然ながらこのタスクを効率的に実行するように設計されている。 特にTFIモデルは、その豊かな物理的挙動とハードウェアに優しい実装により、量子ハードウェアのベンチマークとして人気を博している。

連続時間ダイナミクスをシミュレートする代わりに、密接に関連するキックド・イジング・モデルを採用する。 ダイナミクスは周期的量子回路として正確に表現することができ、各進化ステップは3層の分数2量子ビットゲート RZZ(αZZ)R_{ZZ} (\alpha_{ZZ})、1量子ビットゲート RX(αX)R_X (\alpha_X)RZ(αZ)R_Z (\alpha_Z) の層でインターリーブされている。

私たちは、古典的なシミュレーションとエラー軽減の両方に挑戦的な一般的な角度を使用します。 具体的には、 αZZ=1.0\alpha_{ZZ} = 1.0αX=0.53\alpha_X = 0.53αZ=0.1\alpha_Z = 0.1 を選択し、どの可積分点からも遠い位置にモデルを置いた。

このチュートリアルでは以下のことを行う:

  • QESEMの解析的および経験的な時間見積もり機能を使用して、完全なエラー緩和を行った場合のQPUの予想実行時間を見積もります。
  • ハードウェアにインスパイアされた量子ビットレイアウトとゲート層を使用して、 2D 傾斜磁場イジングモデル回路を構築し、シミュレーションを行う。
  • デバイスの量子ビットの接続性と、実験用に選択したサブグラフを可視化します。
  • オペレータ逆伝播法(OBP) を用いた回路深度の削減手法を実証する。 この手法は、演算者をより多く測定する代償として、回路端からの演算を削減する。
  • QESEMを使用して、複数の観測値に対して同時に不偏の誤差緩和(EM)を実行し、理想的な結果、ノイズのある結果、緩和された結果を比較します。
  • 異なる回路深度にわたって、エラー緩和が磁化に与える影響を分析し、プロットする。

注意: OBPは一般的に、コミットしない観測値の集合を返します。 QESEMは、観測対象が非共約項を含む場合、測定ベースを自動的に最適化します。 いくつかの発見的アルゴリズムを用いて測定基底セットの候補を生成し、異なる基底の数を最小化するセットを選択する。 これは、QESEMが互換性のある観測値を共通のベースにグループ化することで、必要な測定コンフィギュレーションの総数を減らし、効率を向上させることを意味する。


QESEMについて

QESEMは信頼性の高い高精度な特性評価ベースのソフトウェアで、効率的で不偏的な準確率論的誤差緩和を実装しています。 一般的な量子回路のエラーを軽減するように設計されており、アプリケーションを問わない。 これは、 IBM® EagleおよびHeronデバイスでの実用規模の実験を含む、多様なハードウェアプラットフォームで検証されている。 QESEMのワークフローの段階は以下の通りである:

  1. デバイスの特性評価 - ゲートの忠実度をマップし、コヒーレントエラーを特定し、リアルタイムのキャリブレーションデータを提供します。 この段階は、ミティゲーションが利用可能な最も忠実度の高いオペレーションを活用することを保証する。
  2. ノイズを考慮したトランスパイル - 代替の量子ビットマッピング、演算セット、測定ベースを生成して評価し、QPUの推定実行時間を最小化する変種を選択します。
  3. エラー抑制 - ネイティブ・ゲートの再定義、パウリ・ツワーリングの適用、パルス・レベル制御の最適化(サポート対象プラットフォーム)により、忠実度を向上。
  4. 回路特性評価 - テーラーメイドのローカル誤差モデルを構築し、QPU測定値に適合させて残留ノイズを定量化。
  5. エラー・ミティゲーション - マルチタイプの準確率的分解を構築し、ミティゲーションのQPU時間とハードウェアの変動に対する感度を最小化する適応的なプロセスでサンプリングし、大規模な回路ボリュームで高い精度を達成します。

QESEMに関する詳細、およびネイティブ・ヘビー・ヘックス・ジオメトリの103キュービット・高接続性サブグラフを用いた実用規模の実験については、 ibm_marrakesh実用規模の量子回路における信頼性の高い高精度エラー緩和 」を参照してください。

QESEMワークフロー。

要件

ノートブックを実行する前に、以下の Python パッケージをインストールしてください:

  • Qiskit SDK v2.0.0 またはそれ以降 (pip install qiskit)
  • Qiskit Runtime v0.40.0 またはそれ以降 (pip install qiskit-ibm-runtime)
  • Qiskit Functions Catalog v0.8.0 またはそれ以降 ( pip install qiskit-ibm-catalog )
  • Operator Backpropagation Qiskit アドオン v0.3.0 またはそれ以降 ( pip install qiskit-addon-obp )
  • Qiskit Utils アドオン v0.1.1 またはそれ以降 ( pip install qiskit-addon-utils )
  • Qiskit Aer シミュレータ v0.17.1 またはそれ以降 ( pip install qiskit-aer )
  • Matplotlib v3.10.3 またはそれ以降 ( pip install matplotlib )

セットアップ

まず、関連するライブラリをインポートする:

%matplotlib inline

from typing import Sequence

import matplotlib.pyplot as plt
import numpy as np

import qiskit
from qiskit_ibm_runtime import EstimatorV2 as Estimator
from qiskit_ibm_catalog import QiskitFunctionsCatalog
from qiskit_aer import AerSimulator
from qiskit_addon_utils.slicing import combine_slices, slice_by_gate_types
from qiskit_addon_obp import backpropagate
from qiskit_addon_obp.utils.simplify import OperatorBudget
from qiskit.visualization import (
    plot_gate_map,
)

次に、 IBM Quantum Platform のダッシュボードからAPIキーを使用して認証を行ってください。 次に、次のようにQiskit関数を選択します。 (セキュリティ上の理由から、信頼できるマシンを使用している場合は、認証のたびにAPIキーを入力する必要がないよう、 アカウントの認証情報をローカル環境に保存しておくことをお勧めします。)

# Paste here your instance and token strings

instance = "YOUR_INSTANCE"
token = "YOUR_TOKEN"
channel = "ibm_quantum_platform"

catalog = QiskitFunctionsCatalog(
    channel=channel, token=token, instance=instance
)
qesem_function = catalog.load("qedma/qesem")

ステップ1:古典的な入力を量子問題にマッピングする

まず、トロッター回路を作る関数を定義する:

def trotter_circuit_from_layers(
    steps: int,
    theta_x: float,
    theta_z: float,
    theta_zz: float,
    layers: Sequence[Sequence[tuple[int, int]]],
    init_state: str | None = None,
) -> qiskit.QuantumCircuit:
    """
    Generates an ising trotter circuit
    :param steps: trotter steps
    :param theta_x: RX angle
    :param theta_z: RZ angle
    :param theta_zz: RZZ angle
    :param layers: list of layers (can be list of layers in device)
    :param init_state: Initial state to prepare.
     If None, will not prepare any state. If "+", will
     add Hadamard gates to all qubits.
    :return: QuantumCircuit
    """
    qubits = sorted({i for layer in layers for edge in layer for i in edge})
    circ = qiskit.QuantumCircuit(max(qubits) + 1)

    if init_state == "+":
        print("init_state = +")
        for q in qubits:
            circ.h(q)

    for _ in range(steps):
        for q in qubits:
            circ.rx(theta_x, q)
            circ.rz(theta_z, q)

        for layer in layers:
            for edge in layer:
                circ.rzz(theta_zz, *edge)
        circ.barrier(qubits)

    return circ

次に、 AerSimulator を使って理想的な期待値を計算する関数を作る。

なお、大規模な回路(30キュービット以上)については、信念伝播(BP)法を用いたPEPSシミュレーションから得られた事前計算値を使用することを推奨します。 このコードには、 本論文で紹介されたPEPSテンソルネットワークを進化させるためのBPアプローチ(以下、PEPS-BPと呼ぶ)に基づき、テンソルネットワーク Python パッケージ 「quimb 」を用いて、例として35キュービット分の事前計算済み値が含まれています。

def calculate_ideal_evs(circ, obs, num_qubits, step):
    # Predefined results for large circuits - calculated using
    # bppeps for 3, 5, 7, 9 trotter steps
    predefined_35 = [
        0.79537,
        0.78653,
        0.79699,
    ]

    if num_qubits == 35:
        print(
            "Using precalculated ideal values for large circuits calculated "
            "with belief propagation PEPS. Currently only for 35 qubits."
        )
        return predefined_35[step]

    else:
        simulator = AerSimulator()

        # Use Estimator primitive to get expectation value
        estimator = Estimator(simulator)
        sim_result = estimator.run([(circ, [obs])], precision=0.0001).result()

        # Extracting the result
        ideal_values = sim_result[0].data.evs[0]
        return ideal_values

私たちは、Heronデバイスから取得したハードウェアベースの RZZR_{ZZ} レイヤーマッピングを使用しています。そこから、シミュレーションしたい量子ビットの数に応じてレイヤーを切り出します。 10,21,28,35量子ビットのサブグラフを定義し、 2D (お好きなサブグラフに自由に変更してください):

LAYERS_HERON_R2 = [  # the full set of hardware layers for Heron r2
    [
        (2, 3),
        (6, 7),
        (10, 11),
        (14, 15),
        (20, 21),
        (16, 23),
        (24, 25),
        (17, 27),
        (28, 29),
        (18, 31),
        (32, 33),
        (19, 35),
        (36, 41),
        (42, 43),
        (37, 45),
        (46, 47),
        (38, 49),
        (50, 51),
        (39, 53),
        (60, 61),
        (56, 63),
        (64, 65),
        (57, 67),
        (68, 69),
        (58, 71),
        (72, 73),
        (59, 75),
        (76, 81),
        (82, 83),
        (77, 85),
        (86, 87),
        (78, 89),
        (90, 91),
        (79, 93),
        (94, 95),
        (100, 101),
        (96, 103),
        (104, 105),
        (97, 107),
        (108, 109),
        (98, 111),
        (112, 113),
        (99, 115),
        (116, 121),
        (122, 123),
        (117, 125),
        (126, 127),
        (118, 129),
        (130, 131),
        (119, 133),
        (134, 135),
        (140, 141),
        (136, 143),
        (144, 145),
        (137, 147),
        (148, 149),
        (138, 151),
        (152, 153),
        (139, 155),
    ],
    [
        (1, 2),
        (3, 4),
        (5, 6),
        (7, 8),
        (9, 10),
        (11, 12),
        (13, 14),
        (21, 22),
        (23, 24),
        (25, 26),
        (27, 28),
        (29, 30),
        (31, 32),
        (33, 34),
        (40, 41),
        (43, 44),
        (45, 46),
        (47, 48),
        (49, 50),
        (51, 52),
        (53, 54),
        (55, 59),
        (61, 62),
        (63, 64),
        (65, 66),
        (67, 68),
        (69, 70),
        (71, 72),
        (73, 74),
        (80, 81),
        (83, 84),
        (85, 86),
        (87, 88),
        (89, 90),
        (91, 92),
        (93, 94),
        (95, 99),
        (101, 102),
        (103, 104),
        (105, 106),
        (107, 108),
        (109, 110),
        (111, 112),
        (113, 114),
        (120, 121),
        (123, 124),
        (125, 126),
        (127, 128),
        (129, 130),
        (131, 132),
        (133, 134),
        (135, 139),
        (141, 142),
        (143, 144),
        (145, 146),
        (147, 148),
        (149, 150),
        (151, 152),
        (153, 154),
    ],
    [
        (3, 16),
        (7, 17),
        (11, 18),
        (22, 23),
        (26, 27),
        (30, 31),
        (34, 35),
        (21, 36),
        (25, 37),
        (29, 38),
        (33, 39),
        (41, 42),
        (44, 45),
        (48, 49),
        (52, 53),
        (43, 56),
        (47, 57),
        (51, 58),
        (62, 63),
        (66, 67),
        (70, 71),
        (74, 75),
        (61, 76),
        (65, 77),
        (69, 78),
        (73, 79),
        (81, 82),
        (84, 85),
        (88, 89),
        (92, 93),
        (83, 96),
        (87, 97),
        (91, 98),
        (102, 103),
        (106, 107),
        (110, 111),
        (114, 115),
        (101, 116),
        (105, 117),
        (109, 118),
        (113, 119),
        (121, 122),
        (124, 125),
        (128, 129),
        (132, 133),
        (123, 136),
        (127, 137),
        (131, 138),
        (142, 143),
        (146, 147),
        (150, 151),
        (154, 155),
        (0, 1),
        (4, 5),
        (8, 9),
        (12, 13),
        (54, 55),
        (15, 19),
    ],
]

subgraphs = {  # the subgraphs for the different qubit counts such that it's 2D
    10: list(range(22, 29)) + [16, 17, 37],
    21: list(range(3, 12)) + list(range(23, 32)) + [16, 17, 18],
    28: list(range(3, 12))
    + list(range(23, 32))
    + list(range(45, 50))
    + [16, 17, 18, 37, 38],
    35: list(range(3, 12))
    + list(range(21, 32))
    + list(range(41, 50))
    + [16, 17, 18, 36, 37, 38],
    42: list(range(3, 12))
    + list(range(21, 32))
    + list(range(41, 50))
    + list(range(63, 68))
    + [16, 17, 18, 36, 37, 38, 56, 57],
}

n_qubits = 35  # 21, 28, 35, 42
layers = [
    [
        edge
        for edge in layer
        if edge[0] in subgraphs[n_qubits] and edge[1] in subgraphs[n_qubits]
    ]
    for layer in LAYERS_HERON_R2
]

print(layers)

Output:

[[(6, 7), (10, 11), (16, 23), (24, 25), (17, 27), (28, 29), (18, 31), (36, 41), (42, 43), (37, 45), (46, 47), (38, 49)], [(3, 4), (5, 6), (7, 8), (9, 10), (21, 22), (23, 24), (25, 26), (27, 28), (29, 30), (43, 44), (45, 46), (47, 48)], [(3, 16), (7, 17), (11, 18), (22, 23), (26, 27), (30, 31), (21, 36), (25, 37), (29, 38), (41, 42), (44, 45), (48, 49), (4, 5), (8, 9)]]

ここで、選択された部分グラフのHeronデバイス上の量子ビットのレイアウトを可視化する:

catalog = QiskitFunctionsCatalog(
    channel=channel,
    token=token,
    instance=instance,
)
backend = catalog.backend("ibm_fez")  # or any available device

selected_qubits = subgraphs[n_qubits]
num_qubits = backend.configuration().num_qubits
qubit_color = [
    "#ff7f0e" if i in selected_qubits else "#d3d3d3"
    for i in range(num_qubits)
]

plot_gate_map(
    backend=backend,
    figsize=(15, 10),
    qubit_color=qubit_color,
)
plt.show()

Output:

Output of the previous code cell

選択された量子ビットレイアウトの接続性は必ずしも線形ではなく、選択された量子ビット数に応じてHeronデバイスの広い領域をカバーできることに注意してください。

ここで、トロッター回路と、選択した量子ビット数とパラメータに対する平均磁化観測値を生成する:

# Chosen parameters:
theta_x = 0.53
theta_z = 0.1
theta_zz = 1.0
steps = 9

circ = trotter_circuit_from_layers(steps, theta_x, theta_z, theta_zz, layers)
print(
    f"Circuit 2q layers: "
    f"{circ.depth(filter_function=lambda instr: len(instr.qubits) == 2)}"
)
print("\nCircuit structure:")

circ.draw("mpl", scale=0.8, fold=-1, idle_wires=False)
plt.show()

observable = qiskit.quantum_info.SparsePauliOp.from_sparse_list(
    [("Z", [q], 1 / n_qubits) for q in subgraphs[n_qubits]],
    np.max(subgraphs[n_qubits]) + 1,
)  # Average magnetization observable

print(observable)
obs_list = [observable]

Output:

Circuit 2q layers: 27

Circuit structure:
Output of the previous code cell
SparsePauliOp(['IIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIZIII', 'IIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIZIIII', 'IIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIZIIIII', 'IIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIZIIIIII', 'IIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIZIIIIIII', 'IIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIZIIIIIIII', 'IIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIZIIIIIIIII', 'IIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIZIIIIIIIIII', 'IIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIZIIIIIIIIIII', 'IIIIIIIIIIIIIIIIIIIIIIIIIIIIZIIIIIIIIIIIIIIIIIIIII', 'IIIIIIIIIIIIIIIIIIIIIIIIIIIZIIIIIIIIIIIIIIIIIIIIII', 'IIIIIIIIIIIIIIIIIIIIIIIIIIZIIIIIIIIIIIIIIIIIIIIIII', 'IIIIIIIIIIIIIIIIIIIIIIIIIZIIIIIIIIIIIIIIIIIIIIIIII', 'IIIIIIIIIIIIIIIIIIIIIIIIZIIIIIIIIIIIIIIIIIIIIIIIII', 'IIIIIIIIIIIIIIIIIIIIIIIZIIIIIIIIIIIIIIIIIIIIIIIIII', 'IIIIIIIIIIIIIIIIIIIIIIZIIIIIIIIIIIIIIIIIIIIIIIIIII', 'IIIIIIIIIIIIIIIIIIIIIZIIIIIIIIIIIIIIIIIIIIIIIIIIII', 'IIIIIIIIIIIIIIIIIIIIZIIIIIIIIIIIIIIIIIIIIIIIIIIIII', 'IIIIIIIIIIIIIIIIIIIZIIIIIIIIIIIIIIIIIIIIIIIIIIIIII', 'IIIIIIIIIIIIIIIIIIZIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII', 'IIIIIIIIZIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII', 'IIIIIIIZIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII', 'IIIIIIZIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII', 'IIIIIZIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII', 'IIIIZIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII', 'IIIZIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII', 'IIZIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII', 'IZIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII', 'ZIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII', 'IIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIZIIIIIIIIIIIIIIII', 'IIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIZIIIIIIIIIIIIIIIII', 'IIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIZIIIIIIIIIIIIIIIIII', 'IIIIIIIIIIIIIZIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII', 'IIIIIIIIIIIIZIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII', 'IIIIIIIIIIIZIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII'],
              coeffs=[0.02857143+0.j, 0.02857143+0.j, 0.02857143+0.j, 0.02857143+0.j,
 0.02857143+0.j, 0.02857143+0.j, 0.02857143+0.j, 0.02857143+0.j,
 0.02857143+0.j, 0.02857143+0.j, 0.02857143+0.j, 0.02857143+0.j,
 0.02857143+0.j, 0.02857143+0.j, 0.02857143+0.j, 0.02857143+0.j,
 0.02857143+0.j, 0.02857143+0.j, 0.02857143+0.j, 0.02857143+0.j,
 0.02857143+0.j, 0.02857143+0.j, 0.02857143+0.j, 0.02857143+0.j,
 0.02857143+0.j, 0.02857143+0.j, 0.02857143+0.j, 0.02857143+0.j,
 0.02857143+0.j, 0.02857143+0.j, 0.02857143+0.j, 0.02857143+0.j,
 0.02857143+0.j, 0.02857143+0.j, 0.02857143+0.j])

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

OBPの有無によるQPU時間推定

ユーザーは一般に、自分の実験にどれくらいのQPU時間が必要かを知りたがっている。 しかし、これは古典的なコンピューターにとっては難しい問題と考えられている。

QESEMには、実験の実行可能性をユーザーに知らせるための2つの時間推定モードがあります:

  1. 分析的な時間見積もり - 非常に大まかな見積もりが可能で、QPUの時間を必要としない。 これは、トランスパイルパスがQPU時間を短縮する可能性があるかどうかをテストするために使用できる。
  2. 経験的な時間見積もり(ここで実証) - かなり良い見積もりが得られ、QPUの時間を数分使う。

どちらの場合も、QESEMはすべての観測量について必要な精度に到達するまでの推定時間を出力します。

run_on_real_hardware = True

precision = 0.05
if run_on_real_hardware:
    backend_name = "ibm_fez"
else:
    backend_name = "fake_fez"

# Start a job for empirical time estimation
estimation_job_wo_obp = qesem_function.run(
    pubs=[(circ, obs_list)],
    instance=instance,
    backend_name=backend_name,  # E.g. "ibm_brisbane"
    options={
        # "empirical" - gets actual time estimates without running full mitigation
        "estimate_time_only": "empirical",
        "max_execution_time": 120,  # Limits the QPU time, specified in seconds.
        "default_precision": precision,
    },
)
print(estimation_job_wo_obp.job_id)
print(estimation_job_wo_obp.status())

Output:

17d3828e-9fdb-482e-8e9b-392f3eefe313
DONE
# Get the result object (blocking method).
# Use job.status() in a loop for non-blocking.
# This takes 1-3 minutes
result = estimation_job_wo_obp.result()
print(
    f"Empirical time estimation (sec): {result[0].metadata['time_estimation_sec']}"
)

Output:

Empirical time estimation (sec): 1200

次に、オペレータバックプロパゲーション(OBP)を使用します。 (OBP Qiskitアドオンの詳細については、 OBP のドキュメントを参照してください。) バックプロパゲーション用の回路スライスを生成する関数を作成します:

def run_backpropagation(circ_vec, observable, steps_vec, max_qwc_groups=8):
    """
    Runs backpropagation for a list of circuits and observables.
    Returns lists of backpropagated circuits and observables.
    """
    op_budget = OperatorBudget(max_qwc_groups=max_qwc_groups)
    bp_circuit_vec = []
    bp_observable_vec = []

    for i, circ in enumerate(circ_vec):
        slices = slice_by_gate_types(circ)
        bp_observable, remaining_slices, metadata = backpropagate(
            observable,
            slices,
            operator_budget=op_budget,
        )
        bp_circuit = combine_slices(remaining_slices, include_barriers=True)
        bp_circuit_vec.append(bp_circuit)
        bp_observable_vec.append(bp_observable)
        print(f"n.o. steps: {steps_vec[i]}")
        print(f"Backpropagated {metadata.num_backpropagated_slices} slices.")
        print(
            f"New observable has {len(bp_observable.paulis)} terms, "
            f"which can be combined into "
            f"{len(bp_observable.group_commuting(qubit_wise=True))} groups.\n"
            f"After truncation, the error in our observable is bounded by "
            f"{metadata.accumulated_error(0):.3e}"
        )
        print("-----------------")
    return bp_circuit_vec, bp_observable_vec

この関数を呼び出す:

bp_circ_vec, bp_obs_vec = run_backpropagation([circ], observable, [steps])

Output:

n.o. steps: 9
Backpropagated 11 slices.
New observable has 363 terms, which can be combined into 4 groups.
After truncation, the error in our observable is bounded by 0.000e+00
-----------------
print("The remaining circuit after backpropagation looks as follows:")
bp_circ_vec[-1].draw("mpl", scale=0.8, fold=-1, idle_wires=False)
None

Output:

The remaining circuit after backpropagation looks as follows:
Output of the previous code cell

バックプロパゲーションによって、回路の層が2つ減ったことがわかる。 さて、縮小された回路と拡張された観測値ができたので、バックプロパゲートされた回路に対して時間推定を行ってみよう:

# Start a job for empirical time estimation
estimation_job_obp = qesem_function.run(
    pubs=[(bp_circ_vec[-1], [bp_obs_vec[-1]])],
    instance=instance,
    backend_name=backend_name,
    options={
        "estimate_time_only": "empirical",
        "max_execution_time": 120,
        "default_precision": precision,
    },
)
print(estimation_job_obp.job_id)
print(estimation_job_obp.status())

Output:

8bae699d-a16b-4d39-bbd9-d123fbcce55d
DONE
result_obp = estimation_job_obp.result()
print(
    f"Empirical time estimation (sec): {result_obp[0].metadata['time_estimation_sec']}"
)

Output:

Empirical time estimation (sec): 900

OBPは回路の緩和のための時間コストを削減することがわかる。


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

実際のバックエンドで実行する

次に、いくつかのトロッター・ステップで完全な実験を行う。 量子ビット数、要求精度、最大QPU時間は、利用可能なQPUリソースに応じて変更することができる。 最大QPU時間を制限することは、以下の最終プロットでわかるように、最終的な精度に影響することに注意してほしい。

5、7、9トロッター・ステップの4つの回路を 0.05 の精度で解析し、理想的な期待値、ノイズのある期待値、誤差を考慮した期待値を比較する:

steps_vec = [5, 7, 9]

circ_vec = []
for steps in steps_vec:
    circ = trotter_circuit_from_layers(
        steps, theta_x, theta_z, theta_zz, layers
    )
    circ_vec.append(circ)

ここでも、実行時間を短縮するために、各回路でOBPを実行する:

bp_circ_vec_35, bp_obs_vec_35 = run_backpropagation(
    circ_vec, observable, steps_vec
)

Output:

n.o. steps: 5
Backpropagated 11 slices.
New observable has 363 terms, which can be combined into 4 groups.
After truncation, the error in our observable is bounded by 0.000e+00
-----------------
n.o. steps: 7
Backpropagated 11 slices.
New observable has 363 terms, which can be combined into 4 groups.
After truncation, the error in our observable is bounded by 0.000e+00
-----------------
n.o. steps: 9
Backpropagated 11 slices.
New observable has 363 terms, which can be combined into 4 groups.
After truncation, the error in our observable is bounded by 0.000e+00
-----------------

次に、QESEMのフルジョブをバッチで実行する。 QPUバジェットをよりよく制御するために、各ポイントの最大QPUランタイムを制限する。

run_on_real_hardware = True

precision = 0.05
if run_on_real_hardware:
    backend_name = "ibm_marrakesh"
else:
    backend_name = "fake_fez"
# Running full jobs for:
pubs_list = [
    [(bp_circ_vec_35[i], bp_obs_vec_35[i])] for i in range(len(bp_obs_vec_35))
]
# Initiating multiple jobs for different lengths
job_list = []
for pubs in pubs_list:
    job_obp = qesem_function.run(
        pubs=pubs,
        instance=instance,
        backend_name=backend_name,  # E.g. "ibm_brisbane"
        options={
            "max_execution_time": 300,  # Limits the QPU time, specified in seconds.
            "default_precision": 0.05,
        },
    )
    job_list.append(job_obp)

ここでは各ジョブのステータスをチェックする:

for job in job_list:
    print(job.status())

Output:

DONE
DONE
DONE
DONE

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

すべてのジョブの実行が終了したら、ノイズと軽減された期待値を比較することができる。

ideal_values = []
noisy_values = []
error_mitigated_values = []
error_mitigated_stds = []

for i in range(len(job_list)):
    job = job_list[i]
    result = job.result()  # Blocking - takes 3-5 minutes
    noisy_results = result[0].metadata["noisy_results"]

    ideal_val = calculate_ideal_evs(circ_vec[i], observable, n_qubits, i)
    print("---------------------------------")
    print(f"Ideal: {ideal_val}")
    print(f"Noisy: {noisy_results.evs}")
    print(f"QESEM: {result[0].data.evs} \u00b1 {result[0].data.stds}")

    ideal_values.append(ideal_val)
    noisy_values.append(noisy_results.evs)
    error_mitigated_values.append(result[0].data.evs)
    error_mitigated_stds.append(result[0].data.stds)

Output:

Using precalculated ideal values for large circuits calculated with belief propagation PEPS. Currently only for 35 qubits.
---------------------------------
Ideal: 0.79537
Noisy: 0.7039237951821501
QESEM: 0.7828018244130982 ± 0.013257266977728376
Using precalculated ideal values for large circuits calculated with belief propagation PEPS. Currently only for 35 qubits.
---------------------------------
Ideal: 0.78653
Noisy: 0.6478583812958806
QESEM: 0.7875259197423828 ± 0.02703045139248604
Using precalculated ideal values for large circuits calculated with belief propagation PEPS. Currently only for 35 qubits.
---------------------------------
Ideal: 0.79699
Noisy: 0.6171787879868142
QESEM: 0.6918791909168913 ± 0.0740873782039517

最後に、磁化をステップ数に対してプロットすることができる。 これは、QESEM Qiskit Functionをノイズの多い量子デバイスのバイアスフリーエラーミティゲーションに使用することの利点をまとめたものである。

plt.plot(steps_vec, ideal_values, "--", label="ideal")
plt.scatter(steps_vec, noisy_values, label="noisy")
plt.errorbar(
    steps_vec,
    error_mitigated_values,
    yerr=error_mitigated_stds,
    fmt="o",
    capsize=5,
    label="QESEM mitigation",
)
plt.legend()
plt.xlabel("n.o. steps")
plt.ylabel("Magnetization")

Output:

Text(0, 0.5, 'Magnetization')
Output of the previous code cell
Note

QPUの時間を5分に制限したため、第9ステップには大きな統計エラーバーがある。 このステップを15分間実行すると(経験的な時間推定が示唆するように)、エラーバーが小さくなる。 したがって、緩和された値は理想値に近づくことになる。

このページは役に立ちましたか?
バグや誤字の報告、またはコンテンツの要求はGitHubで行ってください。