【最新2024年版】
HHLアルゴリズムで量子コンピュータを使って
とりあえず線型方程式を解いてみた(コード掲載)

IT Quantum Computer

こんにちは、てつをです。今回は、量子コンピュータでHHLアルゴリズムを使って線型方程式を解くまでの流れを、実際のコードと合わせて説明します。私が線型方程式を解く際、様々なサイトを参考にしたのですがどれも情報が古いせいかコードをコピーしても動かせずに苦戦したので、今回は最新版のコードを皆様に紹介したいと思います。

早速HHLアルゴリズムを使って線型方程式を解く前に、少しだけ私がどのような背景知識を持っているかだけ紹介しておきたいと思います。

私が線型方程式を解く前に持っていた知識・スキル

  • 量子コンピュータの本当に基礎の中の基礎の知識
    いちばんやさしい量子コンピュータの教本という本に書いているレベルで、私もこの本は読みました。)
  • 10年前くらいに大学で学んだ量子力学の知識
    (だいぶ忘れてますが、量子力学が人間の直感で理解できない学問であることは知ってます
  • 英語への耐性
    (今回、実際にHHLアルゴリズム解いてみた!みたいな論文を運良く見つけられたので、その論文を読んで自分の力で最後まで解くことが出来ました。)

以上かなと思います。要するに、ほぼ何も知らない状態です、、、

なので、今回の記事では詳細な説明を隅々までするのではなく、とりあえず誰でも量子コンピュータを使ってHHLアルゴリズムで線型方程式が解けるところまで持っていくことをゴールとしたいと思います。とりあえずコードを実行して結果を得るところまでやりたい方は、以下目次の1から3まで読んでいただければと思います。もう少し詳しく知りたい方は、4まで読んでいただけると幸いです。

目次

  1. 今回解く問題の確認
  2. IBM Quantumの新規アカウントを作成
  3. 本サイトのコードをコピー&ペーストして計算実行
  4. おまけ:コードの中身を理解するためには

1.今回解く問題の確認

今回は、

Xベクトル部分を求めることを量子コンピュータで実施してみたいと思います。

問題設定として、以下の様にA, X, bを設定します。

2.IBM Quantumの新規アカウントを作成

IBM Quantumとは?

量子コンピュータを使ってみようと思った場合、既に我々は世の中に普及しているサービスを利用することで量子コンピュータを操作することが出来ます。しかも無料で。(もちろん、無料枠を超えて利用すると有料になります。)量子コンピュータを利用できるサービスは、実機を持ってくるわけにもいかないのですべてクラウドサービスとなります。なので結果としてはAWSAzureGCPがメインになりますが、IBMもクラウドサービスを提供しており、それが「IBM Quantum」です。今回は、こちらのIBM Quantumを利用していきたいと思います。(IBM Quantumの詳細はIBMHPをご覧ください)

IBM Quantum新規アカウントの作成方法

IBMのアカウント作成方法は非常に簡単です。

まず、こちらのIBM Quantumのサイトに入ってください。

そうすると、以下のような画面が表示されるかと思います。

画面左側側の「Create an IBMid」をクリックしてください。

GoogleGithubのアカウントをお持ちの場合は、そのアカウントでIBM Quantumのアカウントと連携が可能です。GoogleGithubアカウントを利用する場合は、「Create an IBMid」の上にあるGoogleGithubのアイコンをクリックしてアカウント連携を進めてください。

するとこのような画面が出てきますので、アカウントをお持ちでない場合の「IBMidの作成」をクリックしてアカウント登録をお願いします。

アカウント登録が完了すると、以下のようなIBM Quantum Platformのトップ画面が見れるかと思います。

では、次の節で中身を簡単に見ていきましょう!!

試しに中身を見てみる

まず、大まかなIBM Quantumの構成を紹介します。

画面右上のハンコ注射みたいなアイコンをクリックすると、以下のようなメニューが出てきます。

IBM Quantumは、このメニューに書いてある「Platform」、「Documentation」、「Learning」に機能を分けることが出来ます。

では、Platformから見ていきましょう。Platformは、ログインした直後の画面です。今度は、画面左上のバーコードみたいなアイコンをクリックすると、以下のようなメニューが出てきます。

Platformは、さらに3つの機能に分けることが出来ます。「Dashboard」は、ログインした時の画面で、自分が投げた計算の状態であったり、自身のAPIキーであったり、その他IBMからのニュースなどを確認できるダッシュボードになっています。

続いて、Compute resourcesをクリックしてみると以下のような画面になります。

こちらの画面では、利用する量子コンピュータのマシン情報や量子コンピュータだけでなくシミュレータの情報がまとまっています。無料枠で利用できるものは限られていますが、実際に計算を実行する際はこのマシンの中から何かしら一つを指定して計算(Job)を投げることになります。

最後に、先ほどのメニューからJobsをクリックしてみてください。まだ計算を1回も実施してない方は何も表示されないかと思いますが、計算を実行していくと以下の画面の様にJobの状態を確認することが出来ます。

画面上に、無料で利用できる枠も表示されていますのでまずはこの枠の範囲内で色々な計算を実行出来たら良いのかなと思います。「Dashboard」の説明については以上です。

続いて、「Documentation」についてですが、こちらは様々な資料が見れるところになりますのでこの場での説明は割愛します。

最後に、「Learning」を見ていきたいと思います。まず、画面右上のハンコ注射アイコンから「Learning」を選択し、続いて画面左上のバーコードアイコンをクリックすることで「Learning」のメニューが出てきます。

Homeについては、学習のホーム画面ということで色々なリソースへ移動できる画面になっています。Catalogについては、学習のコースやチュートリアルを検索できる画面になっています。(こちらは実際に初心者の方がまず触ってみる際に良いコンテンツになっているかと思いますのでご興味ありましたら是非ご覧ください)

そして、3つ目のComposerへ移動してみてください。こちらは、GUI上で量子回路を組むことが出来る画面になっています。自身で画面上に量子回路を組むことで、画面右側にその組んだ回路のPythonコードが自動的に表示されます。(以下は、量子テレポーテーションを学習した際に試しに作ってみた量子回路の参考画面)

こんな感じで見れますので、Pythonやったことなくても量子回路をまず組んでみたいという方がいらっしゃいましたらおすすめの機能になっています。

Learningの最後は、「Lab」となります。実際にクリックしてみてください。すると、画面左にはエクスプローラのようなものがあり、右側には以下のような画面が出ているかと思います。

Labでは、クラウド上で自分が作ったスクリプトを実行して量子コンピュータを動かすことが出来るようになっています。今回のHHLアルゴリズムで線型方程式を解いていくのも、こちらのLab上にスクリプトを準備して計算を実行していくことになります。では、次の章では実際に計算の準備に入っていきましょう。

2.自身のIBM Quantum Labへ本サイトのコードをコピー 

IBM Quantumで新規ファイルを作成 

では、今皆さんはIBM Quantum > Learning > Lab の画面が開かれているかと思いますので、先ほど紹介した画面のNotebookという項目の中に3つアイコンがあると思いますが一番右側のPython3 (ipykernel)をクリックしてください。

クリックすると、数秒後に以下の画面が出てくると思います。

これでコピー&ペーストの準備は完了です。では、次の節から私のコードを1ブロックずつコピーして一部編集も加えて実行していきます。

作成したファイルに本サイトのコードをコピー&ペースト

では、ここからはブロックごとにコピー&ペーストして実行していきます。

まず初めに、以下のコードをコピー&ペーストしてください。


# Imports for Qiskit
from qiskit import QuantumCircuit, execute, Aer, IBMQ
from qiskit.compiler import transpile, assemble
from qiskit.tools.jupyter import *
from qiskit.visualization import *
from qiskit import *
from qiskit import IBMQ, Aer
from qiskit import QuantumCircuit, ClassicalRegister, QuantumRegister
from qiskit.visualization import plot_histogram
import numpy as np

# Various imports 
import numpy as np

from copy import deepcopy
from matplotlib import pyplot as plt

IBMQ.save_account('<ここに自身のAPI Tokenを入れる>')
provider = IBMQ.load_account()
IBMQ.get_provider(hub='ibm-q', group='open', project = 'main')

上記のスクリプトの中で、<ここに自身のAPI Tokenを入れる>というところに自身のAPI Tokenを貼り付けて、スクリプトを実行してください。API Tokenは、IBM Quantumにログインした際の最初の画面の右上にございます。スクリプトの実行については、スクリプトの上側にある再生ボタンのアイコンをクリックするか、Shiftキーを押しながらEnterでも実行できます。

続いて、以下のコードをコピー&ペーストしてそのまま実行してください。

# Create the various registers needed
clock = QuantumRegister(2, name='clock')
input = QuantumRegister(1, name='b')
ancilla = QuantumRegister(1, name='ancilla')
measurement = ClassicalRegister(2, name='c')

# Create an empty circuit with the specified registers
circuit = QuantumCircuit(ancilla, clock, input, measurement)

circuit.barrier()
circuit.draw(output='mpl')

結果として以下のような図が出てきます。

続いて、以下のコードをコピー&ペーストしてそのまま実行してください。

def qft_dagger(circ, q, n):      
    circ.h(clock[1]);
    for j in reversed(range(n)):
      for k in reversed(range(j+1,n)):
        circ.cp(-np.pi/float(2**(k-j)), q[k], q[j]);
    circ.h(clock[0]);
    circ.swap(clock[0], clock[1]);

def qft(circ, q, n):
    circ.swap(clock[0], clock[1]);
    circ.h(clock[0]);
    for j in reversed(range(n)):
      for k in reversed(range(j+1,n)):
        circ.cp(np.pi/float(2**(k-j)), q[k], q[j]);
    circ.h(clock[1]);
def qpe(circ, clock, target):
    circuit.barrier()

    # e^{i*A*t}
    circuit.cu(np.pi/2, -np.pi/2, np.pi/2, 3*np.pi/4, clock[0], input, label='U');
    
    # e^{i*A*t*2}
    circuit.cu(np.pi, np.pi, 0, 0, clock[1], input, label='U2');

    circuit.barrier();
    
    # Perform an inverse QFT on the register holding the eigenvalues
    qft_dagger(circuit, clock, 2)
    
def inv_qpe(circ, clock, target):
    
    # Perform a QFT on the register holding the eigenvalues
    qft(circuit, clock, 2)

    circuit.barrier()

    # e^{i*A*t*2}
    circuit.cu(np.pi, np.pi, 0, 0, clock[1], input, label='U2');

    #circuit.barrier();

    # e^{i*A*t}
    circuit.cu(np.pi/2, np.pi/2, -np.pi/2, -3*np.pi/4, clock[0], input, label='U');

    circuit.barrier()
 
def hhl(circ, ancilla, clock, input, measurement):
    
    qpe(circ, clock, input)

    circuit.barrier()
    
    # This section is to test and implement C = 1
    circuit.cry(np.pi, clock[0], ancilla)
    circuit.cry(np.pi/3, clock[1], ancilla)

    circuit.barrier()
    
    circuit.measure(ancilla, measurement[0])
    circuit.barrier()
    inv_qpe(circ, clock, input)

続いて、以下のコードをコピー&ペーストしてそのまま実行してください。

# State preparation. 
intial_state = [0,1]
circuit.initialize(intial_state, 3)

circuit.barrier()

# Perform a Hadamard Transform
circuit.h(clock)

hhl(circuit, ancilla, clock, input, measurement)

# Perform a Hadamard Transform
circuit.h(clock)

circuit.barrier()


circuit.measure(input, measurement[1])

続いて、以下のコードをコピー&ペーストしてそのまま実行してください。

circuit.draw('mpl',scale=1)

結果として、以下の量子回路が出てきます。(作成した量子回路を描画している)

あと少しです!続いて、以下のコードをコピー&ペーストしてそのまま実行してください。ここで、まずシミュレーターを使った計算を実行しています。

# Execute the circuit using the simulator
#simulator = provider.get_backend('qasm_simulator')
simulator = provider.get_backend('simulator_mps')
job = execute(circuit, backend=simulator, shots=8192)

#Get the result of the execution
result = job.result()

# Get the counts, the frequency of each answer
counts = result.get_counts(circuit)

# Display the results
plot_histogram(counts)

すると、以下のような結果が出てきます。この結果のうち、01と11の結果に着目してください。求めたかった値の二乗と大体同じくらいの比になっていることがご確認いただけますでしょうか?(理想解は01:11=1:9に対して、結果は522:4600⇒1:8.8なので理想解に近い値となっていることがわかる。)

最後に、いよいよ実際の量子コンピュータで計算してみます。

以下のコードをコピー&ペーストしてそのまま実行してください。

provider.backends()

# Choose the backend on which to run the circuit
backend = provider.get_backend('ibm_osaka')

from qiskit.tools.monitor import job_monitor

# Execute the job
job_exp = execute(circuit, backend=backend, shots=8192)

# Monitor the job to know where we are in the queue
job_monitor(job_exp, interval = 2)

すると、今回は量子コンピュータに計算(Job)を投げているので、結果が返ってくるのに少し時間がかかります。しばらくすると、結果として

Job Status: job has successfully run

との出力が出てくるかと思います。

では最後に、以下のコードをコピー&ペーストしてそのまま実行してください。

# Get the results from the computation
results = job_exp.result()

# Get the statistics
answer = results.get_counts(circuit)

# Plot the results
plot_histogram(answer)

結果として、シミュレーターと同じようなヒストグラムが出てきました。しかし、実際の量子コンピュータはまだまだ誤差が大きいせいで、理想的な結果とはかけ離れた値が出力されていることも確認できるかと思います。

量子コンピュータが実用化できるレベルになるにはGoogleの発表によると2029年が目標ということなのであと約5年になります。本当に5年後に量子コンピュータが使えるようになるのかは誰にもわかりませんが、少なくとも量子コンピュータの勉強を進めておいて損はないのかなと思いました!(私はあと50年は生きる予定なので、さすがに50年後には使えるようになってると思います。)

以上