アニーリングマシンを用いたパターン形成
向井かのん

1

(指導教員：工藤和恵)

はじめに2.3

モンテカルロ法を用いて，初期状態から状態遷移を繰り返し，与えられた分布に近づけていく．ある状態を X ，与えられた分布関数を P (X) として，遷移確率を次式のように定める．
)(P (X ′ )′W (X → X ) = min 1,
P (X)

近年，アニーリングマシンの普及が徐々に進んできており，主に最適化問題に関する多くの研究が進められている．これに対して本研究は，アニーリングマシンのさらなる発展を目指し，最適化とは全く別の領域である「時間発展」に応用することを考え，検証したものである．
今回はパターン形成を時間発展の題材に取り上げ，古典的な方法とアニーリング法とでシミュレーション結果を比較した．パターンとは例えばトマトの表面のシワやシマウマの表面の模様のことであり，パターン形成はそのパターンのできるメカニズムを解明する分野である．

2

モデルと方法

2.1

モデル

モンテカルロ法

X は正方格子の白黒 (イジング変数) の組み合わせである．この確率にしたがって状態を遷移させることで与えられた分布に近づけることができる．磁性体のモデルは熱力学的にギブス分布に従うことが知られているので，T を温度，E(X) をエネルギー関数，Z(T ) を分配関数として，P (X) を次のように定める．
P (X; T ) =強磁性体モデルを用い，イジングクラスターを急冷したときのドメイン（白や黒の塊）の大きさの成長を時間をおって観察する．図 1 のように N × N の正方格子を用意する．図 1(a) は初期状態で，それぞれの格子点に白または黒をランダムに割り当てる．図 1(b) は時間経過後のパターンを表す．白はイジング変数 1，黒はイジング変数 −1 に対応する．ハミルトニアン (エネルギー)H は次のようになる．

H=∑

ただし，今回は急冷するという設定なので温度 T は低温で一定とする．

3仮説

イジング模型の連続化に対応する Ginzburg-Landau
自由エネルギーは次のように書ける．
[]∫D
F = d2 r f (m) + (∇m)2(1)2JSi Sj⟨i,j⟩m(r) はイジング変数に対応する秩序パラメタ，f (m)は m = ±1 で最小値をとる関数，D は空間勾配を減らす効果の強さを表す定数である．F の変分は次のようになる．
∫[]
δF= d2 r f ′ (m) − D(∇2 m)(2)δm

J は格子点間の相互作用，Si は格子点 i に対応するイジング変数である．隣り合う格子点 i, j の組みに関して和をとる．
(a)1exp(−E(X)/T )
Z(T )

(b)

ここから時間発展を考えると次のようになる．

∂mδF= −λ∂tδm

図 1: (a) 初期状態．(b) 時間経過後．

λ は変化の速さに対応するパラメタである．式 (3)に 式 (2) を 代 入 す る と ，次 の 時 間 依 存 GinzburgLandau(TDGL) 方程式が得られる．

今回は非保存系 (白と黒の格子点数がそれぞれ変化する) であるため，Glauber Dynamics に従う．この場
1合，ドメインの大きさは時間を t として t 2 に比例することが知られている．

2.2(3)[]∂m= −λ f ′ (m) − D(∇2 m)∂t

シミュレーションの目的と方法(4)

式 (4) を後退差分の形で表すと次のようになる．

次節の仮説「時間発展が各時刻での最適化問題を解くことで計算できる」を実証することが目的である．
そのためにまず既存のシミュレーション方法であるモ1ンテカルロ法を用いて，ドメインの大きさが t 2 に比例することを確認する．次に，これと同じことをアニーリングマシンを用いて再現できれば，アニーリングマ
シンを用いて時間発展問題を解くことができることを示すことができる．

[]
mt − mt−1= −λ f ′ (mt ) − D(∇2 mt )∆t(5)ここで，添字は時間ステップを表す．
ここから，時刻 t での自由エネルギーを Ft とし，次のような量 H を最適化することを考える．
∫∫
112d rmt−1 mt +d2 rm2t (6)H = λFt −
∆t2∆t

17

mt の微小変化に対する H の微小変化は次のようになる．
∫∫
1
12d rmt−1 δmt +d2 rmt δmt
δH = λδFt −∆t∆t(7)δH = 0 が最適化の必要条件となるはずなので，それを式 (7) に代入すると次式が得られる．
{}∫[ ′] mt − mt−122d r λ f (mt ) − D(∇ mt ) +δmt = 0∆t(8)よって，次の式が得られる．
[]
mt − mt−1= −λ f ′ (m) − D(∇2 m)(9)∆t

図 3: 3 MCS おきのスナップショット
シン側にはいくつかのパラメタがあるが，ドメイン成長に大きな影響を与えるのは初期温度であることがわかった．初期温度を調整し，アニーリング法で行った結果が図 4 である．

式 (5) と式 (9) は同一のものである．よって，時間発展方程式を解くことは各時刻で H を最小化することと同じである．H を離散化し，イジング模型の形で表すと次のようになる．
∑∑
H=JSit Sjt −hi Sit(10)⟨i,j⟩iSit は時刻 t における格子点 i に対応するイジング変数で，J = −λD/(∆x)2 , hi = Kt Sit−1 = (−1/∆t)Sit−1とした．

4結果

4.1

モンテカルロ法

図 4: アニーリング法によるドメイン長の推移
ドメイン長の線が青線とほぼ傾きが同じであること
1から，ドメイン長は t 2 に比例して増大することがわ
かった．パターンは図 5 のように遷移した．

既存のシミュレーション方法，モンテカルロ法で行なった場合のドメイン成長は図 2 のようになった．x軸，y 軸共に log スケールで表されており，赤線がド1メインの平均的な大きさ，青線が t 2 である．赤線がほぼ直線となっており，傾きが青線と大幅同じである
1ので，ドメイン長は t 2 に比例して増大していると言える．

図 5: 10 ステップおきのスナップショット

5

まとめと今後の課題

モンテカルロ法とアニーリング法とで同じような結果を得られたことから，前述の仮説，つまり時間発展が各時刻での最適化問題を解くことで計算できることが示された．アニーリングマシンの違いや対象とする問題によって適切なチューニングを行えば，今回のドメイン成長に限らず様々な問題を扱うことが可能であると思われる．
しかしながら，4.2 節の結果は初期温度を調整しないと得られない．その原因を追求することが今後の課題である．
図 2: モンテカルロ法によるドメイン長の推移

参考文献[1] 福島 孝治：“モンテカルロ法の前線―サイコロ振って積分する方法―”，若手研究者・学生向けに最新技術をわかりやすく紹介する講演会「確率的アルゴリズムによる情報処理」講義ノート，(2003).

パターンの状態は図 3 のように遷移した．1 ステッ
プを 1 MCS(モンテカルロステップ，平均して全ての格子点が反転する) を 1 ステップとし，3 ステップおきに表示してある．

4.2

アニーリング法[2] Gary S. Grest, David J. Srolovitz :“Structure and Evolution of Quenched Ising Cluster”，
Phys．Rev．B 30，5150 (1984).

アニーリングマシンに 1 回計算させて結果を受け取るまでを 1 ステップとした．
試行を重ねた結果，J と Kt の比は 1.0 : 2.0 のときのみ適切に成長することがわかった．アニーリングマ

18