專利名稱：模擬方法及程序的制作方法
技術(shù)領(lǐng)域：
本發(fā)明涉及模擬方法及程序，尤其涉及利用了分子動力學(xué)的模擬方法及用于使其在計算機中執(zhí)行的程序。
背景技術(shù)：
利用了分子動力學(xué)的計算機模擬正在不斷施行。分子動力學(xué)中，構(gòu)成作為模擬對象的系統(tǒng)的粒子的運動方程式能夠以數(shù)值方式解出。若作為模擬對象的系統(tǒng)的粒子數(shù)增加，則所需的計算量也增加。以目前計算機的運算能力，例如只能進行到粒子數(shù)為10萬左右的系統(tǒng)的模擬。為了減少模擬所需的計算量，例如在非專利文獻1中，嘗試了在分子動力學(xué)中應(yīng)用重正化方法。另外，本申請所涉及的發(fā)明由本申請發(fā)明人完成，并公開于專利文獻1中?，F(xiàn)有技術(shù)文獻專利文獻專利文獻1 日本特開2006-285866號公報非專利文獻1 稻村，武澤，社本，《關(guān)于基于重正化方法的可變標(biāo)度分子動力學(xué)模擬》(「繰D込 O方法(二基d〈可變7 分子動力學(xué)* ^工>一* 307」)， 日本機械學(xué)會論文集(A篇)，1997年，第63卷，608號，P. 202-207發(fā)明的概要發(fā)明要解決的課題如非專利文獻1中專欄“3.數(shù)值計算法”所記載，在非專利文獻1的模擬方法中無法明確評價溫度。

發(fā)明內(nèi)容
本發(fā)明的目的之一在于提供一種利用分子動力學(xué)的新型模擬方法及用于在計算機中執(zhí)行該方法的程序。本發(fā)明的另一目的在于提供一種可根據(jù)重正化群方法明確評價溫度的模擬方法及用于在計算機中執(zhí)行該方法的程序。用于解決課題的手段根據(jù)本發(fā)明的一種觀點，提供一種模擬方法，其方法具有如下工序(a)工序，考慮包含N個粒子、各粒子質(zhì)量為m、且能夠以表示對粒子間距離的依賴性的無因次函數(shù)f與相互作用系數(shù)ε之積ε f來表示粒子間的相互作用勢能的粒子系統(tǒng)S時，使用配置有該粒子系統(tǒng)S的空間的維數(shù)d確定大于1的第1重正化因子α、0以上d以下的第2重正化因子 Y及0以上的第3重正化因子δ，并通過轉(zhuǎn)換式N' =N/Cid求出重正化的粒子數(shù)N'，通過轉(zhuǎn)換式m' =11^8/(^求出重正化的粒子質(zhì)量111'，通過轉(zhuǎn)換式ε' = ε α γ求出重正化的相互作用系數(shù)ε ‘；及(b)工序，對包含重正化的粒子數(shù)為N'個的粒子、各粒子具有重正化的粒子質(zhì)量m'、且由所述無因次函數(shù)f與重正化的相互作用系數(shù)ε ‘之積ε ‘ f 表示粒子間的相互作用勢能的粒子系統(tǒng)S'執(zhí)行分子動力學(xué)計算。另外，根據(jù)本發(fā)明的其他觀點，提供一種模擬方法，其方法具有如下工序(a)工序，考慮包含N個粒子、各粒子質(zhì)量為m、且能夠以表示對粒子間距離的依賴性的無因次函數(shù)f及相互作用系數(shù)ε之積來表示粒子間的相互作用勢能的粒子系統(tǒng)S時，使用配置有該粒子系統(tǒng)S的空間的維數(shù)d確定大于1的第1重正化因子α和0以上的第3重正化因子S，并通過轉(zhuǎn)換式N' =N/Cid求出重正化的粒子數(shù)N'，通過轉(zhuǎn)換式m' =mas求出重正化的粒子質(zhì)量m'；及(b)工序，對包含重正化的粒子數(shù)N'個的粒子、各粒子具有重正化的粒子質(zhì)量m'、且由所述無因次函數(shù)f與所述相互作用系數(shù)ε之積表示粒子間的相互作用勢能的粒子系統(tǒng)S'執(zhí)行分子動力學(xué)計算。并且，根據(jù)本發(fā)明的其他觀點，提供一種模擬方法，其方法具有如下工序(d)工序，考慮包含N個粒子、各粒子質(zhì)量為m、且能夠以表示對粒子間距離的依賴性的無因次函數(shù)f與相互作用系數(shù)ε之積ε f來表示粒子間的相互作用勢能的粒子系統(tǒng)S，并在配置有所述粒子系統(tǒng)S的空間上考慮XYZ正交坐標(biāo)系時，分別使用大于1的X方向的重正化因子 αχ、Y方向的重正化因子方向的重正化因子αζ，通過轉(zhuǎn)換式N' = N/( α χ α γ α ζ) 求出重正化的粒子數(shù)N'，并且確定0以上的重正化因子δ，通過轉(zhuǎn)換式％' =maxs求出 X方向上的重正化的粒子質(zhì)量％'，通過轉(zhuǎn)換式m/ = ma /求出Y方向上的重正化的粒子質(zhì)量％'，通過轉(zhuǎn)換式mz' = ma /求出Z方向上的重正化的粒子質(zhì)量mz'；及(e)工序， 對包含重正化的粒子數(shù)為N'個的粒子、各粒子具有重正化的粒子質(zhì)量m'、且由所述無因次函數(shù)f與所述相互作用系數(shù)ε之積表示粒子間的相互作用勢能的粒子系統(tǒng)S'，在 X方向的運動方程式中使用X方向上的重正化的粒子質(zhì)量％'執(zhí)行分子動力學(xué)計算，在Y方向的運動方程式中使用Y方向上的重正化的粒子質(zhì)量％'執(zhí)行分子動力學(xué)計算，在Z方向的運動方程式中使用Z方向上的重正化的粒子質(zhì)量mz'執(zhí)行分子動力學(xué)計算。發(fā)明的效果重正化的粒子數(shù)N'少于粒子數(shù)N。因此，能夠以比對粒子系統(tǒng)S執(zhí)行分子動力學(xué)計算時更少的計算量執(zhí)行對粒子系統(tǒng)S'的分子動力學(xué)計算。另外，由于任意地確定0以上的重正化因子δ而進行計算，因此可靈活地執(zhí)行模擬。


圖IA是用于說明標(biāo)度轉(zhuǎn)換法的導(dǎo)出方法的圖。圖IB是用于說明標(biāo)度轉(zhuǎn)換法的導(dǎo)出方法的圖。圖2是表示使用重正化群分子動力學(xué)的對鋁的模擬結(jié)果的圖表。圖3Α是用于說明關(guān)于S-W勢的圖。圖3Β是表示使用重正化群分子動力學(xué)的對硅的模擬結(jié)果的圖表。圖4是表示利用重正化群分子動力學(xué)求出分散關(guān)系的模擬結(jié)果的圖表。符號的說明a -第1重正化因子，N-粒子數(shù)，N'-重正化的粒子數(shù)。
具體實施例方式
首先，對分子動力學(xué)(Molecular Dynamics ；MD)進行簡潔說明。對由N個粒子(例如原子)構(gòu)成、哈密頓算符H表示為如下(1)的粒子系統(tǒng)進行考慮。其中，各粒子質(zhì)量為m，
粒子j的運動量矢量為Pj,粒子j的位置矢量(位置坐標(biāo))為q」，粒子間的相互作用勢能為 φ。
— 2
JV
(1)
N 孖=Σ
；=1
P . Ν
-ζ-+ Σ 椒-q；) 1tjti式。
通過將哈密頓算符H代入哈密頓的正則方程式，得到如下對各粒子的運動方程
dpO^lqi - q .)
=‘ (J = 1，···，Ν)···(2)及
dt ιΦ]
^ =^ = V1 (j = 1”··，Ν)…(3)
dt m
其中，在公式(3)中，\為粒子j的速度矢量。分子動力學(xué)中通過數(shù)值積分并解出由公式(2)及(3)表示的運動方程式來對構(gòu)成粒子系統(tǒng)的各粒子求出各時刻的各粒子的運動量矢量(或者速度矢量)及位置矢量。另外，在數(shù)值積分中多使用Verlet法。關(guān)于 Verlet 法，例如在“J. M. Thi jssen，; Computational Physics，，Cambridge University Press (1999)，，的p. 175中所說明。根據(jù)由分子動力學(xué)計算得到的各粒子的位置、速度能夠算出系統(tǒng)的各種物理量。下面，對應(yīng)用重正化群方法的分子動力學(xué)(將此稱為重正化群分子動力學(xué))進行概括性說明。在重正化群分子動力學(xué)中，使希望求出物理量的期望的系統(tǒng)S對應(yīng)于由少于系統(tǒng) S的粒子數(shù)構(gòu)成的系統(tǒng)(將此稱為重正化的系統(tǒng)S')。其次，對重正化的系統(tǒng)S'執(zhí)行分子動力學(xué)計算。并且，使對重正化的系統(tǒng)S'計算的結(jié)果對應(yīng)于期望的系統(tǒng)S。由此，能夠以比直接對期望的系統(tǒng)S執(zhí)行分子動力學(xué)計算時更少的計算量對期望的系統(tǒng)S算出物理量。 將用于使期望的系統(tǒng)S中的物理量(例如，粒子數(shù)、各粒子質(zhì)量等)與重正化的系統(tǒng)S'中的物理量相互對應(yīng)的轉(zhuǎn)換法稱為標(biāo)度轉(zhuǎn)換法。下面，對本申請發(fā)明人發(fā)現(xiàn)的標(biāo)度轉(zhuǎn)換法進行說明。假定希望算出物理量的期望的粒子系統(tǒng)S由N個粒子構(gòu)成。設(shè)定粒子系統(tǒng)S的所有哈密頓算符H表示為如公式(1)所
7J\ ο
_ 2
孖=Σ
；=1
PN
τ^+ Σ 敝-q；) 1tjti某一粒子j的哈密頓算符&表示為如下。
D 2 N
Η^^Ζ+ΣΦ^-^'·· (4)
^m i = j+\
為了討論重正化，通過表示對粒子間距的依賴性并被無因次化的函數(shù)f、與表示相互作用的強度并具有能量因次的系數(shù)ε (將此稱為相互作用系數(shù)ε)之積，表示相互作用勢能Φ如下。
φ( = ε
/\ q. -q.Φ\^—(\0=^」——L (5)
V σ J例如，對惰性原子適用如下相互作用位勢。其中，r為原子間距離。相互作用系數(shù) ε及常數(shù)ο取對應(yīng)原子種類的值。下面，對從哈密頓算符H導(dǎo)出重正化的哈密頓算符H'的方法進行說明。如以下詳細(xì)所示，重正化的哈密頓算符H'通過如下方法獲得執(zhí)行對粒子系統(tǒng)S的分配函數(shù)Ζ(β ) 的積分的一部分，粗?；茴D算符，接著對積分變量進行再標(biāo)度。對粒子數(shù)固定的正則系綜(Canonical ensemble)的分配函數(shù)Z ( β )相對于粒子系統(tǒng)S表示如下。Ζ(β) = / drNexp(-^H(p, q))··· (7)其中，系數(shù)β通過系統(tǒng)的溫度T與玻爾茲曼常數(shù)、定義為β = l/(kBT)。d「N為相位空間內(nèi)的體積要素，更詳細(xì)地表示為如下。
1 N1..drN = —= —DNpDNq …(8)
VV N 3VyN其中，h3N(h為普朗克常數(shù))。另外，將DpN及D；1分別定義為如下。
NDno = Yidni
P J=I 1
NDf = Yidai
q j=i 1對諧振子以外的位勢的重正化進行討論較為困難。因此，將具有鞍點的相互作用位勢分成3個區(qū)域來考慮。由于粒子間距r —c 與r — 0為固定點，因此重正化轉(zhuǎn)換時位勢不變。因此，假設(shè)鞍點附近的重正化轉(zhuǎn)換對所有粒子間距r都成立。對鞍點的周圍進行泰勒展開，保留至平方波動。若將鞍點的位置設(shè)為Α，相對變位設(shè)為SU，則Su的一次項消失，得到如下公式。φ{(diào)τ0+δη) = φ(τ0) + ]-Su2- (9)
2 or
L」r=r0繼而，關(guān)于粒子j表述鞍點附近的哈密頓算符如下。
P,2 V^H1 =^+ V
]2m其中，Ui及屮分別表示從粒子i及j的鞍點的變位。φ"為關(guān)于粒子間距r的二階微分。下面，對分配函數(shù)Ζ( β)的粒子間相互作用部分的粗粒化進行說明。首先，對配置于1維鏈上的粒子系統(tǒng)的粗?；M行觀察后，向配置于簡單立方晶格上的粒子系統(tǒng)擴張。如圖IA所示，1維鏈上的格點(lattice point)上相互最鄰接配置粒子i與j，且相互最鄰接配置粒子j與k。粒子i與k的中間配置有粒子j。用實線表示最鄰接的粒子彼此的相互作用(最鄰接連接)?？紤]消去并粗?；Ｗ觠。若寫出粒子j所參與的相互
0
1
Xc
\ /
(ro
+
_
y
U
I
UI
\J
(^o
W"
於
6
σ
r
8作用，并對位于粒子i與k之間的粒子j的變位~執(zhí)行積分，則關(guān)于粒子j的影響被重正化后的、粒子i與k的相互作用位勢φ'(化-qk)得出如下公式。
權(quán)利要求
1.一種模擬方法，其特征在于，具有(a)工序，在考慮包含N個粒子、各粒子質(zhì)量為m、且能夠以表示對粒子間距離的依賴性的無因次函數(shù)f與相互作用系數(shù)ε之積來表示粒子間的相互作用勢能的粒子系統(tǒng) S時，使用配置有該粒子系統(tǒng)S的空間的維數(shù)d確定大于1的第1重正化因子α、0以上d 以下的第2重正化因子Y及0以上的第3重正化因子δ，并通過轉(zhuǎn)換式N' = N/ α d求出重正化的粒子數(shù)N'，通過轉(zhuǎn)換式m' =11^8/(^求出重正化的粒子質(zhì)量111'，通過轉(zhuǎn)換式 ε ‘ = ε α γ求出重正化的相互作用系數(shù)ε ‘ ’及(b)工序，對包含重正化的粒子數(shù)為N'個的粒子、各粒子具有重正化的粒子質(zhì)量m'， 且由所述無因次函數(shù)f和重正化的相互作用系數(shù)ε ‘之積ε ‘ f表示粒子間的相互作用勢能的粒子系統(tǒng)S'執(zhí)行分子動力學(xué)計算。
2.如權(quán)利要求1所述的模擬方法，其特征在于，還具有(c)工序，所述(c)工序如下 在將所述工序(b)中用分子動力學(xué)計算求出的所述粒子系統(tǒng)S'中所含的某一粒子的位置矢量表示為q'，將所述工序(b)中用分子動力學(xué)計算求出的所述粒子系統(tǒng)S'中所含的某一粒子的運動量矢量表示為P'，將所述工序(b)中用分子動力學(xué)計算求出的所述粒子系統(tǒng)S'中所含的某一粒子的速度矢量表示為ν'，將所述工序(b)中的分子動力學(xué)計算中的所述粒子系統(tǒng)S'中所含的某一粒子的某一時間間隔表示為t'時，使用在所述工序(a)中確定的所述第1重正化因子α、第2重正化因子Y及第3重正化因子δ中的至少一個， 執(zhí)行通過轉(zhuǎn)換式q = q' α求出所述粒子系統(tǒng)S中所含的某一粒子的位置矢量q的計算、 通過轉(zhuǎn)換式P = P' /α δ/2求出所述粒子系統(tǒng)S中所含的某一粒子的運動量矢量ρ的計算、通過轉(zhuǎn)換式V = V α (- + δ/2)求出所述粒子系統(tǒng)S中所含的某一粒子的速度矢量ν的計算、及通過轉(zhuǎn)換式t = t' α (1+-δ/2)求出假設(shè)對所述粒子系統(tǒng)S執(zhí)行的分子動力學(xué)計算中的某一時間間隔t的計算中的至少1個計算。
3.如權(quán)利要求1或2所述的模擬方法，其特征在于，所述第2重正化因子γ為0。
4.如權(quán)利要求1 3中任一項所述的模擬方法，其特征在于，在將η設(shè)為正整數(shù)時，所述第1重正化因子為2η。
5.如權(quán)利要求1 4中任一項所述的模擬方法，其特征在于，所述維數(shù)d為3。
6.一種模擬方法，其特征在于，具有(a)工序，在考慮包含N個粒子、各粒子質(zhì)量為m、且能夠以表示對粒子間距離的依賴性的無因次函數(shù)f與相互作用系數(shù)ε之積來表示粒子間的相互作用勢能的粒子系統(tǒng)S 時，使用配置有該粒子系統(tǒng)S的空間的維數(shù)d確定大于1的第1重正化因子α和0以上的第3重正化因子δ，并通過轉(zhuǎn)換式N' =N/Cid求出重正化的粒子數(shù)N'，通過轉(zhuǎn)換式m' = Hia s求出重正化的粒子質(zhì)量m'，及(b)工序，對包含重正化的粒子數(shù)為N'個的粒子、各粒子具有重正化的粒子質(zhì)量m'、 且由所述無因次函數(shù)f與所述相互作用系數(shù)ε之積表示粒子間的相互作用勢能的粒子系統(tǒng)S'執(zhí)行分子動力學(xué)計算。
7.如權(quán)利要求6所述的模擬方法，其特征在于，還具有(c)工序，所述(c)工序如下 在將所述工序(b)中用分子動力學(xué)計算求出的所述粒子系統(tǒng)S'中所含的某一粒子的位置矢量表示為q'，將所述工序(b)中用分子動力學(xué)計算求出的所述粒子系統(tǒng)S'中所含的某一粒子的運動量矢量表示為P'，將所述工序(b)中用分子動力學(xué)計算求出的所述粒子系統(tǒng)S'中所含的某一粒子的速度矢量表示為ν'時，使用在所述工序(a)中確定的所述第1 重正化因子α和所述第3重正化因子δ中的至少一方，執(zhí)行通過轉(zhuǎn)換式q = q' α求出所述粒子系統(tǒng)S中所含的某一粒子的位置矢量q的計算、通過轉(zhuǎn)換式P = ρ' /α M求出所述粒子系統(tǒng)S中所含的某一粒子的運動量矢量ρ的計算、及通過轉(zhuǎn)換式ν = ν' α 8/2求出所述粒子系統(tǒng)S中所含的某一粒子的速度矢量ν的計算中的至少1個計算。
8.一種模擬方法，其特征在于，具有(d)工序，在考慮包含N個粒子、各粒子質(zhì)量為m、且能夠以表示對粒子間距離的依賴性的無因次函數(shù)f與相互作用系數(shù)ε之積來表示粒子間的相互作用勢能的粒子系統(tǒng)S、 并且在配置有所述粒子系統(tǒng)S的空間考慮XYZ正交坐標(biāo)系時，分別使用大于1的X方向的重正化因子α x、Y方向的重正化因子α γ及Z方向的重正化因子α ζ，通過轉(zhuǎn)換式N' =N/ (αχαγαζ)求出重正化的粒子數(shù)N'，并且確定0以上的重正化因子δ，通過轉(zhuǎn)換式％'= ma /求出有關(guān)X方向的重正化的粒子質(zhì)量％'，通過轉(zhuǎn)換式m/ =ma/求出有關(guān)Y方向的重正化的粒子質(zhì)量％'，通過轉(zhuǎn)換式mz' = ma zs求出有關(guān)Z方向的重正化的粒子質(zhì)量 mz'；及(e)工序，對包含重正化的粒子數(shù)為N'個的粒子、各粒子具有重正化的粒子質(zhì)量m'、 且由所述無因次函數(shù)f與所述相互作用系數(shù)ε之積表示粒子間的相互作用勢能的粒子系統(tǒng)S'，在X方向的運動方程式中使用有關(guān)X方向的重正化的粒子質(zhì)量mx'來執(zhí)行分子動力學(xué)計算，在Y方向的運動方程式中使用有關(guān)Y方向的重正化的粒子質(zhì)量％'來執(zhí)行分子動力學(xué)計算，在Z方向的運動方程式中使用有關(guān)Z方向的重正化的粒子質(zhì)量mz'來執(zhí)行分子動力學(xué)計算。
9.如權(quán)利要求8所述的模擬方法，其特征在于，還具有(f)工序，所述(f)工序如下 將在所述工序(e)中用分子動力學(xué)計算求出的所述粒子系統(tǒng)S'中所含的某一粒子的位置矢量表示為q'，將該位置矢量q'的X、Y及Z方向各自的成分表示為qx'、q/及Ck'，將在所述工序(e)中用分子動力學(xué)計算求出的所述粒子系統(tǒng)S'中所含的某一粒子的運動量矢量表示為P'，將該運動量矢量P'的X、Y及Z方向各自的成分表示為px'、ργ'及ρζ'， 將在所述工序(e)中用分子動力學(xué)計算求出的所述粒子系統(tǒng)S'中所含的某一粒子的速度矢量表示為ν'，將該速度矢量ν'的X、Y及Z方向各自的成分表示為Vx'、νγ'及νζ'時， 使用所述X方向的重正化因子a x、Y方向的重正化因子a Y、Z方向的重正化因子Ciz及在所述工序(d)中確定的重正化因子δ，執(zhí)行分別通過轉(zhuǎn)換式qx=qx' ax、qY = qY' ^及 qz = qz' az求出所述粒子系統(tǒng)S中所含的某一粒子的位置矢量q的X、Y及Z方向各自的成分％、％及％的計算、分別通過轉(zhuǎn)換式ρχ = ρχ' /αχδ/2,ργ = ργ' / α//2及ρζ = ρζ' / αζδ/2求出所述粒子系統(tǒng)S中所含的某一粒子的運動量矢量P的Χ、Υ及Z方向各自的成分 px、Py及Pz的計算及分別通過轉(zhuǎn)換式Vx = Vx' αχδ/2,νγ = vY' Cm αζδ/2 求出所述粒子系統(tǒng)S中所含的某一粒子的速度矢量ν的X、Y及Z方向各自的成分\、νγ及 νζ的計算中的至少1個計算。
10.如權(quán)利要求8或9所述的模擬方法，其特征在于，所述X方向的重正化因子a χ、γ方向的重正化因子α Z方向的重正化因子CIz中，至少2個互不相同。
11. 一種模擬程序，其特征在于，用于在計算機中執(zhí)行權(quán)利要求1 10中任一項所述的模擬方法。
全文摘要
提供采用了分子動力學(xué)的新的模擬方法。(a)在考慮包含N個粒子、各粒子質(zhì)量為m、且能夠以表示對粒子間距離的依賴性的無因次函數(shù)f與相互作用系數(shù)ε之積εf來表示粒子間的相互作用勢能的粒子系統(tǒng)S時，使用配置有粒子系統(tǒng)S的空間的維數(shù)d確定大于1的第1重正化因子α、0以上d以下的第2重正化因子γ及0以上的第3重正化因子δ，并通過轉(zhuǎn)換式N′＝N/αd求出重正化的粒子數(shù)N′，通過轉(zhuǎn)換式m′＝mαδ/αγ求出重正化的粒子質(zhì)量m′，通過轉(zhuǎn)換式ε′＝εαγ求出重正化的相互作用系數(shù)ε′。(b)對包含重正化的粒子數(shù)為N′個的粒子、各粒子具有重正化的粒子質(zhì)量m′，且由所述無因次函數(shù)f和重正化的相互作用系數(shù)ε′之積ε′f表示粒子間的相互作用勢能的粒子系統(tǒng)S′執(zhí)行分子動力學(xué)計算。
文檔編號G06F19/00GK102257500SQ200980150638
公開日2011年11月23日 申請日期2009年10月29日 優(yōu)先權(quán)日2008年12月19日
發(fā)明者市島大路 申請人:住友重機械工業(yè)株式會社