ブックタイトル日本結晶学会誌Vol62No1

概要

日本結晶学会誌Vol62No1

X線回折実験とシミュレーションのデータ同化による結晶構造解析らかのペナルティー汎関数D［Iexp（θ），Icalc（θ；R）］を導入し，p（I exp（θ）|R）＝exp｛－γD［I exp（θ），I calc（θ；R）］｝のような形で表すことにする（γはある定数）．データ同化においては，I exp（θ）は実験によりすでに確定した量であるので，D［I exp（θ），I calc（θ；R）］は，Rだけの関数であることに注意する．この尤度も凸関数ではなく多数の極大をもつ．以上のように定義した事前確率と尤度をベイズの定理に代入すると，Rに関する事後確率p( R| I (θ)) ? exp{ ?βE( R) ?γD[ I (θ), I (θ; R)]}exp exp calc= exp{ ?β[ E(R) +αND( R )]}（3）が得られる.ここで，E（R）とD（R）（定義は以下で与える）がそれぞれ，O（N），O（1）の量であることを考慮に入れ，重み係数α＝γ/βNを導入した．この事後確率を最大にするには，F（R）≡E（R）＋αND（R）を最小にするRを求めればよい．これを最大事後確率（MAP）推定と呼ぶ．事前分布と尤度は，ともに正しい結晶構造で最大値をとる．それ以外にも両者は多数の極大値をもつが，事前分布（＝全エネルギー）は近距離の構造を，尤度（＝実験との一致度）は長距離の周期構造を反映した量であるため，局所最適点の分布はまったく異なると期待できる．したがって，それらを掛け合わせた事後分布では，正しい結晶構造に対応した大域最適点のみが強調され，より効率的に正しい結晶構造を探索することが可能となる（図1）．次に，実験データとの一致度を表すペナルティー汎関数をどう構成するかについて考えよう．これまで，R因子［式（1）］が広く用いられてきたが，実験データが不完全な場合にもロバストに使える指標として，以下で定義される「結晶化度」λ（R）を導入する．λ( R)=θobs+ ?∑∫Iθcalcθobs ??obs∫Icalc(θ; R)dθ(θ; R)dθ（4）ここで，分母は原子位置Rから計算された回折強度の和，分子はこの回折パターンのうち，実験で観測されているピークと角度が一致するものだけについての和になっている．結晶化度の計算には，実験で得られたピークの位置のみを使い，その強度情報は利用しない．そのため，実験データのノイズやバックグラウンドの影響をうけず，実験データが不完全な場合においても，安定した解析が可能となる．計算により得られたピークの位置が実験のものと完全に一致する場合にはλ（R）＝1，一致しないものがある場合には1よりも小さな非負の値をとる．このλ（R）を用いて，ペナルティー汎関数をD（R）＝1－λ（R）と定義する．図2に，SiO 2多形の1つであるコーサイトに対してデータ同化結晶構造予測手法を適用した結果を示す．コーサイトの空間群はC2/cである．結晶の対称性は非常に低日本結晶学会誌第62巻第1号（2020）Crystallinity1.00.80.60.40.2Optmized with D coesiteOptmized without D coesite0.00 1000 2000 3000 4000Time (fs)図2コーサイトに対する構造最適Structure化．（optimizationfor coesite.）全エネルギーE（R）のみを最適化した場合（黒線）には，結晶化度［式（4）］の値は0.4程度と低い値に留まる．一方で，X線回折ピークの実験との一致度D（R）と組み合わせたハイブリッドコスト関数を最適化すると，結晶化度はほぼ1まで上昇し，実験結果と一致する長周期構造が実現されていることがわかる．く，単位格子に48原子を含む複雑な構造をもっている．計算には96原子を含む，2×1×1のサイズの超格子を用い，温度5000 Kからのシミュレーテッドアニーリングにより最適化を行った．全エネルギーE（R）のみを最適化した場合（黒線）には，結晶化度［式（4）］の値は0.4程度と低い値に留まっており，局所的に正しい構造が得られているものの，正しい周期をもたないアモルファス的な構造に捕らわれてしまっていることが見てとれる．一方で，X線回折ピークの実験との一致度D（R）と組み合わせたハイブリッドコスト関数を最適化すると，結晶化度はほぼ1まで上昇し，最終的に実験結果と一致する長周期構造が得られていることがわかる．実際の構造予測精度（図3）は，重み係数αの選び方に大きく依存しており，αの値が小さすぎても大きすぎても予測精度は低下する．局所最適と大域最適との間のエネルギーバリアを超えるには，コスト関数の第2項αND（R）の大きさが，バリアの高さの数倍程度になるようαの値を調整する必要がある．バリアの高さは，結晶の融点をTmとすると，―3 2 Nk BT m程度と見積もることができる．一方，局所最適構造におけるD（R）＝1－λ（R）の大きさは典型的には2/3程度である（図2）．この両者が釣り合うのは，―4α9kB＝Tmのときである．われわれのシミュレーションでは，―4α9kBをTmの4倍程度に選んだ場合に，予測効率が最大となっており（図3），上記の見積もりが正しいことが確認できる．図4に―4α9kB＝4T mの場合の最適化の様子を示す．実験データとの一致度と組み合わせたハイブリッドコスト関数は，ほぼ単調に減少しているが，全エネルギーはいったん上昇し，その後正しい結晶構造にたどり着いていることがわかる．このように，正しい重み係数を選択することで，局所最適から抜け出し，高い確率で大域最適点にたどり着くことが可能となる．53