3重臨界3状態ポッツ模型の表面臨界現象の共形場理論・テンソルネットワーク解析

平衡統計力学の教科書を紐解くと、「1次元古典系では有限温度での自発的対称性の破れは存在しない」ことの証明を見つけることができる。
つまり、1次元古典系は多体問題としてほとんどつまらない。
では孤立した1次元古典系ではなく、2次元古典系の表面、すなわち2次元のバルクと接触した1次元系を考えてみることにすると、単純な1次元系とは異なる振る舞いが観測できることがある。
とりわけバルクが臨界状態である場合には、その強い相関が表面における秩序の発達を扶翼し、様々な表面臨界現象が観測されることが知られている。
その最たる例は2次元古典3重臨界イジング模型であり、バルクを3重臨界点にfine tuningして表面の温度を変化させると、1次元表面に自発磁化を生じる有限表面温度相転移が観測される[1]。
また、表面に磁場を作用させることでさらなる表面相転移が起きることも知られており、単純な古典1次元系には見られない多彩な相転移現象が展開される。

一方で、2次元臨界古典系の表面相転移現象は境界共形場理論(Boundary conformal fieled theory; BCFT)による記述が有効である。
共形不変性の制限から、可能な境界状態の完全な分類が可能となる場合もあり、さらには固定点の安定性や厳密なスケーリング次元の値を調べることができるため、表面臨界現象を詳細に理解する上で極めて重要である。

ここでは3重臨界イジング模型を対称性\(Z_2\)から\(S_3\)へ拡張した古典スピン模型である、3重臨界3状態ポッツ模型を考える。
この模型も、バルクが3重臨界点である場合に多彩な表面相転移現象を生じることがモンテカルロシミュレーションにより報告されていたが[2]、BCFTによる詳細な解析はなされていなかった。
そこで我々は、ミニマル系列境界共形場理論の\(ADE\)分類理論[3]およびテンソルネットワークくりこみ(Tensor network renormalization; TNR)[4]の手法を用いて、この模型の表面相転移を詳細に調べた[5]。

3状態ポッツ模型の3重臨界現象は中心電荷\(c=6/7\)の\(D\)型ミニマル模型に対応し、そのモジュラー不変な分配関数は\(ADE\)分類理論においてリー代数のペア \((D_4, A_6)\) で分類される。
この共形場理論に対して可能な共形不変な境界状態は12個存在し(ディンキン図形\(D_4\)と\(A_6\)の頂点で特徴づけられる)、特に\(D_4\)のtrialityにより3個の\(Z_3\)対称な境界状態と9個の\(Z_3\)対称性の破れた境界状態が存在することが文献[3]の解析から導かれる。

上記の12個の境界状態をスピン模型の描像で理解するために、TNRによる数値計算を格子上の3状態希釈ポッツ模型で行った。
筆者らが提案したTNRの拡張手法により格子模型から高精度にBCFTのスペクトルを計算することができ、これを通じてBCFTで得られた境界状態と格子模型の示す相図との対応を調べることができる。
数値計算の結果、上記の12個の境界状態のうち\(Z_3\)対称な1つを除いて相図上で対応する境界状態を特定し、その物理的な意味を明らかにした(図参照)。
残った一つの境界状態は12個のうちで最も不安定な固定点に対応し、負のボルツマン重率を用いて実現されるような境界状態に対応すると予想している。
その具体的な格子模型上での実現方法を調べることは今後の研究課題である。
(by 飯野隼平)

図:TNRで得られた3重臨界3状態ポッツ模型の相図。多彩な表面臨界現象が観測される。

References:
[1] I. Affleck, J. Phys. A 33(37), 6473 (2000).
[2] Y. Deng and H. W. J. Blöte, Phys. Rev. E 70, 035107 (2004); Phys. Rev. E 71, 026109 (2005).
[3] R. E. Behrend, P. A. Pearce, V. B. Petkova, and J.-B. Zuber, Nucl. Phys. B 579(3), 707 (2000).
[4] G. Evenbly and G. Vidal, Phys. Rev. Lett. 115, 180405 (2015); S. Iino, S. Morita, and N. Kawashima, Phys. Rev. B 101, 155418 (2020).
[5] S. Iino, arXiv:2007.03182; J. Stat. Phys. 182, 56 (2021).

任意の確率でビットを立てる高速なアルゴリズムの開発

Nビットのビット列について、各ビットを確率 \(p\) でランダムに1にするにはどうすればよいだろうか。単純にそれぞれのビットごとに1にするか0にするかを決めると、 \(N\) 回の乱数生成が必要になる。しかし、アルゴリズムを工夫することで乱数生成回数を減らすことができる。我々はこのようなランダムなビット列を生成するアルゴリズムをいくつか考案し、それを組み合わせることで劇的に計算量を減らすことに成功した[1]。我々が提案したのはBinomial-Shuffle法、Poisson-OR法、有限桁法の三種類と、そのハイブリッド法である。

図1: 三種類のランダムビット生成アルゴリズム

Binomial-Shuffle法は、予めいくつのビットが立つかを計算し、その後でビットをシャッフルするアルゴリズムである。立つビット数の計算にWalker’s Alias法を、ビットのシャッフルにFloydのサンプリング法を使うことで、全体で \(pN+1\) 回の乱数生成でランダムビット列を生成できる。Poisson-OR法は、 \(N\) ビットのうちランダムに1ビットだけビットが立っているビット列を、 \(k\) 個論理和(OR)を取るアルゴリズムである。これはTodo and Suwaのアルゴリズム[2]の特殊な例になっており、ビット列の数 \(k\) をポアソン分布に従うように生成することで、指定の確率でビットが立つようにできる。この手法は \(-N \log(1-p) + 1\) 回の乱数生成を伴う。有限桁法は、確率 \(p\) が二進数表記で有限桁で表せる場合に有効な手法である。確率 \(p\) が二進数表記で \(n\) 桁で表示できる場合、確率 \(1/2\) でランダムにビットが立っているビット列を \(n\) 個生成し、適切に論理和と論理積を取ることで、その確率でビットが1となるようにできる。したがって、確率が \(n\) 桁で表現される場合、有限桁法は \(n\) 回の乱数生成でビット列を生成できる。

Binomial-Shuffle法、Poisson-OR法ともに、確率 \(p\) の値が小さい場合に有効だが、 \(p\) が大きくなると必要となる乱数生成数が増えてしまう。そこで我々は、まず有限桁法で \(p\) に近い確率 \(p’\) を持つビット列を生成し、その差 \(p’-p\) をBinomial-Shuffle法もしくはPoisson-OR法で補正するハイブリッド法を考案した。これにより、32ビットであればたかだか7回、64ビットでもたかだか8回の乱数生成で任意の確率 \(p\) でビットが立っているビット列を生成できること示した。

図2: Directed-Percolationのマルチスピンコーディング実装。臨界点上での時間発展にかかった時間。スカラー実装(Scalar)、有限桁法+Binomial-Shuffle (BS)法による、有限桁法+Poisson-OR (PO)法による補正それぞれについて、32ビット、64ビットの場合を表示。

我々はこのランダムビット列生成アルゴリズムをDirected-Percolationのマルチスピンコーディングに応用し、スカラー実装に比べて14倍の高速化に成功した。このアルゴリズムは他のモデルへの応用も期待されている。

References

  1. H. Watanabe, S. Morita, S. Todo, and N. Kawashima, J. Phys. Soc. Jpn. 88, 024004 (2019)
  2. S. Todo and H. Suwa, J. Phys.: Conf. Ser. 473, 012013 (2013).

乱択特異値分解を用いたテンソルくりこみ群

テンソルネットワーク法は,多体問題のための強力な数値計算手法として研究開発が進められている.テンソルくりこみ群 (Tensor Renormalization Group, TRG) [1] に代表される実空間のくりこみの方法は,粗視化によりテンソルネットワークの縮約を効率良く計算する近似手法である.テンソルの次元(ボンド次元)を大きくすることで近似の精度が向上する.しかしながら,TRGでは計算コストはボンド次元の6乗に比例して増大する.そのため,計算精度を損なわず計算量を削減する手法の開発が求められている.

図1: TRGにおけるテンソルの変形.4階のテンソル生成を回避することで縮約の計算量を削減する.

本研究で我々は,乱択アルゴリズムによる特異値分解 (RSVD) [2] を用いたTRGの新しいアルゴリズムを提案し,TRGの計算コストをボンド次元の5乗に削減することに成功した [3].TRGの主要な演算は「分解」と「縮約」であるが,どちらも計算量がボンド次元の6乗に比例する.RSVDによる分解の計算量削減だけでなく,4階のテンソルを明示的に生成しないことで縮約の計算量も同時に削減する点が提案手法の肝となっている.RSVDは乱数を用いるため統計誤差が含まれるが,オーバーサンプリング変数を増やすことで誤差を減らすことができる.我々は2次元Ising模型におけるベンチマークにより,オーバーサンプリング変数がボンド次元と同じ程度であれば全特異値を求める場合と同じ精度で自由エネルギーが計算できることを示した.

図2: TRGの1ステップ当たりの計算時間.

他のテンソルネットワーク法でも主要な演算は,TRGと同じく「分解」と「縮約」である.そのため,乱択アルゴリズムを用いた手法は,数多くのテンソルネットワーク法に応用されることが期待される.

References

  1. M. Levin and C. P. Nave, Phys. Rev. Lett. 99, 120601 (2007).
  2. N. Halko, P. G. Martinsson, and J. A. Tropp, SIAM Rev. 53, 217 (2011).
  3. S. Morita, R. Igarashi, H.-H. Zhao, and N. Kawashima, Phys. Rev. E 97, 033310 (2018).

テンソルネットワーク法を用いたNa2IrO3有効模型の基底状態解析

トポロジーによって特徴づけられる量子状態の研究が盛んにおこなわれているが,キタエフモデルはそのような状態を基底状態として持ち,かつ厳密に解けるモデルである.キタエフモデルは非常に人工的なモデルであるために,現実の物質でこれに対応するものが存在するかどうかに関心が集まっている.Jackeli と Khaliullin [1] はスピン軌道相互作用が強い物質においては,キタエフモデル的相互作用が重要になることを指摘し,スピン軌道相互作用が大きなイリジウム原子を含んだNa2IrO3などの化合物が候補物質とされた.山地ら [2] は第一原理計算に基づいて,Na2IrO3 の磁性に関する有効モデルを提案した.このモデルは,ハチの巣格子上のキタエフモデルのハミルトニアンに,ハイゼンベルク項や,異方性を表す項などを付け加えたものになっているが,フラストレーションを持つスピン系であるために,通常の量子モンテカルロ法では基底状態の信頼できる計算が不可能である.我々は,このモデルをテンソルネットワーク法などで取り扱った.[3] このモデルでは,相互作用が最近接だけでないために,テンソルの最適化に際して,特別な更新手法を工夫した.(図1)これによる計算の結果,三方晶ゆがみの大きさを自由パラメータとした1次元相図を描いた.(図2)この結果,Na2IrO3に相当する点においては,図1に示されるようなパターンをもった磁化秩序相であることが示された.これは,小さなクラスタの厳密対角化計算による先行研究の結果を支持する.Na2IrO3に相当する点から離れた領域においては新たに3つの相の存在を見出した.第1に,3分の1の六角形において互いに120度の向きにそろった磁化をもつ120度磁気秩序相,第2に,16サイトで1つのユニットセルを作る磁気秩序相,第3に格子の周期とは,非整合な周期をもつ非整合磁気秩序相である.

この結果は,A2IrO3型の磁気化合物が多彩な磁気秩序相を示すであろうことを示唆しており,スピン液体状態の探索にも有用な知見を与えるものである.また,従来は,小規模系の厳密対角化で我慢するしかなかったフラストレート量子スピン系の研究にテンソルネットワーク法が加わったことによって,今後多くの未解決問題に信頼できる解が与えられるものと期待される.

図1:Na2IrO3 のテンソルネットワーク計算に用いられた虚数時間発展の扱い.次近接相互作用があるために特別な取り扱いが必要になる.(文献[3]から転載)
図2 :第一原理計算から導かれたNa2IrO3の有効モデルについてテンソルネットワーク法(iPEPS)などによる計算でえられた相図.(文献[3]から転載)Na2IrO3はΔ=-28meVに相当すると予想されている.

図3  : Na2IrO3の磁気秩序相の模式図

(by 川島 直輝)

[参考文献]

[1] G. Jackeli and G. Khaliullin, Phys. Rev. Lett. 102, 017205 (2009).

[2] Y. Yamaji, Y. Nomura, M. Kurita, R. Arita, and M. Imada, Phys. Rev. Lett. 113, 107201 (2014).

[3] Tsuyoshi Okubo, Kazuya Shinjo, Youhei Yamaji, Naoki Kawashima, Shigetoshi Sota, Takami Tohyama, and Masatoshi Imada, Phys. Rev. B 96, 054434 (2017)

2次元層状ランダムイジング模型の相境界

相互作用が確率的に分布したスピングラス模型は,1975年に提案された Edwards-Anderson 模型を起源とし,ランダムネスとフラストレーションが磁気構造に与える影響について古くから盛んに研究されている.現在では,2次元ランダムボンドイジング模型の相転移点が Kitaev による量子トーラス符号の誤り訂正限界とも関連付けられるなど,統計力学以外の分野からも興味が持たれている.スピンの向きに対する熱平均と相互作用の強さに対する配位平均の2重の平均操作があるため,解析的な計算は困難であり,一般のランダムボンドイジング模型の厳密解は知られていない.一方,McCoy-Wu 模型や Shankar-Murthy 模型など正方格子上で層状に乱れた相互作用を持つイジング模型では,並進対称性より転送行列が簡単に表現できるため厳密解を導出することができる.我々は Shankar-Murthy 模型を拡張した2次元層状ランダムイジング模型を考察し,絶対零度における相転移点を厳密に導出することに成功した.これはランダムボンドイジング模型における非自明な絶対零度相転移点を解析的に導出した初めての例となっている.

図1 : (a) Shankar-Murthy 模型,(b) 周期2の層状乱れがある正方格子,(c) 層状乱れのある六角格子.赤いボンドが反強磁性相互作用,それ以外は強磁性相互作用を表している.

我々が考察した模型は,図1bのように水平方向に2つずらすと元にもどる並進対称性を持っており,水平方向の相互作用のみランダムで垂直方向は強磁性相互作用になっている.この模型は,周期が1の場合に相当する Shankar-Murthy 模型 (図1a) の自明な拡張となっている.垂直方向に進むの転送行列を並進対称性によりブロック対角化し,絶対零度におけるリヤプノフ指数と1次元ランダム磁場イジング模型の基底エネルギーの対応関係に着目することで,リヤプノフ指数の最大値を与えるブロックが変わる点として相転移点を導出した.確率 \(p\) で強磁性相互作用,確率 \(1−p\) で反強磁性相互作用をとる確率分布の場合,得られた強磁性相と常磁性相の絶対零度相転移確率は,正方格子と六角格子 (図1c) でそれぞれ
$$p_c^{(\text{SQ})}=\frac{\sqrt{5}−1}{2}=0.6180⋯$$
$$p_c^{(\text{HEX})}=1−2\sin\frac{\pi}{18}=0.6527⋯$$
である.

図2 : リヤプノフ指数の数値計算により得られた各模型における相境界.横軸が強磁性相互作用の確率,縦軸は温度である.

また,有限温度におけるリヤプノフ指数は,指数関数よりも早く収束するBaiの手法等を用いて数値計算することが可能である.これらの手法を用いて有限温度におけるリヤプノフ指数の交点を計算し,全領域の確率-温度相図 (図2) を高精度に決定した.

[参考文献]

[1] R. Shankar and G. Murthy: Phys. Rev. B 36, 536 (1987).
[2] S. Morita and S. Suzuki: J. Stat. Phys. 162, 123 (2016).

(by 森田 悟史)

複数列表現SU(N)ハイゼンベルグ模型の基底状態

絶対零度でも磁気秩序や対称性の破れを生じない量子スピン液体の探索は現代の物性物理学における大きな研究対象となっている。量子スピン液体を実現する方法の一つとして、スピンの高い対称性により量子揺らぎを増大することが考えられる。SU(2)対称性を持ったスピンが相互作用する通常の反強磁性ハイゼンベルグ模型において、スピンの対称性をSU(N)に拡張したSU(N)ハイゼンベルグ模型では、Nが十分に大きな極限では、バレンスボンド固体(VBS)と呼ばれる、磁気秩序が無いまま格子の並進対称性や回転対称性が破れた状態が実現することが知られており[2]、Nが小さい状況での反強磁性秩序と、Nが大きいVBS秩序との間に、スピン液体状態が存在する可能性について議論が行われている[3,4]。

我々は、このSU(N)反強磁性ハイゼンベルグ模型において、SU(N)のNの値や表現を変えた場合の基底状態について考察した。ここで、「表現を変える」ことは、通常のSU(2)スピンにおいて、スピンの大きさSを変えることに相当しており、SU(N)表現の列の数を n とすると、S=n/2 の関係が成立している。 Nが大きい領域での計算(1/N展開)では、VBS状態のパターンは、nを4で割った余りで決まり、余りが1,3の場合には、Columnar VBS(図1(a))と呼ばれる秩序、余りが2の場合には、Nematic VBS(図1(b))と呼ばれる秩序、余りが0の場合には、S=1の一次元スピン系のハルデーン相と同様に、対称性の破れがないVBS状態が実現すると期待されている[2]。Nがほどほどの大きさでの基底状態については、n=1の場合に量子モンテカルロ法による数値計算から、Nが4以下で反強磁性状態、5以上でColumnar VBS状態となって中間相がないとが確認されている[3,4]一方で、nが2以上の場合については、反強磁性秩序とVBS秩序の相境界については明確な結論が得られていなかった。

図1 : バレンスボンド固体(VBS)状態の模式図。破線、実線、濃い線の順に、スピン間の相関が強い。(a) Columnar VBS、(b) Nematic VBS。[1]より。
本研究では、並列化した量子モンテカルロ法を用いてこれまでよりも大きな系での計算を行うことで、SU(N)表現の列の数 n が2の場合には、Nが9以下では反強磁性秩序、Nが10以上ではNematic VBS秩序が基底状態となり、中間相が存在しないことを明確に示した(図2)。一方、n=3の場合には、Nが15以上で反強磁性秩序が消失する一方で、明確なVBS秩序は観測されず、中間相の存在を否定する結果は得られなかった。ただし、VBS秩序パラメタのサイズ依存性と 1/N 展開との比較から、今回の計算したシステムサイズ L128では強い有限サイズ効果に隠れて秩序が観測されないことが示唆されるため、今回の結果は、直ちに中間相の存在を意味するものではない。将来、より大きな系( L500)での計算が可能になれば、n が3以上での中間相の有無がはっきりすると期待される。

図2 : 複数列表現のSU(N)ハイゼンベルグ模型の基底状態相図。横軸はSU(N)のNの大きさ。縦軸は、表現の列数(SU(2)でのスピンの大きさSに相当)。[1]より。

(by 大久保 毅)

[参考文献]

[1] T. Okubo, K. Harada, J. Lou and N. Kawashima: Phys. Rev. B 92, 134404 (2015).
[2] N. Read and S. Sachdev: Phys. Rev. B 42, 4568 (1990), N. Read and S. Sachdev: Nucl. Phys. B 316, 609 (1989).
[3] K. Harada, N. Kawashima and M. Troyer: Phys. Rev. Lett. 90, 117203 (2003).
[4] N. Kawashima and Y. Tanabe: Phys. Rev. Lett. 98, 057202 (2007).

危険なイレレバント場に対するスケーリング関係式

連続相転移で現れる臨界現象では,通常,「イレレバント」な性質は大きな役割を果たさずに「レレバント」な性質により臨界指数などの特徴が支配される.そのため,臨界現象は物理系の詳細に依存せずに,物理系の対称性や空間の次元等の少数の性質だけで,その振る舞いの特徴が決まっている.しかし,特殊な場合には,臨界現象に重要な役割を果たすイレレバント変数が存在し,それは「危険なイレレバント場」と呼ばれている.我々は,そのような危険なイレレバント変数が物理系の対称性を低下させる場合について考察し,臨界現象を特徴付ける臨界指数の間に新しいスケーリング関係式が成立することを示した.

図1 : 危険なイレレバント場 \(\lambda\) が存在する場合の一般的なくり込みの図.[1]より.
我々は図1のような二つの固定点(臨界点を支配する臨界固定点 \(X\) と,秩序状態を支配する固定点 \(Y\))を含む一般的なくり込み群の流れを考察した.このくり込みの流れでは,対称性を破る場 \(\lambda\) は,臨界点近傍では,ある長さスケール \(\xi\) まででほとんどゼロになる一方で,もう少し長いスケール \(\xi^\prime\) では,その値が有限に「復活」する.そのため,二つの長さスケールの中間の大きさ, \(\xi \ll L \ll \xi^\prime\) では, 対称性を破る場 \(\lambda\) があたかも無い(\(\lambda=0\))ように感じられ,新しい対称性が現れたように見えることになる.

このようなくり込みの流れが実現している古典的な例は,3次元のq状態クロック模型と呼ばれる,q 個の向やすい方向を持った,平面内で自由な方向を向けるベクトルスピンが相互作用する模型である.この模型において,短い相関長 \(\xi\) を特徴付ける臨界指数 \(\nu\) と長い相関長 \(\xi^\prime\) を特徴付ける臨界指数 \(\nu^\prime\) との間に成立するスケーリング関係式については,これまでにドメインウォールの存在を仮定し,その自由エネルギーが系の体積や断面積に比例するすると考える議論が行われていた[2,3,4].我々は,くり込みの流れ(図1)を元にしたより一般的な議論により,二つの固定点 \(X, Y\) での \(\lambda\) のスケーリング次元,\(y_\lambda < 0\) と \(y_\lambda^\prime > 0\) を用いて,\(\nu\) と \(\nu^\prime\) との間に,$$\frac{\nu^\prime}{\nu} = 1 – \frac{y_\lambda}{y_\lambda^\prime}$$という関係式が成立することを示した.

図2 : 平均磁化の方向に依存したスピン揺らぎのフーリエ変換 \(\tilde{S}_k\) の温度依存性の両対数プロット.(a) 3次元 6状態 クロック模型 (b) 4次元 6状態 クロック模型.[1]より.
また,提案したスケーリング関係式が成立していることを確認するために,平均磁化の方向に依存したスピン揺らぎのフーリエ変換 \(\tilde{S}_k\)を定義しその振る舞いをモンテカルロシミュレーションにより計算した.この \(\tilde{S}_k\) は有限サイズスケーリング $$\tilde{S}_{k} \sim L^\mu g \left [(T-T_c)L^{\nu^\prime}\right]$$ に従うことが期待される.図2 に示したように,数値計算結果は,期待される有限サイズスケーリングによく従っており,我々が提案したスケーリング関係式が確かに成立していることが確認できた.

(by 大久保 毅)

参考文献

  1. T. Okubo, K. Oshikawa, H. Watanabe and N. Kawashima: Phys. Rev. B 91, 174417 (2015).
  2. Y. Ueno and K. Mitsubo: Phys. Rev. B 43, 8654 (1991).
  3. M. Oshikawa: Phys. Rev. B 61, 3430 (2000).
  4. J. Low, A. W. Sandvik, and L. Balents: Phys. Rev. Lett. 99, 207203 (2007).

気泡のOstwald成長過程におけるスケーリング則

シャンパンや炭酸飲料の栓をあけるとたくさんの泡が出るが、その後、大きい気泡がより大きく、小さい気泡がより小さくなるOstwald成長という現象が起きる。 Ostwald成長は、合金系や液滴生成など、一次転移を引き起こす系において広く見られる普遍的な現象であり、その成長則は1960年代にLifshitz, Slyozov, Wagnerらによって完成され、LSW理論という古典論としてまとめられている[1,2]。 しかし、気泡生成現象、特に液相と気相が同じ物質である単一成分系のOstwald成長においては実験が難しく、 数値計算においても液滴生成に比べて気泡生成は要求される計算量が莫大になることからこれまで研究が進んでいなかった。 特にLSW理論で仮定されている自己スケーリング則と平均場的描像については、気泡生成系で成り立つかどうかは非自明である。

図1:減圧シミュレーションにおける気泡の可視化。理研の稲岡氏に依る。

そこで我々は分子動力学法を用いた急減圧シミュレーションにより気泡分布関数の時間発展を調べた。 これまでの分子動力学法を用いた先行研究においては、泡は数個から10個程度しか表現できず、 分布関数の評価は困難であったが、京コンピュータ上で数億粒子規模の計算を行うことで、気泡分布関数を直接評価することが可能となった。 その結果、気泡分布関数が確かに自己スケーリング則に従っていること、スケーリング指数が低温では3/2、 高温では1になることがわかった。これは、気泡の成長過程が低温では界面律速であるが、高温では拡散律速にクロスオーバーすることに対応している。本研究により、気泡生成系においても自己スケーリング則が成立しており、かつ指数が古典論による予測と 一致することがわかった。これは、気泡生成過程において圧力の緩和が十分に早いため、個々の気泡が感じる圧力が一様であり、 結果として平均場描像が正当化されることを意味する。

図2:気泡数の時間発展の温度依存性。いずれも急減圧の直後に増えた後、ベキ的に減衰するが、温度が低い場合の指数は3/2、高い場合には1となる。

LSW理論は、準安定領域における核生成率を予想する古典核生成論と多くの仮定を共有するが、 古典核生成論が気泡生成率の予測に失敗するのに対し、LSW理論は気泡成長を良く記述することがわかった。 古典核生成論が気泡生成率を記述できない理由については、今後の研究課題として残されている。

(by 渡辺宙志)

参考文献

[1] I. Lifshitz and V. Slyozov, J. Phys. Chem. Solids 19, 35 (1961).
[2] C. Wagner, Z. Elektrochem. 65, 581 (1961).
[3] H. Watanabe, M. Suzuki, H. Inaoka, and N. Ito, J. Chem. Phys. 141 234703 (2014).

AIPによるプレスリリース: How the Physics of Champagne and Soda Bubbles May Help Address the World’s Future Energy Needs

非局所更新による並列化量子モンテカルロアルゴリズム

量子多体系の大規模計算によって解明が期待できる問題は数多く残されている。 近年のハイパフォーマンス・コンピュータは、京コンピューターに代表されるように、多数コアによる大規模並列により計算性能を稼ぐものが主流であり、大規模計算を行う場合、アルゴリズムの並列化が有効な手段として挙げられる。ファインマンの経路積分表示に基づく世界線量子モンテカルロ法は比較的大きなサイズの量子多体系を統計誤差の範囲内で精密に計算できる数値計算手法であるが、中でもワームアルゴリズム[1,2]は、大域的な更新ができる、汎用性の高い優れたアルゴリズムである。 しかし、ワームアルゴリズムでは、空間の1点に着目し、そこに挿入されたワームと呼ばれるオペレーターがループを描きながら世界線の状態を更新していく、イベント駆動型であることから、並列化が非自明である。

我々は向き付きループアルゴリズム(DLA)[2]に基づき、複数のワームを配位空間中に挿入でき、配位空間をドメイン分割できる並列計算向きのアルゴリズム(Parallelized Multi-Worm Algorithm = PMWA)を考案した(図1)。複数ワームはモデルハミルトニアンにソース場ηを付加し調節することにより導入され、測定された物理量をη=0の極限に外挿することによってソース場のない時の結果が得られる。ワームを複数扱う際には、これまでのワームアルゴリズムと同じ手続きでは時間反転対称性を満たすように状態更新がでないが、我々のアルゴリズムはこれを満たすようにできている。また、ドメイン間通信によりワームのドメイン間移動とドメイン境界の状態更新が有効的に果たされ、エルゴード性も保たれている。この際の情報転送量はドメイン表面積にしか依存しないため、プロセッサ並列数が増えても高い並列化効率を維持することができる。 このアルゴリズムはこれまでのワームアルゴリズム同様、ソフトコア・ボーズ粒子系や、不符号問題のない量子スピン系に適用できる。

図1:PMWAでのワームがいるの配位空間。

我々はこのアルゴリズムを正方格子上拡張ハードコア・ボーズ・ハバードモデル

に適用し、ベンチマーク計算を行った。 ここでiサイトのボゾンの消滅(生成)演算子を bi (bi†)とし、 tは隣接サイト間の粒子のホッピングエネルギー、Vは隣接サイトの粒子間相互作用、μは化学ポテンシャルをそれぞれ表している。我々は物理量のηに関する外挿則を導き、DLAの結果と比較することで、物理量がそれに従い振る舞うことを確認した。 図2(a)は超流動相でのボーズ・アインシュタイン凝縮秩序変数Qをηの関数としてプロットしたものである。ここではLを一方向あたりの全格子数、βを逆温度としたとき、最大でL×L×β =10,240×10,240×16のサイズの計算を3,200プロセッサ並列により実現した。 これはシングルプロセッサの計算で実行できるサイズを大きく上回っている。 また、我々はアルゴリズムの精度を評価するために、標準誤差のドメイン分割数Nの依存性を調べた。 分割による緩和速度の低下は非常に小さいことがわかった(図2(b))。 また、ストロングスケーリングでの並列化効率が非常に良いために、同じ実時間内で従来のアルゴリズムと比較した場合、 N>8で従来のアルゴリズムの精度を上回り、Nを上げるにつれさらに精度が向上することを確認できた(図2(c))。

図2:(a)様々なシステムサイズにおけるボーズ・アインシュタイン凝縮秩序変数のη依存性。 (b)モンテカルロステップ数を固定したときの標準誤差の分割数依存性。 (c)実時間を固定したときの標準誤差の分割数依存性。[3]
(by 正木 晶子)

参考文献

[1] N. Prokof’ev, B. Svistunov and I. Tupitsyn, Sov. Phys. JETP 87, 310 (1998) 等.
[2] O. F. Syljuasen and A. W. Sandvik, Phys. Rev. E 66, 046701 (2012).
[3] A. Masaki-Kato, T. Suzuki, K. Harada, S. Todo and N. Kawashima, Phys. Rev. Lett. 112, 140603(2014).