量子コンピュータの応用 - VQEの基本とこれを応用した組合せ最適化問題のプログラミング

量子コンピュータを応用したアルゴリズムの一つであるVQE(Variational Quantum Eigensolver)について、その基本からblueqatを使用したプログラミングまでの解説を行います。


VQEとは


VQEは Variational Quantum Eigensolver の略で、変分量子固有値ソルバーと言います。ソルバーとは、数式や数値解析を用いて最適な解を求める機能またはソフトウェアのことです。

よってVQEとは、変分原理を用いて量子力学の理論から最適な固有値を求めるもの、という意味になります。

変分原理と量子力学および固有値の関係について順を追って解説します。

VQEアルゴリズムを使用することで、化学分野では材料科学計算や量子化学計算へ、また、組合せ最適化問題の分野への応用などに期待されています。

それでは、ここから変分量子固有値ソルバーという言葉を手綱に、基本からしっかり解説していきましょう。


固有値とは何か?


まず、固有値について解説します。固有値は数学の分野の話になります。

固有値の数学的な定義は以下です(難しく感じる場合はこの定義部分を一旦飛ばして読んで下さい)。

◆固有値・固有ベクトル・固有空間

Fn次元ベクトル空間Vの線形変換、An次正方行列とするとき、

1. F(x)=λx かつ x0 なる xV が存在するとき、1. F(x)=λx かつ x≠0 なる x∈V が存在するとき、

λFの固有値、xを固有値λに対するFの固有ベクトル、λをFの固有値、xを固有値λに対するFの固有ベクトル、

W(λ)=xF(x)=λxW(λ)={x|F(x)=λx}

を固有値λに対するFの固有空間という。を固有値λに対するFの固有空間という。

2. Ax=λx かつ x0 なる xCn が存在するとき、2. Ax=λx かつ x≠0 なる x∈C^n が存在するとき、

λAの固有値、xを固有値λに対するAの固有ベクトル、λをAの固有値、xを固有値λに対するAの固有ベクトル、

W(λ)=xAx=λxW(λ)={x|Ax=λx}

を固有値λに対するAの固有空間という。を固有値λに対するAの固有空間という。


平たく説明すると、A=(ajk)をn x n行列、xx を列ベクトルとして、以下の一次方程式

Ax=λxAx=λx

について、ベクトルxx に行列Aを作用させても、もとのベクトルxxλλ 倍になる特別なベクトルxx を固有ベクトル、その根となるλλ を固有値と言います。

固有値は一般に複素数です。また、固有値は重根となる場合があり、その場合を縮退していると言います。

量子力学について次の節で詳しく説明しますが、量子力学においてはしばしば、固有ベクトルを固有関数、または固有状態と言います。

量子力学では、固有値と固有関数を以下のように記述します。

Aを物理量とするとき、Aを物理量とするとき、

Aφ(x)=αφ(x)  ...(1)Aφ(x)=αφ(x)  ...式(1)

ならば、αAの固有値、φ(x)を固有関数(固有状態)という。ならば、αをAの固有値、φ(x)を固有関数(固有状態)という。



量子力学とVQE


量子コンピュータは原理的に量子力学の理論を応用したコンピュータと言えます。

量子力学とは、物質を構成する原子や電子、イオン化されたこれらの粒子、さらにこれらを構成するさらに小さい各種の素粒子(まとめて量子と言います)の挙動(振る舞い)を説明する理論です。

これらの量子の振る舞いを記述する最も基本的な式がシュレーディンガーの方程式です。

以下に時間依存しないシュレーディンガー方程式を紹介します。

22md2φ(x)dx2+V(x)φ(x)=Eφ(x)-\frac{\hbar^2}{2m}\frac{d^2\varphi(x)}{dx^2}+V(x)\varphi(x) = E\varphi(x)

:プランク定数\hbar:プランク定数

V:ポテンシャルエネルギーV:ポテンシャルエネルギー

E:粒子のエネルギーE:粒子のエネルギー

上記の式について、22md2dx2+V(x)-\frac{\hbar^2}{2m}\frac{d^2}{dx^2}+V(x) をハミルトニアンHH とおくと以下のように簡潔に記述することができます。

Hφ(x)=Eφ(x)  ...(2)H\varphi(x)=E\varphi(x)  ...式(2)

上式を見ると、前節で説明した固有値・固有ベクトルの関係である式(1)と同じ形をしていることがわかります。シュレーディンガー方程式とはハミルトニアンHH の固有関数φ(x)\varphi(x) 、固有値EE を求める式であると言えます。

シュレーディンガーの方程式【式(2)】を応用して粒子の状態である固有値(エネルギー状態)を求めることがVQEの理論の裏付けとなります。


期待値の公式


変分原理についての説明に入る前に、期待値の求め方について説明します。期待値はVQEでは以下のように書かれますが、なぜそのように書くのかを確認します(この部分の説明は飛ばしても構いません)。

期待値=ψAψ〉  ...(3)期待値 = 〈\psi|A\psi〉  ...式(3)

Aは物理量Aは物理量

一般に物理量とは、位置や運動エネルギーなどの測定可能な量のことです。物理量の測定値は実数となります。

ここで以下の量子力学の公理を確認します。

量子系ではHilbert空間(ヒルベルト空間)上で考察を行います。

物理量はHilbert空間上のエルミート演算子で表現され、測定値はその固有値とします。量子状態ψ\psi の下で物理量Aを測定するとき、測定値が固有値aa となる確率PrPr は以下の式で与えられます。

Pr(A=aψ)=ψPaψ〉  ...(4)Pr(A=a | |\psi〉)=〈\psi|P_a\psi〉  ...式(4)

PaAの固有値aに対する固有射影演算子P_aはAの固有値aに対する固有射影演算子

式(4)は確率を求める式です。この式をBornの確率規則と言います。

上式は、エルミート演算子Aに対応する物理量の測定を行った時のAの固有値のうち一つが得られる確率を意味しています。

ここでは有限次元を考えますので、Hilbert空間は一般のベクトル空間と考えて良いです。

固有値aa は実数値です。また、PaP_a は固有値aa に対する固有射影演算子です。

固有射影演算子とは、固有値を与える射影演算子で、さらに射影演算子とは以下を満たす演算子です。

Pa=Pa=Pa2  ...(5)P_a=P_a†=P_a^2  ...式(5)

aPa=I  ...(6)\sum_{a}^{}P_a=\mathbb{I}  ...式(6)

式(5)を用いて、式(4)について以下が言えます。

これは、PrPr が0以上を取ることを意味します。

ψPaψ=ψPa2ψ=PaψPaψ0  ...(7)〈\psi|P_a\psi〉=〈\psi|P_a^2\psi〉=〈P_a\psi|P_a\psi〉≧0  ...式(7)

同様に、式(6)を用いて式(4)について以下が言えます。

これは、各確率の和を取ると1になることを意味します。

aψPaψ=\sum_{a}^{}〈\psi|P_a\psi〉=

ψ(aPa)ψ=ψIψ=ψ2=1  ...(8)〈\psi|(\sum_{a}^{}P_a)\psi〉=〈\psi|\mathbb{I}\psi〉=||\psi||^2 = 1  ...式(8)

式(7)式(8)より、式(4)は確率分布となることが分かります。

上記の公理を用いて期待値を算出することができます。

状態ψ|\psi〉 のもとで物理量Aを測定するときの測定値の期待値EE は固有値aa と確率PrPr の積aPaaP_a の和を取り、

Eψ[A]=aaPr(A=aψ)E_\psi[A] = \sum_{a}^{}aPr(A=a||\psi〉)

=aaψPaψ=\sum_{a}^{}a〈\psi|P_a\psi〉

=ψ(aaPa)ψ=〈\psi|(\sum_{a}^{}aP_a)\psi〉

=ψAψ〉  ...(9)=〈\psi|A\psi〉  ...式(9)

上式において、以下のスペクトル分解定理を使用しました。

任意の正規演算子は以下のように表すことができる。任意の正規演算子は以下のように表すことができる。

A=aσ(A)aPaA=\sum_{a∊\sigma(A)}^{}aP_a

よって、式(3)の期待値を求めることができました。

これを期待値の公式と呼びます。


変分原理を理解する


変分原理と言っても古典力学や電磁気学、量子力学などいくつかの分野で定義された異なる理論がありますが、ここでは量子力学での変分原理であるリッツの変分原理を取り上げます。

期待値の公式【式(3)】より、期待値を状態ψ|\psi〉 と物理量Aで求めことができることがわかりました。

すでに説明した通り、物理量Aは位置やエネルギーなどさまざまな量を表す演算子です。ここで、その物理量の一つでありエネルギー量を表す演算子のHH (ハミルトニアン)をAの代わりに用います。

ハミルトニアンHH と量子状態ψ|\psi〉 にあるときの期待値は、ある固有値(固有エネルギー)E0E_0 を用いて、以下のように書き表すことができます。

ψHψE0  ...(10)〈\psi|H|\psi〉≧E_0  ...式(10)

これは、量子力学の観点から以下のように説明できます。

すなわち、さまざまな状態において、電子などの粒子は量子力学の法則に従いいくつかのエネルギー状態を取ります。

最もエネルギーの低い状態を基底状態と呼びます。また、エネルギーが高くなるにつれて、順に第一励起状態、第二励起状態、さらに上位に状態が変化します(下図)。

よって、以下が成り立ちます。

適当な境界条件を持つ任意の状態ψ\psi に対するハミルトニアンHの期待値は、基底状態のエネルギーE0E_0 よりも常に大きいか等しい。

これは式(10)が成立することを意味します。これが変分原理です。


ψ\psi へのパラメータθ\theta の導入


上記の変分原理を応用し、エネルギーが最小となる固有値E0E_0 を求めるためにはどうすればよいでしょうか。この答えは単純で、様々な異なる状態のを与えたときに、期待値が最も小さくなるときのψ\psi を探し出せば良いわけです。そのときの期待値が固有値E0E_0 となります。

様々な異なるψ\psi を与えるために、新たに角度のパラメータθ\theta を用意します。そして、ψ\psi がパラメータθ\theta に依存する一種の関数とみなします。

期待値は以下のように書き直すことができます。

ψ(θ)Hψ(θ)〈\psi(\theta)|H|\psi(\theta)〉


VQEアルゴリズムを用いたプログラミングの実践


変分原理を理解した上で、VQEを用いたプログラミングに入りたいと思います。

プログラミングには、量子コンピューターシミュレータとして評価の高い、blueqat株式会社(東京・本郷)が提供するblueqat(ブルーキャット)を使用します。

また、解く問題としては組合せ最適化問題の一つである頂点被覆問題を取り扱います。


(1)頂点被覆問題について


頂点被覆問題(Vertex Cover)とは、グラフ・ネットワーク問題の一つで、複数の頂点とそれを結ぶ辺があり、どの辺も少なくともその始点または終点のどちらかが、二つにクラス分けした頂点群のうちの一方のクラスに属する頂点(のいずれか)と接するもので、そのクラスを最小にする問題です。

数学的に説明すると以下となります。

頂点の部分集合XVに対して、どの枝も少なくとも一方の端点がXに含まれるとき、すなわち、任意の(i,j)Eに対してXi,j0が成り立つとき、Xを頂点被覆という。頂点の部分集合X⊆Vに対して、どの枝も少なくとも一方の端点がXに含まれるとき、すなわち、任意の(i,j)∈Eに対してX∩i,j≠0が成り立つとき、Xを頂点被覆という。

以下に図を示します。

ここで、v1,v2,v3,v4,v5は頂点を、e1,e2,e3,e4,e5は辺を表します。図1の通り、どの辺についてもその端点の一方が赤丸の頂点に属しています。赤丸の数が最小になるケースを求める問題が頂点被覆問題です。

例えば、下の図2のような頂点を選択した場合、頂点の数は3個なので最適解にはなりません。

この例題を解く際には、イジング模型に基づくハミルトニアンの定式化を行います。今回の例題では、ハミルトニアンHH を以下のように定義します。

H=532x132x232x352x412x5+x1x2+x1x4+x2x3+x3x4+x4x5H=5-\frac{3}{2}x_1-\frac{3}{2}x_2-\frac{3}{2}x_3-\frac{5}{2}x_4-\frac{1}{2}x_5+x_1x_2+x_1x_4+x_2x_3+x_3x_4+x_4x_5

上記の定式化についての詳細は、以下の技術ブログをご参照下さい。

■Blueqatで簡単な頂点被覆問題を解く

https://qiita.com/ttabata/items/b4ab222ef9c5ee99233c


(2)変分原理によるVQEアルゴリズムの計算手順


変分原理を用いたVQEアルゴリズムは以下の手順で計算を行います。

①パラメータθ\theta を用いた量子回路を準備(ansatz:アンザッツ)

②量子ビットを測定し、期待値を算出する

③古典最適化計算により、期待値を最小にするパラメータθ\theta を探索する

VQEでは、量子コンピュータによる期待値の測定と、古典最適化計算(古典コンピュータの計算)によるパラメータθ\theta の探索処理の両方を用いて解を算出します。

処理①~③を模式図で示すと以下となります。こちらはVQEの論文からの引用です。

引用:(https://arxiv.org/pdf/1304.3061.pdf)

処理①に記載のansatzとは、量子計算に必要な量子状態ψ|\psi〉 を生成するための量子回路のことで、任意回転ゲートと入力パラメータθ\theta を用いて量子状態を生成します。この量子状態を用いて期待値を計算します。


(3)blueqatによるVQEアルゴリズム計算のプログラミング


blueqatには上記の計算手順①~③を内部関数で実行してくれる機能が備わっており、プログラマーはパラメータや期待値計算、最小値の探索に手をかけなくても良くなっています。これはblueqatが優れている点の一つと言えます。

プログラマーがすべきことは、パラメータとハミルトニアンHH を定義することのみです。

では、さっそくパラメータ及びハミルトニアンHH の定義を行いましょう。前掲の例題である頂点被覆問題から、ハミルトニアンを以下とします。

H=532q032q132q252q312q4+q0q1+q0q3+q1q2+q2q3+q3q4H=5-\frac{3}{2}q_0-\frac{3}{2}q_1-\frac{3}{2}q_2-\frac{5}{2}q_3-\frac{1}{2}q_4+q_0q_1+q_0q_3+q_1q_2+q_2q_3+q_3q_4

例題では変数をxx としていますが、ここでは量子ビット(quantum bit)の意味で、アルファベットqq を使用しました。また、インデックス番号を0始まりとしています。

基本的にはこのハミルトン式を組み込めば良いです。プログラムのソースコードは以下となります。

from blueqat import Circuit
from blueqat.pauli import X, Y, Z, I
from blueqat.pauli import qubo_bit as q
from blueqat.vqe import AnsatzBase, Vqe


class QubitAnsatz(AnsatzBase):
    def __init__(self, hamiltonian):
        super().__init__(hamiltonian, 10)
        self.step = 1
    
    def get_circuit(self, params):
        a, b, c, d, e, f, g, h, i, j = params   # パラメータの定義1
        return Circuit().ry(a)[0].rz(b)[0].ry(c)[1].rz(d)[1].ry(e)[2].rz(f)[2].ry(g)[3].rz(h)[3].ry(i)[4].rz(j)[4]   # パラメータの定義2


H = 5 - 3*q(0)/2 - 3*q(1)/2 - 3*q(2)/2 - 5*q(3)/2 - 1*q(4)/2 + q(0)*q(1) + q(0)*q(3) + q(1)*q(2) + q(2)*q(3) + q(3)*q(4)   # ハミルトニアンの定義
H = H.to_expr().simplify()


runner = Vqe(QubitAnsatz(H))
result = runner.run()


print('Result by VQE')
print(result.most_common())


ソースコード中で「パラメータの定義1」と「パラメータの定義2」の記述を行いました。

今回の例題では、ハミルトニアン内に含まれる量子ビット(q)がq0~q4まで5個あります。それぞれの量子ビットに対して回転ゲートry,rzを適用するため、回転ゲートの総数は10個となり、それぞれにパラメータを設定するため、パラメータ数も10個です。

よって、「パラメータの定義1」に10個分のパラメータa~jを宣言しています。

「パラメータの定義2」では上記に定義した10個のパラメータを用いて、回転ゲートを宣言します。量子ビット(q)がq0~q4まで5個分あるので、回転ゲートry,rzもそれぞれ、ry(0)~ry(4)の5個分と、rz(0)~rz(4)の5個分を宣言します。

「ハミルトニアンの定義」にハミルトニアンを定義します。

プログラミングは以上です。


(4)実行結果の確認


実行結果は以下となりました。結果の見方は例題と同じです。

q0~q4は以下の値を取ります。

ここでは、0,1はクラス分けを意味しています。0を取る頂点のクラスと、1を取る頂点のクラスに分かれます。

このクラス分けが今回の問題の解答となります。


つまり、q1とq3が同じクラス、q0,q2,q4が別の同じクラスに属することがわかります。

q1,q3は図1の赤丸の頂点、q0,q2,q4はその他の丸の頂点に分けられます。

この分け方が問題の解となります。


以上でプログラミングの説明を終わります。

Tetsuro Tabata
Comments
Tetsuro Tabata
Related posts

blueqat Inc.

Shibuya Scramble Square 39F 2-24-12, Shibuya, Shibuya-ku, Tokyo
Contact: info@blueqat.com