補遺: 適応的量子モンテカルロ法

量子アニーリングを古典計算機で調べるとき、もっとも基本的な道具の1つが量子モンテカルロ法です。横磁場を持つ量子Ising模型を、鈴木--Trotter分解によって「虚時間方向に並んだ古典スピンの模型」へ写します。2日目の講義で出てきたSQA、Simulated Quantum Annealingは、この見方に基づいています。

標準的な横磁場なら、この変換はきれいに働きます。しかし、量子アニーリングの性能を高めようとして、より凝った量子揺らぎを入れると、すぐに符号問題が顔を出します。非対角成分から負の確率重みが生じ、モンテカルロ法で「確率としてサンプルする」という足場が崩れるからです。

適応的量子モンテカルロ法は、この困難をすべてのnon-stoquasticハミルトニアンに対して解決する万能法ではありません。狙いはもっと具体的です。横方向の量子揺らぎが、全スピンの横磁化

$${\hat{m}}_x=\frac{1}{N}\sum _{i=1}^N{\hat{\sigma }}_i^x$$

の関数として書ける場合、その非線形な量子揺らぎを、横磁化に応じて値が変わる有効横磁場へ置き換えます。すると、各ステップで解くべき問題は通常の横磁場Ising模型になります。横磁場そのものは、シミュレーションの途中で測定された横磁化に合わせて更新されます。この「有効横磁場を自分で調整しながら量子モンテカルロを走らせる」というのが、適応的という名前の中身です。

1. 標準的な横磁場QMCを思い出す

出発点として、標準的な横磁場付きIsing ハミルトニアンを置きます。

$$\hat{H}=H_0({\hat{\mathit{\sigma }}}^z)-\mathrm{\Gamma }\sum _{i=1}^N{\hat{\sigma }}_i^x$$

$H_0$ は最小化したい古典Ising ハミルトニアンです。${\hat{\sigma }}_i^z$ の固有値を ${\sigma }_i=\pm 1$ と見れば、$H_0$ は古典スピン配置のエネルギーを返します。一方、${\hat{\sigma }}_i^x$ は ${\sigma }_i^z$ の基底では非対角です。古典スピンを反転させる方向の量子揺らぎだと思えばよいです。

分配関数は

$$Z=\mathrm{Tr}\exp (-\beta \hat{H})$$

です。$H_0$ と横磁場項は一般に可換ではないので、そのまま古典確率分布には見えません。そこで逆温度方向を $M$ 個に分け、鈴木--Trotter分解を使います。

$$Z\simeq \sum _{\{{\sigma }_{i,t}\}}\exp [-\frac{\beta }{M}\sum _{t=1}^M{H}_0({\mathit{\sigma }}_t)+B(\mathrm{\Gamma })\sum _{i=1}^N\sum _{t=1}^M{\sigma }_{i,t}{\sigma }_{i,t+1}]$$

ここで $t$ はTrotter方向の添字で、周期境界条件 ${\sigma }_{i,M+1}={\sigma }_{i,1}$ を入れます。横磁場は、Trotter方向に隣り合うコピーどうしの結合

$$B(\mathrm{\Gamma })=\frac{1}{2}\log \coth \left(\frac{\beta \mathrm{\Gamma }}{M}\right)$$

へ変わります。$\mathrm{\Gamma }$ が小さいと $B(\mathrm{\Gamma })$ は大きくなり、Trotter方向のコピーはそろいやすくなります。$\mathrm{\Gamma }$ が大きいとコピーは揺らぎやすくなります。量子揺らぎを、もう1つの古典方向の相互作用として読むわけです。

この形なら重みは正です。したがって、メトロポリス法、熱浴法、交換モンテカルロ法など、普通のマルコフ連鎖モンテカルロでサンプルできます。

2. non-stoquasticな量子揺らぎ

量子アニーリングでは、横磁場だけでなく、横方向の相互作用を加えたくなることがあります。典型例は、反強磁性的なXX相互作用です。横方向にもスピンどうしを相互作用させることで、一次相転移を避け、エネルギーギャップの閉じ方を緩やかにできる場合があるからです。

ここで扱う量子揺らぎは、次の形にまとめられます。

$$\hat{H}=H_0({\hat{\mathit{\sigma }}}^z)-Nf({\hat{m}}_x)$$

$f$ は1変数関数です。標準的な横磁場は

$$f(m_x)=\mathrm{\Gamma }{m}_x$$

に対応します。なぜなら $-Nf({\hat{m}}_x)=-\mathrm{\Gamma }\sum _i{\hat{\sigma }}_i^x$ だからです。

反強磁性的なXX相互作用を含む場合は、たとえば

$$f(m_x)=\mathrm{\Gamma }{m}_x-\frac{\gamma }{2}{m}_x^2$$

と書けます。このとき量子揺らぎの項は

$$-Nf({\hat{m}}_x)=-\mathrm{\Gamma }\sum _i{\hat{\sigma }}_i^x+\frac{\gamma }{2N}{\left(\sum _i{\hat{\sigma }}_i^x\right)}^2$$

です。$\gamma >0$ なら、横方向には反強磁性的な相互作用が入ります。この項があると、${\hat{\sigma }}^z$ 基底で見た非対角行列要素の符号が標準横磁場の場合のようにそろわなくなり、素朴な鈴木--Trotter分解では符号問題が出ます。

ここで大事なのは、難しさの源が「横方向にある非線形関数 $f({\hat{m}}_x)$」だという点です。$H_0$ の形にはかなり自由度があります。全結合スピングラスでも、有限次元格子のスピングラスでも、制約を含む最適化ハミルトニアンでもかまいません。一方で、量子揺らぎの側は ${\hat{m}}_x$ の関数として書ける、という構造を使います。

3. 補助変数で横磁化を取り出す

$f({\hat{m}}_x)$ が非線形だから困るなら、${\hat{m}}_x$ をいったん独立した変数として取り出せばよいです。これは前の補遺で見たハバード・ストラトノビッチ変換と同じ匂いのする操作です。多体相互作用を、補助変数との結合へほどく。

Trotter分解した各スライス $t$ で、横磁化に対応する補助変数 $m_{x,t}$ を入れます。形式的には、恒等式

$$1=\int d{m}_{x,t}\delta (m_{x,t}-\frac{1}{N}\sum _{i=1}^N{\hat{\sigma }}_{i,t}^x)$$

を挿入します。$\delta$ 関数は、「$m_{x,t}$ は、そのスライスの横磁化そのものである」という制約を課しています。

さらに $\delta$ 関数をFourier積分で表します。規格化定数を省けば、

$$\delta (y)={\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}\frac{d\tilde{m}}{2\pi }\exp \left(i\tilde{m}y\right)$$

のように書けます。このあと積分経路を取り直し、係数を吸収して、実数の共役場として $\tilde{m}$ を扱います。$\tilde{m}$ は $m$ と共役な補助変数です。前の補遺の $\nu$ が制約のずれに対応するラグランジュ乗数だったように、ここでの $\tilde{m}$ は横磁化を固定するための場です。

この操作を使うと、指数の中に現れる $f(m_{x,t})$ と ${\tilde{m}}_{x,t}$ の関係は、概念的には

$$\frac{\beta N}{M}[f(m_{x,t})-{\tilde{m}}_{x,t}{m}_{x,t}]+\frac{\beta {\tilde{m}}_{x,t}}{M}\sum _i{\hat{\sigma }}_{i,t}^x$$

となります。最後の項に注目してください。これは、強さ ${\tilde{m}}_{x,t}$ の横磁場と同じ形です。つまり、もともとの非線形な量子揺らぎ $Nf({\hat{m}}_x)$ が、補助変数を固定したもとでは「普通の横磁場」に見えるようになります。

熱力学極限では、$m_{x,t}$ と ${\tilde{m}}_{x,t}$ の積分は鞍点で評価できます。$m_{x,t}$ で微分すると、

$${\tilde{m}}_{x,t}=f^{\prime }(m_{x,t})$$

が得られます。$f^{\prime }(m)$ は $f$ の微分です。この式が適応的量子モンテカルロ法の心臓部です。横磁化 $m_x$ が決まると、有効横磁場 ${\tilde{m}}_x$ が $f^{\prime }(m_x)$ として決まります。

静的近似では、Trotterスライスに依存しない

$$m_{x,t}=m_x,{\tilde{m}}_{x,t}={\tilde{m}}_x$$

を仮定します。この近似のもとでは、解くべき量子系は

$${\hat{H}}_{\mathrm{e}\mathrm{f}\mathrm{f}}(m_x)=H_0({\hat{\mathit{\sigma }}}^z)-f^{\prime }(m_x)\sum _{i=1}^N{\hat{\sigma }}_i^x$$

です。右辺は、標準的な横磁場付きIsing模型です。non-stoquasticな項は、有効横磁場の値をどう決めるかという外側の問題へ移りました。

4. 条件付き分布として普通のQMCを走らせる

$m_x$ を固定すると、有効横磁場は

$${\mathrm{\Gamma }}_{\mathrm{e}\mathrm{f}\mathrm{f}}=f^{\prime }(m_x)$$

です。したがってTrotter方向の結合は

$$B(m_x)=\frac{1}{2}\log \coth \left(\frac{\beta f^{\prime }(m_x)}{M}\right)$$

になります。固定した $m_x$ のもとでの条件付き分布は

$$P(\{{\mathit{\sigma }}_t\}\mid m_x)\propto \exp [-\frac{\beta }{M}\sum _{t=1}^M{H}_0({\mathit{\sigma }}_t)+B(m_x)\sum _{i,t}{\sigma }_{i,t}{\sigma }_{i,t+1}]$$

です。これは標準的なSQAと同じ形です。違うのは、横磁場の値が固定パラメータ $\mathrm{\Gamma }$ ではなく、いま仮に置いている $m_x$ から $f^{\prime }(m_x)$ として決まることだけです。

ここでは ${\mathrm{\Gamma }}_{\mathrm{e}\mathrm{f}\mathrm{f}}>0$ として書いています。通常の横磁場QMCとして扱える範囲に有効横磁場を置くことが、符号問題を避けるための実装上の前提になります。

補助変数 $m_x$ 自身の分布も書けます。規格化を除けば、

$$P(m_x)\propto Z(m_x)\exp [-\beta N\{m_x{f}^{\prime }(m_x)-f(m_x)\}]$$

です。$Z(m_x)$ は、条件付き分布 $P(\{{\mathit{\sigma }}_t\}\mid m_x)$ の規格化因子です。この式は少し抽象的に見えますが、意味は単純です。補助変数 $m_x$ のよさは、固定した $m_x$ のもとでスピン配置がどれだけ多くの重みを持つか、そして $f(m_x)$ とその接線 $m_x{f}^{\prime }(m_x)$ の差で決まります。

大きな系では、$P(m_x)$ は鞍点の近くに鋭く集中します。したがって、$m_x$ をランジュバン過程で揺らしながらサンプルすることもできますが、多くの場合は鞍点条件を直接使うほうが見通しがよいです。

5. 横磁化を測り、有効横磁場を更新する

固定した $m_x$ のもとでQMCを走らせたら、横磁化を推定します。Trotter表示では、横磁化は隣り合うTrotterスライスの相関から読めます。1つの便利な推定量は

$${\hat{m}}_x=\frac{1}{NM}\sum _{i=1}^N\sum _{t=1}^M{[\tanh \left(\frac{\beta {\mathrm{\Gamma }}_{\mathrm{e}\mathrm{f}\mathrm{f}}}{M}\right)]}^{{\sigma }_{i,t}{\sigma }_{i,t+1}}$$

です。指数に ${\sigma }_{i,t}{\sigma }_{i,t+1}$ が入っているので、一瞬ぎょっとしますが、意味は素朴です。隣同士が同じなら $\tanh (\beta {\mathrm{\Gamma }}_{\mathrm{e}\mathrm{f}\mathrm{f}}/M)$、反対ならその逆数が入ります。Trotter方向の「切れ目」の出方から、横方向にどれだけスピンが揺らいだかを見ているわけです。

鞍点条件は

$$m_x={⟨{\hat{m}}_x⟩}_{P(\mathit{\sigma }\mid m_x)}$$

です。左辺は有効横磁場を決めるために仮に置いた横磁化、右辺はその有効横磁場で実際にQMCを走らせて得られる横磁化です。この2つが一致したところが、もとの非線形な量子揺らぎを正しく反映する点になります。

したがって、適応的量子モンテカルロ法の手順は次のように書けます。

1. 横磁化の仮値 $m_x^{(s)}$ を置きます。
2. ${\mathrm{\Gamma }}_{\mathrm{e}\mathrm{f}\mathrm{f}}^{(s)}=f^{\prime }(m_x^{(s)})$ を計算します。
3. ${\mathrm{\Gamma }}_{\mathrm{e}\mathrm{f}\mathrm{f}}^{(s)}$ を持つ標準的な横磁場QMCを走らせる。
4. サンプルから ${\hat{m}}_x^{(s)}$ を見積もる。
5. $m_x$ を更新し、物理量が収束するまで繰り返します。

数値的には、更新を少しなめらかにするために

$$m_x^{(s+1)}=(1-\alpha )m_x^{(s)}+\alpha {\hat{m}}_x^{(s)},0<\alpha \le 1$$

のように混合してもよいです。$\alpha =1$ なら、推定した横磁化をそのまま次の仮値にします。$\alpha$ を小さくすると、振動しやすい場合にゆっくり進められます。

この方法の面白いところは、non-stoquasticなハミルトニアンを直接サンプルしていない点です。各ステップでサンプルしているのは、符号問題のない標準横磁場の模型です。けれども横磁場の強さを、サンプルから得た横磁化に合わせて変えることで、もとの非線形な量子揺らぎの効果を取り込んでいきます。

6. 反強磁性XX相互作用で見る

反強磁性XX相互作用を含む関数

$$f(m_x)=\mathrm{\Gamma }{m}_x-\frac{\gamma }{2}{m}_x^2$$

では、

$$f^{\prime }(m_x)=\mathrm{\Gamma }-\gamma m_x$$

です。したがって有効横磁場は

$${\mathrm{\Gamma }}_{\mathrm{e}\mathrm{f}\mathrm{f}}=\mathrm{\Gamma }-\gamma m_x$$

になります。横磁化が大きいほど、有効横磁場は弱くなります。これは反強磁性的なXX相互作用の直観と合っています。全員が同じ横方向にそろうことを、$\gamma$ の項が押し戻すからです。

検証によく使われる簡単な古典部分は、無限レンジ強磁性模型です。

$$H_0(\mathit{\sigma })=-N{\left(\frac{1}{N}\sum _{i=1}^N{\sigma }_i^z\right)}^2-h\sum _{i=1}^N{\sigma }_i^z$$

この模型は平均場解析で厳密に調べやすいので、適応的量子モンテカルロ法が正しい物理量に近づくかを確認する舞台として使いやすいです。典型的な数値実験では、$N=4,8,16,32$ の系を、$M=128$ のTrotter数、$\beta =50$ の低温で調べ、内部エネルギー、縦磁化、横磁化を比較します。有限サイズ効果は残りますが、反強磁性XX相互作用を入れた結果に近づくことが確認できます。

もちろん、ここで大切なのは「この模型だけが解けた」という話ではありません。平均場で答えが分かる模型を物差しにして、符号問題を避ける置き換えが機能するかを確かめる。いったんその道筋が見えれば、平均場だけでは扱いにくい $H_0$ にも量子モンテカルロの装置を使える可能性が出てきます。

7. データ解析としての読み方

同じ考え方には、もう1つの実装があります。適応的に横磁場を更新する代わりに、普通の横磁場QMCをたくさん走らせておき、あとから交点を読む方法です。

まず、さまざまな横磁場 ${\mathrm{\Gamma }}_{\mathrm{e}\mathrm{f}\mathrm{f}}$ について標準的なQMCを実行し、横磁化

$$m_x=m_x({\mathrm{\Gamma }}_{\mathrm{e}\mathrm{f}\mathrm{f}})$$

を測ります。これは「普通の横磁場をかけたとき、系がどれだけ横方向に磁化するか」という応答曲線です。

一方、もとの非線形な量子揺らぎでは、鞍点条件から

$${\mathrm{\Gamma }}_{\mathrm{e}\mathrm{f}\mathrm{f}}=f^{\prime }(m_x)$$

が必要です。逆関数で書けば

$$m_x=(f^{\prime }{)}^{-1}({\mathrm{\Gamma }}_{\mathrm{e}\mathrm{f}\mathrm{f}})$$

です。したがって、標準QMCから得た曲線 $m_x({\mathrm{\Gamma }}_{\mathrm{e}\mathrm{f}\mathrm{f}})$ と、関数 $f$ から決まる曲線 $(f^{\prime }{)}^{-1}({\mathrm{\Gamma }}_{\mathrm{e}\mathrm{f}\mathrm{f}})$ の交点を探せばよいです。交点では、標準QMCで得られる横磁化と、非線形な量子揺らぎが要求する横磁化が一致しています。

反強磁性XX相互作用の例では、

$$m_x=\frac{\mathrm{\Gamma }-{\mathrm{\Gamma }}_{\mathrm{e}\mathrm{f}\mathrm{f}}}{\gamma }$$

です。これは直線です。標準QMCの結果を横磁場ごとに並べ、その曲線とこの直線が交わる場所を読みます。そこで得られた横磁化や内部エネルギーを、non-stoquasticハミルトニアンの物理量として再構成します。

この方法は、横磁場を逐次更新しないので、計算の見通しがよい場合があります。すでに標準QMCのデータがあるなら、あとから別の $\mathrm{\Gamma },\gamma$ に対する結果を読み替えることもできます。一方で、交点が複数ある場合には、自由エネルギーを比べて熱力学的に安定な解を選ばなければなりません。

8. 何を可能にし、何を残すのか

適応的量子モンテカルロ法は、符号問題そのものを一般に消す方法ではありません。適用しやすいのは、量子揺らぎが少数の横方向の集団変数で表せる場合です。ここで扱ったもっとも単純な場合では、集団変数は全体の横磁化 $m_x$ でした。部分集合ごとの横磁化に分ければ、不均一なXX相互作用へ拡張する道もありますが、補助変数の数が増えれば計算量も増えます。

また、鞍点を使うため、系が十分大きいほど理屈は鋭くなります。小さい系では有限サイズ効果が残ります。横磁化を推定しながら有効横磁場を更新するので、通常の横磁場QMCより収束が遅くなることもあります。反復の初期値、平均を取る長さ、混合率 $\alpha$、交換モンテカルロを使うかどうか、といった実装上の選択も効いてきます。

それでも、この方法にははっきりした価値があります。non-stoquasticな量子揺らぎは、量子アニーリングを単なる横磁場から広げるための重要な候補です。しかし、すぐに符号問題が出るため、古典計算機で大きな系を調べることが難しいです。適応的量子モンテカルロ法は、その一部を通常の横磁場QMCの反復として調べるための足場を与えます。

横磁場だけなら古典計算機で模擬できるからつまらない、と急いで切ってしまう必要はありません。どの量子揺らぎなら古典的に追えるのか、どこから先で本当に難しくなるのでしょうか。その境目を調べること自体が、量子アニーリングの研究です。適応的量子モンテカルロ法は、その境目に立って、non-stoquasticハミルトニアンの一部を手で触れる形にする方法だと読めます。

参考文献: [1]。

参考文献

1. Masayuki Ohzeki, “Quantum Monte Carlo simulation of a particular class of non-stoquastic Hamiltonians in quantum annealing,” Scientific Reports 7, 41186 (2017). 出典