calibrating-qubits-pulse.ipynb

Kernel: Python 3

Qiskit Pulseで量子ビットを較正する

Qiskitは量子コンピューティングのオープンソース化されたフレームワークです(参照 1)。Qiskitを使用して、量子回路を構築し、シミュレーションし、量子デバイス上で実行することが可能です。

Qiskit Pulseは、特定のハードウェアに依存しない一般の量子デバイスに対してパルスレベルの制御(例: 入力シグナルの連続時間ダイナミクスの制御)を指定する言語を提供します(参照 2)。

このチュートリアルでは、QiskitとQiskit Pulseを使った単一量子ビットのキャリブレーションと特性評価実験の実装について説明します。これらは通常、デバイスが製造されシステムに取り付けられた直後にラボで行われる最初の実験です。このチュートリアルは教育的であり、学生は二準位系のダイナミクスを実験的に調べることに使えます。すべての単位は標準SI系（つまり、Hz、秒など）です。

それぞれの実験でシステムに関する詳細情報を取得し、通常、それを後続の実験で使用します。そのため、このnotebookは基本的に順番に実行する必要があります。

目次

1. はじめに

2. 周波数スイープを使った量子ビット周波数の探索

3. 較正と$\pi$パルスの利用 3.1 ラビ実験を使った$\pi$パルスの較正 3.2 0か1の決定 3.3 反転回復法を使った$T_1$の測定

4. 量子ビットコヒーレンスの決定 4.1 ラムゼイ実験を使った正確な量子ビット周波数の測定 4.2 ハーンエコーを使った$T_2$の測定

5. 参考文献

1. はじめに

最初の準備として基本的な依存関係を設定します。このnotebookでのキャリブレーション実験には現実のノイズの多いデバイスを使用するため、IBMQアカウントをロードして適切なバックエンドを設定する必要があります。

from qiskit.tools.jupyter import *

from qiskit import IBMQ
IBMQ.load_account()
provider = IBMQ.get_provider(hub='ibm-q', group='open', project='main')
backend = provider.get_backend('ibmq_armonk')

/home/dsche/.local/lib/python3.8/site-packages/qiskit/providers/ibmq/ibmqfactory.py:192: UserWarning: Timestamps in IBMQ backend properties, jobs, and job results are all now in local time instead of UTC.
  warnings.warn('Timestamps in IBMQ backend properties, jobs, and job results '

バックエンドがPulse機能をサポートしていることを確認するため、バックエンドの構成をチェックします。この構成を確認することで、バックエンドのセットアップ構造に関する一般情報を知ることができます。

backend_config = backend.configuration()
assert backend_config.open_pulse, "Backend doesn't support Pulse"

例えば、バックエンド構成内のバックエンド・パルスのサンプリング時間を見つけることができます。これは、キャリブレーションルーチンを構築して実行するときに非常に有用な値になります。

dt = backend_config.dt
print(f"Sampling time: {dt*1e9} ns")    # dtを秒の単位で返すので
                                        # 1e9をかけてナノ秒を得ます

Sampling time: 0.2222222222222222 ns

バックエンドの使い方を学ぶにはバックエンドのデフォルト値を知ることから始まります。バックエンドのデフォルト値には基本的な量子演算子を実行するための量子ビット周波数とデフォルトプログラムの推定値が含まれています。次の方法でそれらにアクセスできます：

backend_defaults = backend.defaults()

2. 周波数スイープを使った量子ビット周波数の探索

まずは量子ビット周波数を探すことから始めます。量子ビット周波数は基底状態と励起状態間のエネルギーの差であり、それぞれ$\vert0\rangle$と$\vert1\rangle$の状態としてラベル付けします。この周波数は、量子ビットに特定の量子演算子を実行するためのパルスを作成するのに重要で、この周波数を求めることがキャリブレーションの最終目標です！

超伝導量子ビットでは、より高いエネルギー準位も利用できますが、どの遷移を励起するかを制御できるようにシステムを非調和に作ります。このようにして、2つのエネルギー準位を分離し各量子ビットを基本的な二準位系として扱う事で、より高いエネルギー状態を無視することができます。

一般的な実験室の設定は、ネットワーク・アナライザーと呼ばれるツールを使用して特定の範囲の周波数をスイープし、吸収のサインを探すことで量子ビット周波数を見つけることができます。この測定により、量子ビット周波数の大まかな推定値が得られます。後半で、ラムゼイ・パルスシーケンスを使用してより正確な測定を行う方法について説明します。

まず、量子ビットを探すためスイープする周波数範囲を定義します。この範囲は任意に広げられるので、backend_defaultsを使って推定量子ビット周波数の周りの40MHzの幅に制限します。そして1MHz単位で周波数を変更します。

import numpy as np

# 単位変換係数 -> バックエンドプロパティはSI系（Hz、秒）で返されるため
GHz = 1.0e9 # Gigahertz
MHz = 1.0e6 # Megahertz
us = 1.0e-6 # Microseconds
ns = 1.0e-9 # Nanoseconds

# 以下の量子ビットについて量子ビット周波数を探索する
qubit = 0

# デフォルトの推定周波数を中心値としてスイープを行う
center_frequency_Hz = backend_defaults.qubit_freq_est[qubit]        # デフォルト周波数はHzで与えられる
                                                                    # 注意：これは将来のリリースで変更予定
print(f"Qubit {qubit} has an estimated frequency of {center_frequency_Hz / GHz} GHz.")

# データを見やすくするスケーリング
scale_factor = 1e-14

# 推定周波数を中心に40MHzのスパンでスイープする
frequency_span_Hz = 40 * MHz
# 1MHzの単位で
frequency_step_Hz = 1 * MHz

# 推定周波数より20 MHz高く、20 MHz低くスイープ
frequency_min = center_frequency_Hz - frequency_span_Hz / 2
frequency_max = center_frequency_Hz + frequency_span_Hz / 2
# 実験のために周波数のnp配列をGHz単位で作成
frequencies_GHz = np.arange(frequency_min / GHz, 
                            frequency_max / GHz, 
                            frequency_step_Hz / GHz)

print(f"The sweep will go from {frequency_min / GHz} GHz to {frequency_max / GHz} GHz \
in steps of {frequency_step_Hz / MHz} MHz.")

Qubit 0 has an estimated frequency of 4.974445646137273 GHz.
The sweep will go from 4.954445646137273 GHz to 4.994445646137273 GHz in steps of 1.0 MHz.

次に、実験に使用するパルスを定義します。ガウシアン・パルスである駆動パルスから始めます。

先述のdtの値を覚えていますか？パルスの持続時間はdtで与えられます。次のセルでは、駆動パルスの長さをdtで定義します。

# サンプル数は16の倍数である必要があります
def get_closest_multiple_of_16(num):
    return int(num + 8 ) - (int(num + 8 ) % 16)

from qiskit import pulse            # これで、Pulseのすべての機能にアクセスできます。
from qiskit.pulse import Play
from qiskit.pulse import pulse_lib  # このPulseモジュールは、一般的なパルス形状のサンプルパルスを構築するのに役立ちます


# 駆動パルスのパラメーター (us = マイクロ秒)
drive_sigma_us = 0.075                     # ガウシアンの実際の幅を決定
drive_samples_us = drive_sigma_us*8        # 切り捨てパラメーター。sigmaの８倍の長さでカットする。
                                           # ガウシアンには自然な有限長がないためです。

drive_sigma = get_closest_multiple_of_16(drive_sigma_us * us /dt)       # dt単位でのガウシアンの幅
drive_samples = get_closest_multiple_of_16(drive_samples_us * us /dt)   # dt単位での切り捨てパラメーター
drive_amp = 0.3
# 駆動パルスを構築
drive_pulse = pulse_lib.gaussian(duration=drive_samples,
                                 sigma=drive_sigma,
                                 amp=drive_amp,
                                 name='freq_sweep_excitation_pulse')

量子ビットを適切に測定するには、測定マップを確認する必要があります。ハードウェアの制約の確認です。1つの量子ビットに対して取得されると、残りの量子ビットに対しても行われます。 Pulseでプログラムをビルドするときは、この制約を尊重する必要があります。私たちの量子ビットが入っている量子ビットのグループをチェックしましょう：

# この量子ビットに必要な測定マップインデックスを見つける
meas_map_idx = None
for i, measure_group in enumerate(backend_config.meas_map):
    if qubit in measure_group:
        meas_map_idx = i
        break
assert meas_map_idx is not None, f"Couldn't find qubit {qubit} in the meas_map!"

次に、測定パルスを定義します。パルスをハードコーディングするのではなく、バックエンドのデフォルトの命令スケジュールマップからキャリブレーション済みの測定パルスを取得できます。頻繁にキャリブレーションされるため、測定パルスを自分で定義するよりも正確です。この測定パルスには取得チャネルも含まれるため、取得チャネルは手動で追加する必要はありません。

inst_sched_map = backend_defaults.instruction_schedule_map
measure = inst_sched_map.get('measure', qubits=backend_config.meas_map[meas_map_idx])

最後に、パルスを適用するチャネルを指定します。駆動、測定、取得チャネルは、量子ビットのIDによってインデックス化されます。

### 必要なチャネルを集める
drive_chan = pulse.DriveChannel(qubit)
meas_chan = pulse.MeasureChannel(qubit)
acq_chan = pulse.AcquireChannel(qubit)

パルスパラメータが定義され、実験用のパルス形状が作成されたので、パルススケジュールの作成に進むことができます。

各周波数で、その周波数の駆動パルスを量子ビットに送信し、そのパルスの直後に測定します。 パルスエンベロープは周波数に依存しないため、再利用可能なscheduleを作成し、周波数構成配列でドライブパルス周波数を指定します。

# 基本のスケジュールを作成
# 駆動パルスを駆動チャネルで動かすことから開始
schedule = pulse.Schedule(name='Frequency sweep')
schedule += Play(drive_pulse, drive_chan)
# `<<` は、スケジュールの開始時刻をある期間だけシフトすることを意味する特別な構文です
schedule += measure << schedule.duration

# スイープする周波数を設定（Hz単位で）
frequencies_Hz = frequencies_GHz*GHz
schedule_frequencies = [{drive_chan: freq} for freq in frequencies_Hz]

健全性チェックとして、常にパルススケジュールを確認することをお勧めします。 これは、以下のようにschedule.draw()を使って行われます。

schedule.draw(label=True)

上記のschedulesとschedule_frequenciesをQobjと呼ばれるプログラムオブジェクトにアセンブルし、量子デバイスに送信します。量子ビット応答の適切な推定値を得るために、各スケジュール（周波数スイープの各ポイント）をnum_shots_per_frequency回繰り返すように要求します。

測定設定も指定します。meas_level=0は生データ（ショットごとの複素数値の配列）を返し、meas_level=1はカーネルデータ（ショットごとに1つの複素数値）を返し、meas_level=2は分類されたデータ（ショットごとに0または1のビット）を返します。ラボにいて、0と1を分類するための識別器(discriminator)をまだキャリブレーションしていないと仮定して、meas_level=1を選択し、作業しているものを複製します。個々のショットではなく、結果の'avg'を求めます。

from qiskit import assemble

num_shots_per_frequency = 1024
frequency_sweep_program = assemble(schedule,
                                   backend=backend, 
                                   meas_level=1,
                                   meas_return='avg',
                                   shots=num_shots_per_frequency,
                                   schedule_los=schedule_frequencies)

/home/dsche/.local/lib/python3.8/site-packages/qiskit/compiler/assemble.py:320: RuntimeWarning: Dynamic rep rates not supported on this backend. rep_time will be used instead of rep_delay.
  warnings.warn(

別のユニット変更警告が表示される場合がありますが、これは無視してかまいません。最後に、以下を使用して、アセンブルされたプログラムをバックエンドで実行できます：

job = backend.run(frequency_sweep_program)

後で結果を取り出すためにjob_idをプリントし、job_monitor()でジョブのステータスをモニターすると良いです。

# print(job.job_id())
from qiskit.tools.monitor import job_monitor
job_monitor(job)

Job Status: job has successfully run

ジョブが実行が完了したら、以下を使って結果を取得できます：

frequency_sweep_results = job.result(timeout=120) # タイムアウトのパラメーターは120秒にセット

結果を抽出し、matplotlibを使ってプロットします：

import matplotlib.pyplot as plt

sweep_values = []
for i in range(len(frequency_sweep_results.results)):
    # i番目の実験から結果を得る
    res = frequency_sweep_results.get_memory(i)*scale_factor
    # この実験から`qubit`の結果を得る
    sweep_values.append(res[qubit])

plt.scatter(frequencies_GHz, np.real(sweep_values), color='black') # スイープ値の実数部分をプロット
plt.xlim([min(frequencies_GHz), max(frequencies_GHz)])
plt.xlabel("Frequency [GHz]")
plt.ylabel("Measured signal [a.u.]")
plt.show()

上に見られるように、中心付近のピークは量子ビット周波数の位置に対応しています。信号はパワーブロードニングを示しています。これは、中心周波数に近づくと、量子ビットを非共振に駆動できることを示しています。ピーク周波数の値を取得するために、値を共鳴応答曲線に適合させます。これは通常、ローレンツ分布の形です。

from scipy.optimize import curve_fit

def fit_function(x_values, y_values, function, init_params):
    fitparams, conv = curve_fit(function, x_values, y_values, init_params)
    y_fit = function(x_values, *fitparams)
    
    return fitparams, y_fit

fit_params, y_fit = fit_function(frequencies_GHz,
                                 np.real(sweep_values), 
                                 lambda x, A, q_freq, B, C: (A / np.pi) * (B / ((x - q_freq)**2 + B**2)) + C,
                                 [-5, 4.975, 1, 3] # フィッティングの初期パラメーター
                                )

plt.scatter(frequencies_GHz, np.real(sweep_values), color='black')
plt.plot(frequencies_GHz, y_fit, color='red')
plt.xlim([min(frequencies_GHz), max(frequencies_GHz)])

plt.xlabel("Frequency [GHz]")
plt.ylabel("Measured Signal [a.u.]")
plt.show()

A, rough_qubit_frequency, B, C = fit_params
rough_qubit_frequency = rough_qubit_frequency*GHz # 周波数の単位をHzからGHzへ
print(f"We've updated our qubit frequency estimate from "
      f"{round(backend_defaults.qubit_freq_est[qubit] / GHz, 5)} GHz to {round(rough_qubit_frequency/GHz, 5)} GHz.")

We've updated our qubit frequency estimate from 4.97445 GHz to 4.97453 GHz.

3. 較正と$\pi$パルスの利用

3.1 ラビ実験を使った$\pi$パルスの較正

量子ビットの周波数がわかったら、次のステップは、$\pi$パルスの強度を決定することです。量子ビットを2レベルのシステムとして厳密に言うと、$\pi$パルスとは、量子ビットを$\vert0\rangle$から$\vert1\rangle$に、またはその逆に変換するパルスです。これは、$X$ または$X180$ゲート、ビットフリップ演算子とも呼ばれます。前の周波数スイープ実験からこの遷移を駆動するためのマイクロ波周波数が既に分かっているため、$\vert0\rangle$から$\vert1\rangle$への$\pi$回転を実現するのに必要な振幅を求めます。望ましい回転は、下の図のブロッホ球に示されています -- $\pi$パルスは、ブロッホ球をスイープする角度から名前が付けられていることがわかります。

駆動パルスの振幅を少しずつ変化させ、量子ビットの状態を毎回測定します。量子ビットが$\vert0\rangle$から$\vert1\rangle$へ行ったり来たりすること、つまりラビ振動と一般に呼ばれれいる振動が見られることが予想されます。

# この実験では以下について前回の実験から得た値を用います：
    # `qubit`,
    # `measure`, and
    # `rough_qubit_frequency`.

# ラビ実験のパラメーター
num_rabi_points = 50

# 駆動パルスの振幅値：0から0.75の間で等間隔に配置された50個の振幅
drive_amp_min = 0
drive_amp_max = 0.75
drive_amps = np.linspace(drive_amp_min, drive_amp_max, num_rabi_points)

# ラビ実験を構築：
#    量子ビット周波数において、毎回駆動パルスの振幅を変化させながら
#    パルスを駆動し、その後測定する
rabi_schedules = []
for drive_amp in drive_amps:
    rabi_pulse = pulse_lib.gaussian(duration=drive_samples, amp=drive_amp, 
                                    sigma=drive_sigma, name=f"Rabi drive amplitude = {drive_amp}")
    this_schedule = pulse.Schedule(name=f"Rabi drive amplitude = {drive_amp}")
    this_schedule += Play(rabi_pulse, drive_chan)
    # 周波数スイープの実験で使った測定命令を再利用する
    this_schedule += measure << this_schedule.duration
    rabi_schedules.append(this_schedule)

周波数スイープ実験の時と基本的に同じスケジュールに見えます。違いは、周波数を変化させているのではなく、駆動パルスの振幅を変化させる実験を実行していることです。

rabi_schedules[-1].draw(label=True)

# Qobjにスケジュールをアセンブルします
num_shots_per_point = 1024

rabi_experiment_program = assemble(rabi_schedules,
                                   backend=backend,
                                   meas_level=1,
                                   meas_return='avg',
                                   shots=num_shots_per_point,
                                   schedule_los=[{drive_chan: rough_qubit_frequency}]
                                                * num_rabi_points)

/home/dsche/.local/lib/python3.8/site-packages/qiskit/compiler/assemble.py:320: RuntimeWarning: Dynamic rep rates not supported on this backend. rep_time will be used instead of rep_delay.
  warnings.warn(

print(job.job_id())
job = backend.run(rabi_experiment_program)
job_monitor(job)

Job Status: job has successfully run

rabi_results = job.result(timeout=120)

結果が得られたので、それらを抽出し、正弦曲線にフィットさせます。私たちが選んだ駆動振幅の範囲で、量子ビットがブロッホ球の上で$|0\rangle$から始めて完全に何回か回転すると予想されます。この正弦波の振幅は、$|1\rangle$状態を生成したラビ駆動振幅でのショットの割合を示します。信号が最大（すべて $|0\rangle$の状態）から最小（すべて$|1\rangle$の状態）に振動するために必要な駆動振幅を求めます。 -- これによって、$\pi$パルスを実行する振幅がキャルブレーションされます。

# 0を中心にデータを合わせる
def baseline_remove(values):
    return np.array(values) - np.mean(values)

rabi_values = []
for i in range(num_rabi_points):
    # i番目の実験から`qubit`の結果を得る
    rabi_values.append(rabi_results.get_memory(i)[qubit]*scale_factor)

rabi_values = np.real(baseline_remove(rabi_values))

plt.xlabel("Drive amp [a.u.]")
plt.ylabel("Measured signal [a.u.]")
plt.scatter(drive_amps, rabi_values, color='black') # ラビ実験結果の実部のみプロット
plt.show()

fit_params, y_fit = fit_function(drive_amps,
                                 rabi_values, 
                                 lambda x, A, B, drive_period, phi: (A*np.cos(2*np.pi*x/drive_period - phi) + B),
                                 [3, 0.1, 0.5, 0])

plt.scatter(drive_amps, rabi_values, color='black')
plt.plot(drive_amps, y_fit, color='red')

drive_period = fit_params[2] # 3個目のパラメーターがラビ振動の周期

plt.axvline(drive_period/2, color='red', linestyle='--')
plt.axvline(drive_period, color='red', linestyle='--')
plt.annotate("", xy=(drive_period, 0), xytext=(drive_period/2,0), arrowprops=dict(arrowstyle="<->", color='red'))
plt.annotate("$\pi$", xy=(drive_period/2-0.03, 0.1), color='red')

plt.xlabel("Drive amp [a.u.]", fontsize=15)
plt.ylabel("Measured signal [a.u.]", fontsize=15)
plt.show()

pi_amp = abs(drive_period / 2)
print(f"Pi Amplitude = {pi_amp}")

Pi Amplitude = 0.23400325324347926

我々の$\pi$パルス！

今、決定した振幅で我々のパルスを定義して、後の実験で使えるようにします。

pi_pulse = pulse_lib.gaussian(duration=drive_samples,
                              amp=pi_amp, 
                              sigma=drive_sigma,
                              name='pi_pulse')

3.2 0か1の決定

我々の$\pi$パルスがキャリブレーションされたので、適切な確率で$\vert1\rangle$の状態を作ることができます。これを使って、$\vert0\rangle$ と $\vert1\rangle$を繰り返し測定してその測定されたシグナルをプロットすることで、測定でどのように見えるかを確認できます。 これは、識別器（discriminator）を構築するために使います。識別器とは、測定されカーネル化された複素数値（meas_level=1）を受け取り、それを0または1（meas_level=2）として分類する関数のことです。

# 2つのスケジュールを作成

# 基底状態のスケジュール
gnd_schedule = pulse.Schedule(name="ground state")
gnd_schedule += measure

# 励起状態のスケジュール
exc_schedule = pulse.Schedule(name="excited state")
exc_schedule += Play(pi_pulse, drive_chan)  # 上記2Aで見つけています
exc_schedule += measure << exc_schedule.duration

gnd_schedule.draw(label=True)

exc_schedule.draw(label=True)

基底状態と励起状態の準備スケジュールを1つのQobjにアセンブルします。これらはそれぞれnum_shots回実行されます。 今回は、結果が$|0\rangle$または$|1\rangle$に分類されていないため、meas_level=1を選択します。代わりに、カーネル化されたデータが必要です：つまり、カーネル関数を通過してショットごとに1つの複素数を生成した生の取得データです。（カーネルは、生の測定データに適用される内積と考えることができます。） exc_scheduleでのみ使用されますが、両方のスケジュールに同じ周波数を渡します。

# 実行の設定
num_shots = 1024

gnd_exc_program = assemble([gnd_schedule, exc_schedule],
                           backend=backend,
                           meas_level=1,
                           meas_return='single',
                           shots=num_shots,
                           schedule_los=[{drive_chan: rough_qubit_frequency}] * 2)

# print(job.job_id())
job = backend.run(gnd_exc_program)
job_monitor(job)

Job Status: job has successfully run

gnd_exc_results = job.result(timeout=120)

これで結果が得られたので、準備した2つの母集団を単純な散布図で視覚化できます。基底状態プログラムの結果は青で、励起状態準備プログラムの結果は赤で表示されます。 注：母集団が不規則な形（ほぼ円形ではない）の場合は、ノートブックを再実行してみてください。

gnd_results = gnd_exc_results.get_memory(0)[:, qubit]*scale_factor
exc_results = gnd_exc_results.get_memory(1)[:, qubit]*scale_factor

plt.figure(figsize=[4,4])
# 全ての結果をプロット
# gnd_scheduleの結果は青でプロット
plt.scatter(np.real(gnd_results), np.imag(gnd_results), 
                s=5, cmap='viridis', c='blue', alpha=0.5, label='state_0')
# exc_scheduleの結果は赤でプロット
plt.scatter(np.real(exc_results), np.imag(exc_results), 
                s=5, cmap='viridis', c='red', alpha=0.5, label='state_1')

# 0と1の状態の結果の平均は大きなドットでプロット
mean_gnd = np.mean(gnd_results) # 実部と虚部両方の平均をとる
mean_exc = np.mean(exc_results)
plt.scatter(np.real(mean_gnd), np.imag(mean_gnd), 
            s=200, cmap='viridis', c='black',alpha=1.0, label='state_0_mean')
plt.scatter(np.real(mean_exc), np.imag(mean_exc), 
            s=200, cmap='viridis', c='black',alpha=1.0, label='state_1_mean')

plt.ylabel('I [a.u.]', fontsize=15)
plt.xlabel('Q [a.u.]', fontsize=15)
plt.title("0-1 discrimination", fontsize=15)

plt.show()

$|0\rangle$と$|1\rangle$の2つの集団が独自のクラスターを形成していることがはっきりとわかります。カーネル化された測定結果（meas_level=1から）は、これら2つのクラスターを最適に分離する識別器を適用することによって（meas_level=2に）分類されます。最適な分離は、IQ平面の直線であり、上の図で大きなドットでプロットした平均結果から等距離にあり、2つのドットを結ぶ直線に垂直です。

与えられた点が基底状態の平均に近い場合は0を返し、励起状態の平均に近い場合は1を返すことにより、迅速な分類関数を設定できます。

import math

def classify(point: complex):
    """与えられた状態を|0>か|1>に分類する"""
    def distance(a, b):
        return math.sqrt((np.real(a) - np.real(b))**2 + (np.imag(a) - np.imag(b))**2)
    return int(distance(point, mean_exc) < distance(point, mean_gnd))

3.3 反転回復法を使った$T_1$の測定

量子ビットの$T_1$時間は、量子ビットが励起状態から基底状態に減衰するのにかかる時間です。量子コンピューターで実行する意味のあるプログラムの実行時間として、重要です。

$T_1$の測定は以前の実験と同様であり、キャリブレーションした$\pi$パルスを使用します。再び、単一の駆動パルス、$\pi$パルスを適用し、次に測定パルスを適用します。ただし、今回はすぐには測定を行いません。遅延を挿入し、その遅延を実験間で変化させます。遅延時間に対して測定されたシグナルをプロットすると、量子ビットがエネルギーで緩和するにつれて、シグナルが指数関数的に減衰することが分かります。減衰時間は、量子ビットの$T_1$または緩和時間です！

# T1実験のパラメーター
time_max_us = 450
time_step_us = 6
times_us = np.arange(1, time_max_us, time_step_us)
# dtの単位へ変換
delay_times_dt = times_us * us / dt
# 以前キャリブレーションし使ったのと同じ`pi_pulse`と量子ビット周波数を使います

# 実験のスケジュールを作ります
t1_schedules = []
for delay in delay_times_dt:
    this_schedule = pulse.Schedule(name=f"T1 delay = {delay * dt/us} us")
    this_schedule += Play(pi_pulse, drive_chan)
    this_schedule |= measure << int(delay)
    t1_schedules.append(this_schedule)

T1のスケジュールも確認できます。この実験を実際に理解するには、sched_idxの値を変えて、次のセルを複数回実行して、いくつかのスケジュールを確認してみてください。sched_idxを増やすと、すこし後で測定パルスが開始するのがわかります。

sched_idx = 0
t1_schedules[sched_idx].draw(label=True)

# 実行の設定
num_shots = 256

t1_experiment = assemble(t1_schedules,
                         backend=backend, 
                         meas_level=1,
                         meas_return='avg',
                         shots=num_shots,
                         schedule_los=[{drive_chan: rough_qubit_frequency}] * len(t1_schedules))

job = backend.run(t1_experiment)
# print(job.job_id())
job_monitor(job)

Job Status: job has successfully run

t1_results = job.result(timeout=120)

t1_values = []
for i in range(len(times_us)):
    t1_values.append(t1_results.get_memory(i)[qubit]*scale_factor)
t1_values = np.real(t1_values)

plt.scatter(times_us, t1_values, color='black') 
plt.title("$T_1$ Experiment", fontsize=15)
plt.xlabel('Delay before measurement [$\mu$s]', fontsize=15)
plt.ylabel('Signal [a.u.]', fontsize=15)
plt.show()

We can then fit the data to a decaying exponential, giving us T1!

# データをフィッティング
fit_params, y_fit = fit_function(times_us, t1_values, 
            lambda x, A, C, T1: (A * np.exp(-x / T1) + C),
            [-3, 3, 100]
            )

_, _, T1 = fit_params

plt.scatter(times_us, t1_values, color='black')
plt.plot(times_us, y_fit, color='red', label=f"T1 = {T1:.2f} us")
plt.xlim(0, np.max(times_us))
plt.title("$T_1$ Experiment", fontsize=15)
plt.xlabel('Delay before measurement [$\mu$s]', fontsize=15)
plt.ylabel('Signal [a.u.]', fontsize=15)
plt.legend()
plt.show()

4. 量子ビットコヒーレンスの決定

4.1 ラムゼー実験を使った正確な量子ビット周波数の測定

次に、量子ビット周波数をより正確に決定します。ラムゼー・パルスシークエンスを使って行います。このパルスシークエンスでは、最初に$\pi/2$（"pi over two"）パルスを適用し、$\Delta t$時間待ち、また別の$\pi/2$パルスを適用します。パルスと同じ周波数で量子ビットからのシグナルを測定するので、適用したパルスと量子ビットの間の周波数における差で振動を観測する必要があります。

# ラムゼイ実験のパラメーター
time_max_us = 1.8
time_step_us = 0.025
times_us = np.arange(0.1, time_max_us, time_step_us)
# dtの単位へ変換
delay_times_dt = times_us * us / dt

# 駆動パラメーター
# pi/2の駆動振幅はpiパルスの振幅の半分です
drive_amp = pi_amp / 2
# x_90はpi_over_2の簡潔な方法です：つまりX軸周りの90度回転
x90_pulse = pulse_lib.gaussian(duration=drive_samples,
                               amp=drive_amp, 
                               sigma=drive_sigma,
                               name='x90_pulse')

# ラムゼー実験のスケジュールを作成
ramsey_schedules = []

for delay in delay_times_dt:
    this_schedule = pulse.Schedule(name=f"Ramsey delay = {delay * dt / us} us")
    this_schedule |= Play(x90_pulse, drive_chan)
    this_schedule |= Play(x90_pulse, drive_chan) << int(this_schedule.duration + delay)
    this_schedule |= measure << int(this_schedule.duration)

    ramsey_schedules.append(this_schedule)

$T_1$スケジュールと同じように、作ったスケジュールのいくつかを検査するために次のセルを何回か実行して明らかにできます。ramsey_schedulesが増えると、２つの$\pi/2$パルスの間の遅延が増加します。

ramsey_schedules[0].draw(label=True)

ここで、一般に使用される実験的なトリックを適用します。既知の量だけパルスをオフレゾナンス（非共振）に駆動します。これをdetuning_MHzと呼びます。測定されたラムゼー・シグナルは、detuning_MHz近くの周波数で少しのオフセットで共振を示すはずです。この小さなオフセットが、rough_qubit_frequencyが量子ビット周波数からどれだけ離れているかを正確に示しています。

# 実行設定
num_shots = 256

detuning_MHz = 2 
ramsey_frequency = round(rough_qubit_frequency + detuning_MHz * MHz, 6) # need ramsey freq in Hz
ramsey_program = assemble(ramsey_schedules,
                             backend=backend,
                             meas_level=1,
                             meas_return='avg',
                             shots=num_shots,
                             schedule_los=[{drive_chan: ramsey_frequency}]*len(ramsey_schedules)
                            )

job = backend.run(ramsey_program)
# print(job.job_id())
job_monitor(job)

Job Status: job has successfully run

ramsey_results = job.result(timeout=120)

ramsey_values = []
for i in range(len(times_us)):
    ramsey_values.append(ramsey_results.get_memory(i)[qubit]*scale_factor)
    
plt.scatter(times_us, np.real(ramsey_values), color='black')
plt.xlim(0, np.max(times_us))
plt.title("Ramsey Experiment", fontsize=15)
plt.xlabel('Delay between X90 pulses [$\mu$s]', fontsize=15)
plt.ylabel('Measured Signal [a.u.]', fontsize=15)
plt.show()

データを正弦波に当てはめ、我々が関心のある情報-- つまり$\Delta f$ -- を抽出します。

fit_params, y_fit = fit_function(times_us, np.real(ramsey_values),
                                 lambda x, A, del_f_MHz, C, B: (
                                          A * np.cos(2*np.pi*del_f_MHz*x - C) + B
                                         ),
                                 [5, 1./0.4, 0, 0.25]
                                )

# オフレゾナス コンポーネント
_, del_f_MHz, _, _, = fit_params # 時間がusなので周波数はMHz

plt.scatter(times_us, np.real(ramsey_values), color='black')
plt.plot(times_us, y_fit, color='red', label=f"df = {del_f_MHz:.2f} MHz")
plt.xlim(0, np.max(times_us))
plt.xlabel('Delay between X90 pulses [$\mu$s]', fontsize=15)
plt.ylabel('Measured Signal [a.u.]', fontsize=15)
plt.title('Ramsey Experiment', fontsize=15)
plt.legend()
plt.show()

これで、del_f_MHzが分かったので、量子ビット周波数の推定値をアップデートできます。

precise_qubit_freq = rough_qubit_frequency + (del_f_MHz - detuning_MHz) * MHz # get new freq in Hz
print(f"Our updated qubit frequency is now {round(precise_qubit_freq/GHz, 6)} GHz. "
      f"It used to be {round(rough_qubit_frequency / GHz, 6)} GHz")

Our updated qubit frequency is now 4.974591 GHz. It used to be 4.974528 GHz

4.2 ハーン・エコーを使った$T_2$の測定

次に、我々の量子ビットのコヒーレンスタイム$T_2$を測定できます。この実験を行うために使われるパルスシークエンスは、ハーン・エコーとして知られていて、NMRコミュニティーに由来する用語です。ハーン・エコーの実験は上記のラムゼー実験にとてもよく似ていて、２つの$\pi/2$パルスの間に$\pi$パルスを追加します。時間$\tau$における$\pi$パルスは位相の累積を反転させ、時間$2\tau$にエコーを生成します。時間$2\tau$に最後の$\pi/2$パルスを適用して測定を行います。

ハーン・エコー実験における遅延時間がヒコーレンスタイム$T_2$です。

# T2実験のパラメーター
tau_max_us = 200
tau_step_us = 4
taus_us = np.arange(2, tau_max_us, tau_step_us)
# dtの単位へ変換
delay_times_dt = taus_us * us / dt

# 以前の実験から得たpi_pulseとx90_pulse を使います

t2_schedules = []
for tau in delay_times_dt:
    this_schedule = pulse.Schedule(name=f"T2 delay = {tau *dt/us} us")
    this_schedule |= Play(x90_pulse, drive_chan)
    this_schedule |= Play(pi_pulse, drive_chan) << int(this_schedule.duration + tau)
    this_schedule |= Play(x90_pulse, drive_chan) << int(this_schedule.duration + tau)
    this_schedule |= measure << int(this_schedule.duration)
    
    t2_schedules.append(this_schedule)

t2_schedules[0].draw(label=True)

# 実行の設定
num_shots_per_point = 512

t2_experiment = assemble(t2_schedules,
                         backend=backend,
                         meas_level=1,
                         meas_return='avg',
                         shots=num_shots_per_point,
                         schedule_los=[{drive_chan: precise_qubit_freq}]
                                      * len(t2_schedules))

/home/dsche/.local/lib/python3.8/site-packages/qiskit/compiler/assemble.py:320: RuntimeWarning: Dynamic rep rates not supported on this backend. rep_time will be used instead of rep_delay.
  warnings.warn(

job = backend.run(t2_experiment)
# print(job.job_id())
job_monitor(job)

Job Status: job has successfully run

t2_results = job.result(timeout=120)

t2_values = []
for i in range(len(taus_us)):
    t2_values.append(t2_results.get_memory(i)[qubit]*scale_factor)

plt.scatter(2*taus_us, np.real(t2_values), color='black')
plt.xlabel('Delay between X90 pulse and $\pi$ pulse [$\mu$s]', fontsize=15)
plt.ylabel('Measured Signal [a.u.]', fontsize=15)
plt.title('Hahn Echo Experiment', fontsize=15)
plt.show()

fit_params, y_fit = fit_function(2*taus_us, np.real(t2_values),
             lambda x, A, B, T2: (A * np.exp(-x / T2) + B),
             [-3, 0, 100])

_, _, T2 = fit_params
print()

plt.scatter(2*taus_us, np.real(t2_values), color='black')
plt.plot(2*taus_us, y_fit, color='red', label=f"T2 = {T2:.2f} us")
plt.xlim(0, np.max(2*taus_us))
plt.xlabel('Delay between X90 pulse and $\pi$ pulse [$\mu$s]', fontsize=15)
plt.ylabel('Measured Signal [a.u.]', fontsize=15)
plt.title('Hahn Echo Experiment', fontsize=15)
plt.legend()
plt.show()

C. 動的デカップリング

単一の$\pi$パルスは位相の累積の反転による準静的ノイズを除去できます。このコンセプトは、$\pi$パルスを何回も続けて適用することで準静的とは近似できないノイズに拡張することです。このテクニックは、動的デッカプリングとして一般に知られており、ノイズのさまざまな周波数をキャンセルし、量子ビットからより長いコヒーレンス時間を抽出するために使われます。

# 動的デカップリング実験のパラメーター
tau_us_min = 1
tau_us_max = 40
tau_step_us = 1.5
taus_us = np.arange(tau_us_min, tau_us_max, tau_step_us)
# dtの単位への変換
taus_dt = taus_us * us / dt
num_pi_pulses = 6 # 6つのpiパルスを適用
print(f"Total time ranges from {2.*num_pi_pulses*taus_us[0]} to {2.*num_pi_pulses*taus_us[-1]} us")

Total time ranges from 12.0 to 462.0 us

T2DD_schedules = []
for delay in taus_dt:
    this_schedule = pulse.Schedule(name=f"T2DD delay = {delay * dt/us} us")
    this_schedule |= Play(x90_pulse, drive_chan)
    this_schedule |= Play(pi_pulse, drive_chan) << int(this_schedule.duration + delay)

    for _ in range(num_pi_pulses - 1):
        this_schedule |= Play(pi_pulse, drive_chan) << int(this_schedule.duration + 2*delay)

    this_schedule |= Play(x90_pulse, drive_chan) << int(this_schedule.duration + delay)
    this_schedule |= measure << int(this_schedule.duration)
    
    T2DD_schedules.append(this_schedule)

T2DD_schedules[0].draw(label=True)

num_shots_per_point = 1024

T2DD_experiment = assemble(T2DD_schedules,
                             backend=backend,
                             meas_level=1,
                             meas_return='avg',
                             shots=num_shots_per_point,
                             schedule_los=[{drive_chan: precise_qubit_freq}]
                                          * len(T2DD_schedules))

job = backend.run(T2DD_experiment)
print(job.job_id())
job_monitor(job)

5f85cfa300f7d2001a1f7466
Job Status: job has successfully run

T2DD_results = job.result(timeout=120)

times_us = 2.*num_pi_pulses*taus_us
DD_values = []
for i in range(len(taus_us)):
    DD_values.append(T2DD_results.get_memory(i)[qubit]*scale_factor)

plt.scatter(times_us, np.real(DD_values), color='black')
plt.xlim(0, np.max(times_us))
plt.xlabel('Total time before measurement [$\mu$s]', fontsize=15)
plt.ylabel('Measured Signal [a.u.]', fontsize=15)
plt.title('Dynamical Decoupling Experiment', fontsize=15)
plt.show()

# データをフィッテイング
fit_func = lambda x, A, B, T2DD: (A * np.exp(-x / T2DD) + B)
fitparams, conv = curve_fit(fit_func, times_us, np.real(DD_values), [3.5, 0.8, 150])

_, _, T2DD = fitparams
plt.scatter(times_us, np.real(DD_values), color='black')
plt.plot(times_us, fit_func(times_us, *fitparams), color='red', label=f"T2DD = {T2DD:.2f} us")
plt.xlim([0, np.max(times_us)])
plt.xlabel('Total time before measurement [$\mu$s]', fontsize=15)
plt.ylabel('Measured Signal [a.u.]', fontsize=15)
plt.title('Dynamical Decoupling Experiment', fontsize=15)
plt.legend()
plt.show()

5. References

1. H. Abraham, I. Y. Akhalwaya, G. Aleksandrowicz, T. Alexander, G. Alexandrowics, E. Arbel, A. Asfaw, C. Azaustre, P. Barkoutsos, G. Barron, L. Bello, Y. Ben-Haim, L. S. Bishop, S. Bosch, D. Bucher, CZ, F. Cabrera, P. Calpin, L. Capelluto, J. Carballo, C.-F. Chen, A. Chen, R. Chen, J. M. Chow, C. Claus, A. W. Cross, A. J. Cross, J. Cruz- Benito, C. Culver, A. D. C ́orcoles-Gonzales, S. Dague, M. Dartiailh, A. R. Davila, D. Ding, E. Dumitrescu, K. Dumon, I. Duran, P. Eendebak, D. Egger, M. Everitt, P. M. Fern ́andez, A. Frisch, A. Fuhrer, J. Gacon, Gadi, B. G. Gago, J. M. Gambetta, L. Garcia, S. Garion, Gawel-Kus, L. Gil, J. Gomez-Mosquera, S. de la Puente Gonz ́alez, D. Green- berg,J.A.Gunnels,I.Haide,I.Hamamura,V.Havlicek,J.Hellmers,L􏰀.Herok,H.Horii, C. Howington, W. Hu, S. Hu, H. Imai, T. Imamichi, R. Iten, T. Itoko, A. Javadi-Abhari, Jessica, K. Johns, N. Kanazawa, A. Karazeev, P. Kassebaum, V. Krishnan, K. Kr- sulich, G. Kus, R. LaRose, R. Lambert, J. Latone, S. Lawrence, P. Liu, P. B. Z. Mac, Y. Maeng, A. Malyshev, J. Marecek, M. Marques, D. Mathews, A. Matsuo, D. T. Mc- Clure, C. McGarry, D. McKay, S. Meesala, A. Mezzacapo, R. Midha, Z. Minev, P. Mu- rali, J. Mu ̈ggenburg, D. Nadlinger, G. Nannicini, P. Nation, Y. Naveh, Nick-Singstock, P. Niroula, H. Norlen, L. J. O’Riordan, S. Oud, D. Padilha, H. Paik, S. Perriello, A. Phan, M. Pistoia, A. Pozas-iKerstjens, V. Prutyanov, J. P ́erez, Quintiii, R. Raymond, R. M.-C. Redondo, M. Reuter, D. M. Rodr ́ıguez, M. Ryu, M. Sandberg, N. Sathaye, B. Schmitt, C. Schnabel, T. L. Scholten, E. Schoute, I. F. Sertage, Y. Shi, A. Silva, Y. Siraichi, S. Sivarajah, J. A. Smolin, M. Soeken, D. Steenken, M. Stypulkoski, H. Takahashi, C. Taylor, P. Taylour, S. Thomas, M. Tillet, M. Tod, E. de la Torre, K. Trabing, M. Treinish, TrishaPe, W. Turner, Y. Vaknin, C. R. Valcarce, F. Varchon, D. Vogt- Lee, C. Vuillot, J. Weaver, R. Wieczorek, J. A. Wildstrom, R. Wille, E. Winston, J. J. Woehr, S. Woerner, R. Woo, C. J. Wood, R. Wood, S. Wood, J. Wootton, D. Yeralin, J. Yu, L. Zdanski, Zoufalc, azulehner, drholmie, fanizzamarco, kanejess, klinvill, merav aharoni, ordmoj, tigerjack, yang.luh, and yotamvakninibm, “Qiskit: An open-source framework for quantum computing,” 2019.

2. D. C. McKay, T. Alexander, L. Bello, M. J. Biercuk, L. Bishop, J. Chen, J. M. Chow, A. D. C ́orcoles, D. Egger, S. Filipp, J. Gomez, M. Hush, A. Javadi-Abhari, D. Moreda, P. Nation, B. Paulovicks, E. Winston, C. J. Wood, J. Wootton, and J. M. Gambetta, “Qiskit backend specifications for OpenQASM and OpenPulse experiments,” 2018.

Note: 'Qiskit Pulse' was formerly known as 'OpenPulse'.

import qiskit
qiskit.__qiskit_version__

{'qiskit-terra': '0.15.1',
 'qiskit-aer': '0.6.1',
 'qiskit-ignis': '0.4.0',
 'qiskit-ibmq-provider': '0.8.0',
 'qiskit-aqua': '0.7.5',
 'qiskit': '0.20.0'}