【課題】量子計算

\(\newcommand{\ket}[1]{|#1\rangle}\)

問題1: 足し算を実機で行う

実習の内容の延長です。ここでは並列足し算回路を実機で実行します。

実習のおさらいをすると、まずもともとの足し算のロジックをそのまま踏襲した回路を作り、それではゲート数が多すぎるので効率化した回路を作成しました。

実は効率化した回路でもまだゲートの数が多すぎて、4ビット+4ビットの計算では答えがスクランブルされてしまいます。回路が小規模になればそれだけ成功確率も上がりますので、\((n_1, n_2)\)の値として(4, 4)以外に(3, 3)、(2, 2)、(1, 1)も同時に試すことにしましょう。

# まずは全てインポート
import numpy as np
import matplotlib.pyplot as plt
from qiskit import QuantumRegister, ClassicalRegister, QuantumCircuit, IBMQ, Aer, transpile, execute
from qiskit.providers.ibmq import least_busy
from qiskit.visualization import plot_histogram
from qiskit.tools.monitor import job_monitor
from utils.optimized_additions import optimized_additions, get_initial_layout

print('notebook ready')
IBMQ.enable_account('__paste_your_token_here__')

provider = IBMQ.get_provider(hub='ibm-q', group='open', project='main')

実習を参考に、次のセルで効率化前の回路を返すmake_original_circuit関数を実装してください。

def make_original_circuit(n1, n2):
    """A function to define a circuit with the original implementation of additions given n1 and n2
    """
    n3 = np.ceil(np.log2((2 ** n1) + (2 ** n2) - 1)).astype(int)

    reg1 = QuantumRegister(n1, 'r1')
    reg2 = QuantumRegister(n2, 'r2')
    reg3 = QuantumRegister(n3, 'r3')
    creg = ClassicalRegister(n1 + n2 + n3)

    # QuantumCircuit can be instantiated from multiple registers
    circuit = QuantumCircuit(reg1, reg2, reg3, creg)
    
    ##################
    ### EDIT BELOW ###
    ##################

    #circuit.?

    ##################
    ### EDIT ABOVE ###
    ##################

    return circuit


# 4 + 4 + 5 = 13量子ビット以上のマシンで ibm-q/open/main からアクセス可能なのはibmq_16_melbourneのみ(2021年3月現在)
backend = provider.get_backend('ibmq_16_melbourne')

# List of circuits
circuits = []

for n1, n2 in [(4, 4), (3, 3), (2, 2), (1, 1)]:
    print('Original circuit with n1, n2 = {}, {}'.format(n1, n2))
    circuit = make_original_circuit(n1, n2)
    
    print('  Transpiling..')
    circuit = transpile(circuit, backend=backend, optimization_level=3)
    print('  Done. Ops: N(Rz)={rz}, N(SX)={sx}, N(CNOT)={cx}'.format(**circuit.count_ops()))
    circuit.name = 'original_{}_{}'.format(n1, n2)
    circuits.append(circuit)

    print('Optimized circuit with n1, n2 = {}, {}'.format(n1, n2))    
    circuit = optimized_additions(n1, n2)

    print('  Transpiling..')
    initial_layout = get_initial_layout(backend, n1, n2)
    circuit = transpile(circuit, backend=backend, routing_method='basic', initial_layout=initial_layout, optimization_level=3)
    print('  Done. Ops: N(Rz)={rz}, N(SX)={sx}, N(CNOT)={cx}'.format(**circuit.count_ops()))
    circuit.name = 'optimized_{}_{}'.format(n1, n2)
    circuits.append(circuit)

shots = 8192    
    
print('Submitting {} circuits to {}'.format(len(circuits), backend.name()))
job = backend.run(circuits, shots=shots)

job_monitor(job, interval=2)

counts_list = job.result().get_counts()

ジョブが返ってきたら、8192ショットのうち正しい足し算を表しているものがどれだけあるか調べてみましょう。

def count_correct_additions(counts_list, n1, n2, shots):
    for counts, ctype in zip(counts_list, ['Original ', 'Optimized']):
        correct_additions = 0

        for key, value in counts.items():
            x1 = int(key[-n1:], 2)
            x2 = int(key[-n1 - n2:-n1], 2)
            x3 = int(key[:-n1 - n2], 2)

            if x1 + x2 == x3:
                correct_additions += value
            
        print('{} circuit ({}, {}): {} / 8192 = {}'.format(ctype, n1, n2, correct_additions, correct_additions / shots))

icirc = 0
for n1, n2 in [(4, 4), (3, 3), (2, 2), (1, 1)]:
    count_correct_additions(counts_list[icirc:icirc + 2], n1, n2, shots)
    icirc += 2

回路が均一にランダムに\(0\)から\(2^{n_1 + n_2 + n_3} - 1\)までの数を返す場合、レジスタ1と2のそれぞれの値の組み合わせに対して正しいレジスタ3の値が一つあるので、正答率は\(2^{n_1 + n_2} / 2^{n_1 + n_2 + n_3} = 2^{-n_3}\)となります。実機の結果と比べてみると、(4, 4), (3, 3), (2, 2)では残念ながら元の回路も効率化した回路も正答率がほとんどこの値に近くなっているのがわかります。かろうじて(1, 1)において正答率がランダムなケース(\(2^{-2}\))より高くなっています。

Quantum Volume

実は上のテストをIBMの提供するより高性能なマシンで実行すると、元の回路と効率化した回路とで正答率が異なることが見えます。\((n_1, n_2)\) = (1, 1)に関しては効率化回路のほうがCNOTが少しだけ多くなるケースがあり、その場合実際に元の回路の正答率のほうが高くなります。例としてibmq_torontoというマシンで実行した結果は以下のとおりです。

Original circuit with n1, n2 = 4, 4
  Transpiling..
  Done. Ops: N(Rz)=370, N(SX)=313, N(CNOT)=290
Optimized circuit with n1, n2 = 4, 4
  Transpiling..
  Done. Ops: N(Rz)=142, N(SX)=71, N(CNOT)=142
Original circuit with n1, n2 = 3, 3
  Transpiling..
  Done. Ops: N(Rz)=203, N(SX)=164, N(CNOT)=183
Optimized circuit with n1, n2 = 3, 3
  Transpiling..
  Done. Ops: N(Rz)=95, N(SX)=53, N(CNOT)=84
Original circuit with n1, n2 = 2, 2
  Transpiling..
  Done. Ops: N(Rz)=104, N(SX)=74, N(CNOT)=81
Optimized circuit with n1, n2 = 2, 2
  Transpiling..
  Done. Ops: N(Rz)=60, N(SX)=35, N(CNOT)=41
Original circuit with n1, n2 = 1, 1
  Transpiling..
  Done. Ops: N(Rz)=36, N(SX)=30, N(CNOT)=12
Optimized circuit with n1, n2 = 1, 1
  Transpiling..
  Done. Ops: N(Rz)=25, N(SX)=14, N(CNOT)=13
Original  circuit (4, 4): 289 / 8192 = 0.0352783203125
Optimized circuit (4, 4): 340 / 8192 = 0.04150390625
Original  circuit (3, 3): 531 / 8192 = 0.0648193359375
Optimized circuit (3, 3): 960 / 8192 = 0.1171875
Original  circuit (2, 2): 1250 / 8192 = 0.152587890625
Optimized circuit (2, 2): 1984 / 8192 = 0.2421875
Original  circuit (1, 1): 6114 / 8192 = 0.746337890625
Optimized circuit (1, 1): 5288 / 8192 = 0.6455078125

効率化した回路では(3, 3)の段階(回路3)で、すでに正答率が\(2^{-4}\)を大きく上回っているのがわかります。

IBMQでは量子コンピュータ一つ一つをQuantum Volume(QV、量子体積)[CBS+19]という指標で評価しています1。QVは簡単に言えば「量子コンピュータ上である特定の形を持った回路を安定的に実行できる量子ビット数と回路の長さ」を測っていて、QVの値が大きいマシンほど高性能と言えます。QVは全体の量子ビット数だけではなく、量子ビット間の接続度合、(特に制御ゲートの)エラー率、トランスパイラの性能などを総合的に評価した指標です。

2021年3月現在、IBMQのマシンでQVの最大値は128ですが、ibmq_16_melbourneのQVは8です。QVの値は2の実効量子ビット数乗と考えればよく、これはつまり(かなり乱暴に言うと)15ビットマシンibmq_16_melbourneの上で信頼度高く実行できるのは3量子ビット回路までだということを意味します。

上で利用したibmq_torontoのQVは32でした。量子ビット数が少ないマシンであればopenのアカウントでもQV=32のマシンが複数存在するので、そのどれかを使って\((n_1, n_2)\) = (1, 1)の足し算を走らせてみましょう。

backend_filter = lambda b: (not b.configuration().simulator) and b.configuration().n_qubits >= 4 and b.configuration().quantum_volume >= 32 and b.status().operational
backend = least_busy(provider.backends(filters=backend_filter))

print('Using backend', backend.name())
n1 = n2 = 1

circuits = []

print('Original circuit with n1, n2 = {}, {}'.format(n1, n2))
circuit = make_original_circuit(n1, n2)
print('  Transpiling..')
circuit = transpile(circuit, backend=backend, optimization_level=3)
print('  Done. Ops: N(Rz)={rz}, N(SX)={sx}, N(CNOT)={cx}'.format(**circuit.count_ops()))
circuit.name = 'original_{}_{}'.format(n1, n2)
circuits.append(circuit)

print('Optimized circuit with n1, n2 = {}, {}'.format(n1, n2))    
circuit = optimized_additions(n1, n2)
print('  Transpiling..')
initial_layout = get_initial_layout(backend, n1, n2)
circuit = transpile(circuit, backend=backend, routing_method='basic', initial_layout=initial_layout, optimization_level=3)
print('  Done. Ops: N(Rz)={rz}, N(SX)={sx}, N(CNOT)={cx}'.format(**circuit.count_ops()))
circuit.name = 'optimized_{}_{}'.format(n1, n2)
circuits.append(circuit)

print('Submitting circuits')
job = backend.run(circuits, shots=shots)

job_monitor(job, interval=2)

counts_list = job.result().get_counts()

count_correct_additions(counts_list, n1, n2, shots)

ibmq_16_melbourneを使った場合と比べて、正答率はどうでしょうか?

完全にランダム(正答率0.25)よりはいい結果が得られるはずなので、公開されているバックエンドのエラー率から正答率を予想し、結果と照らし合わせてみましょう。

IBMQの各バックエンドに対して、次の3種類のエラー率が公開されています。

  • 1量子ビットゲート(\(X\)とSX)の量子ビットごとのエラー率

  • CNOTの量子ビット間接続ごとのエラー率

  • 測定の量子ビットごとのエラー率

これまで触れて来ませんでしたが、測定自体でもエラーが起こる余地があります。測定のエラーとは、具体的には状態\(\ket{0}\)\(\ket{1}\)を測定したときにどのくらいの割合で0や1でない結果を得るかということを表した値です。

バックエンドのエラー率はbackend.properties()メソッドから取り出します。以下のコードを実行して、get_initial_layout関数から得られるマッピングのリストに含まれる量子ビットとそれらの間の接続について、SX、CNOT、測定のエラー率を取得しましょう。2021年3月現在、公表されている\(X\)とSXのエラー率は同じなので(IBMが別々の測定をしていません)、SXのエラー率を持って1量子ビットゲート全般のエラー率と考えます。

def get_error_rates(backend, qubits):
    links = set()
    for qubit1 in qubits:
        for qubit2 in qubits:
            if qubit1 != qubit2:
                links.add((qubit1, qubit2))
                
    qubits = set(qubits)

    properties = backend.properties()

    sx_error = dict()
    cx_error = dict()
    meas_error = dict()

    for gate in properties.gates:
        if gate.gate == 'sx':
            qubit = gate.qubits[0]
            if qubit not in qubits:
                continue
            
            for param in gate.parameters:
                if param.name == 'gate_error':
                    sx_error[qubit] = param.value
                
        elif gate.gate == 'cx':
            link = tuple(gate.qubits)
            if link not in links:
                continue

            for param in gate.parameters:
                if param.name == 'gate_error':
                    cx_error[link] = param.value

    for qubit in qubits:
        qprops = properties.qubits[qubit]

        for prop in qprops:
            if prop.name == 'readout_error':
                meas_error[qubit] = prop.value
                
    return sx_error, cx_error, meas_error


qubits = get_initial_layout(backend, n1, n2)
print('Qubits:', qubits)
sx_error, cx_error, meas_error = get_error_rates(backend, qubits)
print('Single-qubit gate error rates:')
for qubit in sorted(sx_error):
    print(' {}: {:.2e}'.format(qubit, sx_error[qubit]))
print('CNOT error rates:')
for link in sorted(cx_error):
    print(' {}: {:.2e}'.format(link, cx_error[link]))
print('Measurement error rates')
for qubit in sorted(meas_error):
    print(' {}: {:.2e}'.format(qubit, meas_error[qubit]))

CNOTのエラー率がSXのそれよりも一桁高く、さらに測定のエラー率がCNOTのそれの数倍になっているのがわかります。また、量子ビットごとや接続ごとにエラー率が大きく異なることも見て取れます。

それでは、これらの数字を元に、足し算回路の正答率を予想してみましょう。単純化するために以下の仮定をおきます。

  • 1量子ビットゲートのエラーは無視できる

  • CNOTのエラー率については平均値を取り、回路中すべてのCNOTゲートが同じ確率\(\epsilon\)でエラーを起こす

  • 回路中2箇所以上でエラーが起きた場合は結果がランダムになり、確率\(2^{-n_3}\)で正しい答えに帰着する

すると、\(n\)ビット回路にCNOTが\(t\)個あり、ビット\(j\)の測定エラー率が\(\delta_j\)のとき、エラーが一切起こらない確率が

\[ P_0 = (1 - \epsilon)^t \prod_{j=0}^{n-1} (1 - \delta_j), \]

エラーが一度だけ起こる確率が

\[ P_1 = t \epsilon (1 - \epsilon)^{t-1} \prod_{j=0}^{n-1} (1 - \delta_j) + (1 - \epsilon)^t \sum_{j=0}^{n-1} \delta_j \prod_{k \neq j} (1 - \delta_k) \]

です。エラーが二度以上起こる確率は\(1 - P_0 - P_1\)です。エラーが一度だけ起きた場合は絶対に正しい答えが得られないことを考えると、正答率\(F\)

\[ F = P_0 + 2^{-n_3} (1 - P_0 - P_1) \]

です。

それでは、cx_errormeas_error、 CNOTの数n_cxと、randomized_rate = \(2^{-n_3}\)から上の計算で\(F\)を返す関数を書いて、結果を実機でみられた正答率と比較してください。

def predict_correct_fraction(cx_error, meas_error, n_cx, randomized_rate):
    ##################
    ### EDIT BELOW ###
    ##################

    ##################
    ### EDIT ABOVE ###
    ##################
    pass

print(predict_correct_fraction(cx_error, meas_error, 12, np.power(2., -2)))

提出するもの

  • make_original_circuit関数のコード

  • 実機で各回路を実行した結果

  • predict_correct_fraction関数のコードと得られる正答率の予想

問題2: 符号が反転している基底を見つける

実習の前半で出てきたequal superposition状態

\[ H^{\otimes n} \ket{0} = \frac{1}{\sqrt{2^n}} \sum_{k=0}^{2^n-1} \ket{k} \]

をそのまま測定すると、全ての整数\(k\)に対応するビット列が等しい確率で現れます。測定でビット列が現れる確率はそのビット列に対応する計算基底の振幅の絶対値自乗で決まるので、重ね合わせにおいてある整数\(\tilde{k}\)の符号だけ逆転している以下の状態でもやはり全ての整数が確率\(1/2^n\)で得られます。

\[ \frac{1}{\sqrt{2^n}} \left( \sum_{k \neq \tilde{k}} \ket{k} - \ket{\tilde{k}} \right) \]

(一般には、全ての計算基底にバラバラに位相因子\(e^{i\theta_{k}}\)がかかっていても確率は同じです。)

さて、後の実習で登場するグローバー探索というアルゴリズムは、上のように一つの計算基底の符号を逆転させるブラックボックス演算子(どの基底かは事前に知られていない)が与えられたときに、符号の反転が起こっている計算基底を効率よく見つけ出すための手法です。グローバー探索を利用すると、例えば\(N\)個のエントリーのあるデータベースから特定のエントリーを探し出すのに、\(\mathcal{O}(\sqrt{N})\)回データベースを参照すればいいということがわかっています。

今から考えるのはそのような効率的な方法ではなく、同じようにブラックボックス演算子が与えられたときに、原理的には符号の反転が起こっている基底を見つけることができる、という手法です。そのために振幅の干渉を利用します。

まずは具体性のために\(n=3\)として、ブラックボックスは\(k=5\)の符号を反転させるとします。ここでブラックボックスの中身が完全に明かされてしまっていますが、これは実装上の都合で、重要なのは検索アルゴリズムが中身(5)を一切参照しないということです。

後で便利なように、まずはブラックボックスを単体の回路として定義します。

num_qubits = 3
needle = 5

haystack_register = QuantumRegister(num_qubits, name='haystack') # ビット数を指定してレジスタを作る
blackbox_circuit = QuantumCircuit(haystack_register, name='blackbox') # これまでの例と異なり、ビット数ではなくレジスタを指定して回路を作る

for i in range(num_qubits):
    if ((needle >> i) & 1) == 0:
        blackbox_circuit.x(haystack_register[i]) # Xゲートをレジスタのi番目のビットにかける(これまでの例との引数の違いに注意)
        
# レジスタの(0番から)最後から二番目のビットまでで制御し、最後のビットを標的にする
blackbox_circuit.mcp(np.pi, haystack_register[:-1], haystack_register[-1])

# 後片付け
for i in range(num_qubits):
    if ((needle >> i) & 1) == 0:
        blackbox_circuit.x(haystack_register[i])
        
blackbox_circuit.draw('mpl')
_images/quantum_computation_22_0.png

ここまでは単純な量子回路をゼロから書くの問題5と同じです。

Qiskitでは、QuantumCircuitオブジェクト全体を一つのゲートのようにみなして、それから制御ゲートを派生させることができます。

blackbox = blackbox_circuit.to_gate()
cblackbox = blackbox.control(1)

新たに定義したcblackboxゲートは、\(n+1\)ビット回路に

circuit.append(cblackbox, qargs=range(circuit.num_qubits))

のようにして組み込むことができます。この場合、0番ビットが制御ビットになります。

それでは、この制御ブラックボックスゲートを利用して、equal superpositionにあるhaystackレジスタで干渉を起こして、観測でneedleが識別できるような回路を書いてください。

test_register = QuantumRegister(1, 'test')
test = test_register[0]
circuit = QuantumCircuit(test_register, haystack_register)

# equal superpositionを作る(このようにゲート操作のメソッドにレジスタを渡すと、レジスタの各ビットにゲートがかかります。)
circuit.h(haystack_register)

# 干渉を起こす方法のヒント
circuit.h(test)

##################
### EDIT BELOW ###
##################

#circuit.?

##################
### EDIT ABOVE ###
##################

circuit.measure_all()

circuit.draw('mpl')

回路が完成したら、qasm_simulatorで実行し、ヒストグラムをプロットしてください。

simulator = Aer.get_backend('qasm_simulator')
sim_job = execute(circuit, backend=simulator, shots=10000)
sim_result = sim_job.result()
plot_histogram(sim_result.get_counts(circuit), figsize=(16, 4))
_images/quantum_computation_28_0.png

実機ではどんな結果が得られるでしょうか。

backend_filter = lambda b: (not b.configuration().simulator) and (b.configuration().n_qubits >= circuit.num_qubits) and b.status().operational
backend = least_busy(provider.backends(filters=backend_filter))

print('Jobs will run on', backend.name())

circuit_tr = transpile(circuit, backend=backend)
job = backend.run(circuit_tr, shots=8192)

job_monitor(job, interval=2)
result = job.result()
plot_histogram(result.get_counts(circuit), figsize=(16, 4))

シミュレータでの結果とずれていたら、トランスパイルした回路circuit_trを見て原因を考えてみてください。作ったのは一見単純な回路ですが、実際に実行されるものはどうなっているでしょうか。

circuit_tr.draw('mpl')

提出するもの

  • 完成した回路のコードと実行結果(実機とシミュレータ)のヒストグラム

  • ヒストグラムから何が読み取れるかの考察

  • haystackレジスタが一般の\(n\)ビットであるとき、この方法でneedleを探すことの問題点(実行時間の観点から)に関する考察

参考文献

CBS+19

Andrew W. Cross, Lev S. Bishop, Sarah Sheldon, Paul D. Nation, and Jay M. Gambetta. Validating quantum computers using randomized model circuits. Phys. Rev. A, 100:032328, Sep 2019. URL: https://link.aps.org/doi/10.1103/PhysRevA.100.032328, doi:10.1103/PhysRevA.100.032328.

1

QVはハードウェアの詳細に依存しないように定義されているので、量子ビット型の量子コンピュータであればIBMのマシンに限らずすべてQVで評価できます。実際、業界で徐々にQVを標準ベンチマークとして使う動きが広がってきているようです。