专利名称：模拟方法及程序的制作方法
技术领域：
本发明涉及模拟方法及程序，尤其涉及利用了分子动力学的模拟方法及用于使其在计算机中执行的程序。
背景技术：
利用了分子动力学的计算机模拟正在不断施行。分子动力学中，构成作为模拟对象的系统的粒子的运动方程式能够以数值方式解出。若作为模拟对象的系统的粒子数增加，则所需的计算量也增加。以目前计算机的运算能力，例如只能进行到粒子数为10万左右的系统的模拟。为了减少模拟所需的计算量，例如在非专利文献1中，尝试了在分子动力学中应用重正化方法。另外，本申请所涉及的发明由本申请发明人完成，并公开于专利文献1中。现有技术文献专利文献专利文献1 日本特开2006-285866号公报非专利文献1 稻村，武泽，社本，《关于基于重正化方法的可变标度分子动力学模拟》(「繰D込 O方法(二基d〈可变7 分子動力学* ^工>一* 307」)， 日本机械学会论文集(A篇)，1997年，第63卷，608号，P. 202-207发明的概要发明要解决的课题如非专利文献1中专栏“3.数值计算法”所记载，在非专利文献1的模拟方法中无法明确评价温度。
发明内容
本发明的目的之一在于提供一种利用分子动力学的新型模拟方法及用于在计算机中执行该方法的程序。本发明的另一目的在于提供一种可根据重正化群方法明确评价温度的模拟方法及用于在计算机中执行该方法的程序。用于解决课题的手段根据本发明的一种观点，提供一种模拟方法，其方法具有如下工序(a)工序，考虑包含N个粒子、各粒子质量为m、且能够以表示对粒子间距离的依赖性的无因次函数f与相互作用系数ε之积ε f来表示粒子间的相互作用势能的粒子系统S时，使用配置有该粒子系统S的空间的维数d确定大于1的第1重正化因子α、0以上d以下的第2重正化因子 Y及0以上的第3重正化因子δ，并通过转换式N' =N/Cid求出重正化的粒子数N'，通过转换式m' =11^8/(^求出重正化的粒子质量111'，通过转换式ε' = ε α γ求出重正化的相互作用系数ε ‘；及(b)工序，对包含重正化的粒子数为N'个的粒子、各粒子具有重正化的粒子质量m'、且由所述无因次函数f与重正化的相互作用系数ε ‘之积ε ‘ f 表示粒子间的相互作用势能的粒子系统S'执行分子动力学计算。另外，根据本发明的其他观点，提供一种模拟方法，其方法具有如下工序(a)工序，考虑包含N个粒子、各粒子质量为m、且能够以表示对粒子间距离的依赖性的无因次函数f及相互作用系数ε之积来表示粒子间的相互作用势能的粒子系统S时，使用配置有该粒子系统S的空间的维数d确定大于1的第1重正化因子α和0以上的第3重正化因子S，并通过转换式N' =N/Cid求出重正化的粒子数N'，通过转换式m' =mas求出重正化的粒子质量m'；及(b)工序，对包含重正化的粒子数N'个的粒子、各粒子具有重正化的粒子质量m'、且由所述无因次函数f与所述相互作用系数ε之积表示粒子间的相互作用势能的粒子系统S'执行分子动力学计算。并且，根据本发明的其他观点，提供一种模拟方法，其方法具有如下工序(d)工序，考虑包含N个粒子、各粒子质量为m、且能够以表示对粒子间距离的依赖性的无因次函数f与相互作用系数ε之积ε f来表示粒子间的相互作用势能的粒子系统S，并在配置有所述粒子系统S的空间上考虑XYZ正交坐标系时，分别使用大于1的X方向的重正化因子 αχ、Y方向的重正化因子方向的重正化因子αζ，通过转换式N' = N/( α χ α γ α ζ) 求出重正化的粒子数N'，并且确定0以上的重正化因子δ，通过转换式％' =maxs求出 X方向上的重正化的粒子质量％'，通过转换式m/ = ma /求出Y方向上的重正化的粒子质量％'，通过转换式mz' = ma /求出Z方向上的重正化的粒子质量mz'；及(e)工序， 对包含重正化的粒子数为N'个的粒子、各粒子具有重正化的粒子质量m'、且由所述无因次函数f与所述相互作用系数ε之积表示粒子间的相互作用势能的粒子系统S'，在 X方向的运动方程式中使用X方向上的重正化的粒子质量％'执行分子动力学计算，在Y方向的运动方程式中使用Y方向上的重正化的粒子质量％'执行分子动力学计算，在Z方向的运动方程式中使用Z方向上的重正化的粒子质量mz'执行分子动力学计算。发明的效果重正化的粒子数N'少于粒子数N。因此，能够以比对粒子系统S执行分子动力学计算时更少的计算量执行对粒子系统S'的分子动力学计算。另外，由于任意地确定0以上的重正化因子δ而进行计算，因此可灵活地执行模拟。
图IA是用于说明标度转换法的导出方法的图。图IB是用于说明标度转换法的导出方法的图。图2是表示使用重正化群分子动力学的对铝的模拟结果的图表。图3Α是用于说明关于S-W势的图。图3Β是表示使用重正化群分子动力学的对硅的模拟结果的图表。图4是表示利用重正化群分子动力学求出分散关系的模拟结果的图表。符号的说明a -第1重正化因子，N-粒子数，N'-重正化的粒子数。
具体实施例方式
首先，对分子动力学(Molecular Dynamics ；MD)进行简洁说明。对由N个粒子(例如原子)构成、哈密顿算符H表示为如下(1)的粒子系统进行考虑。其中，各粒子质量为m，
粒子j的运动量矢量为Pj,粒子j的位置矢量(位置坐标)为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的速度矢量。分子动力学中通过数值积分并解出由公式(2)及(3)表示的运动方程式来对构成粒子系统的各粒子求出各时刻的各粒子的运动量矢量(或者速度矢量)及位置矢量。另外，在数值积分中多使用Verlet法。关于 Verlet 法，例如在“J. M. Thi jssen，; Computational Physics，，Cambridge University Press (1999)，，的p. 175中所说明。根据由分子动力学计算得到的各粒子的位置、速度能够算出系统的各种物理量。下面，对应用重正化群方法的分子动力学(将此称为重正化群分子动力学)进行概括性说明。在重正化群分子动力学中，使希望求出物理量的期望的系统S对应于由少于系统 S的粒子数构成的系统(将此称为重正化的系统S')。其次，对重正化的系统S'执行分子动力学计算。并且，使对重正化的系统S'计算的结果对应于期望的系统S。由此，能够以比直接对期望的系统S执行分子动力学计算时更少的计算量对期望的系统S算出物理量。 将用于使期望的系统S中的物理量(例如，粒子数、各粒子质量等)与重正化的系统S'中的物理量相互对应的转换法称为标度转换法。下面，对本申请发明人发现的标度转换法进行说明。假定希望算出物理量的期望的粒子系统S由N个粒子构成。设定粒子系统S的所有哈密顿算符H表示为如公式(1)所
7J\ ο
_ 2
孖=Σ
；=1
PN
τ^+ Σ 敝-q；) 1tjti某一粒子j的哈密顿算符&表示为如下。
D 2 N
Η^^Ζ+ΣΦ^-^'·· (4)
^m i = j+\
为了讨论重正化，通过表示对粒子间距的依赖性并被无因次化的函数f、与表示相互作用的强度并具有能量因次的系数ε (将此称为相互作用系数ε)之积，表示相互作用势能Φ如下。
φ( = ε
/\ q. -q.Φ\^—(\0=^」——L (5)
V σ J例如，对惰性原子适用如下相互作用位势。其中，r为原子间距离。相互作用系数 ε及常数ο取对应原子种类的值。下面，对从哈密顿算符H导出重正化的哈密顿算符H'的方法进行说明。如以下详细所示，重正化的哈密顿算符H'通过如下方法获得执行对粒子系统S的分配函数Ζ(β ) 的积分的一部分，粗粒化哈密顿算符，接着对积分变量进行再标度。对粒子数固定的正则系综(Canonical ensemble)的分配函数Z ( β )相对于粒子系统S表示如下。Ζ(β) = / drNexp(-^H(p, q))··· (7)其中，系数β通过系统的温度T与玻尔兹曼常数、定义为β = l/(kBT)。d「N为相位空间内的体积要素，更详细地表示为如下。
1 N1..drN = —= —DNpDNq …(8)
VV N 3VyN其中，h3N(h为普朗克常数)。另外，将DpN及D；1分别定义为如下。
NDno = Yidni
P J=I 1
NDf = Yidai
q j=i 1对谐振子以外的位势的重正化进行讨论较为困难。因此，将具有鞍点的相互作用位势分成3个区域来考虑。由于粒子间距r —c 与r — 0为固定点，因此重正化转换时位势不变。因此，假设鞍点附近的重正化转换对所有粒子间距r都成立。对鞍点的周围进行泰勒展开，保留至平方波动。若将鞍点的位置设为Α，相对变位设为SU，则Su的一次项消失，得到如下公式。φ{τ0+δη) = φ(τ0) + ]-Su2- (9)
2 or
L」r=r0继而，关于粒子j表述鞍点附近的哈密顿算符如下。
P,2 V^H1 =^+ V
]2m其中，Ui及屮分别表示从粒子i及j的鞍点的变位。φ"为关于粒子间距r的二阶微分。下面，对分配函数Ζ( β)的粒子间相互作用部分的粗粒化进行说明。首先，对配置于1维链上的粒子系统的粗粒化进行观察后，向配置于简单立方晶格上的粒子系统扩张。如图IA所示，1维链上的格点(lattice point)上相互最邻接配置粒子i与j，且相互最邻接配置粒子j与k。粒子i与k的中间配置有粒子j。用实线表示最邻接的粒子彼此的相互作用(最邻接连接)。考虑消去并粗粒化粒子j。若写出粒子j所参与的相互
0
1
Xc
\ /
(ro
+
_
y
U
I
UI
\J
(^o
W"
於
6
σ
r
8作用，并对位于粒子i与k之间的粒子j的变位~执行积分，则关于粒子j的影响被重正化后的、粒子i与k的相互作用位势φ'(化-qk)得出如下公式。
权利要求
1.一种模拟方法，其特征在于，具有(a)工序，在考虑包含N个粒子、各粒子质量为m、且能够以表示对粒子间距离的依赖性的无因次函数f与相互作用系数ε之积来表示粒子间的相互作用势能的粒子系统 S时，使用配置有该粒子系统S的空间的维数d确定大于1的第1重正化因子α、0以上d 以下的第2重正化因子Y及0以上的第3重正化因子δ，并通过转换式N' = N/ α d求出重正化的粒子数N'，通过转换式m' =11^8/(^求出重正化的粒子质量111'，通过转换式 ε ‘ = ε α γ求出重正化的相互作用系数ε ‘ ’及(b)工序，对包含重正化的粒子数为N'个的粒子、各粒子具有重正化的粒子质量m'， 且由所述无因次函数f和重正化的相互作用系数ε ‘之积ε ‘ f表示粒子间的相互作用势能的粒子系统S'执行分子动力学计算。
2.如权利要求1所述的模拟方法，其特征在于，还具有(c)工序，所述(c)工序如下 在将所述工序(b)中用分子动力学计算求出的所述粒子系统S'中所含的某一粒子的位置矢量表示为q'，将所述工序(b)中用分子动力学计算求出的所述粒子系统S'中所含的某一粒子的运动量矢量表示为P'，将所述工序(b)中用分子动力学计算求出的所述粒子系统S'中所含的某一粒子的速度矢量表示为ν'，将所述工序(b)中的分子动力学计算中的所述粒子系统S'中所含的某一粒子的某一时间间隔表示为t'时，使用在所述工序(a)中确定的所述第1重正化因子α、第2重正化因子Y及第3重正化因子δ中的至少一个， 执行通过转换式q = q' α求出所述粒子系统S中所含的某一粒子的位置矢量q的计算、 通过转换式P = P' /α δ/2求出所述粒子系统S中所含的某一粒子的运动量矢量ρ的计算、通过转换式V = V α (- + δ/2)求出所述粒子系统S中所含的某一粒子的速度矢量ν的计算、及通过转换式t = t' α (1+-δ/2)求出假设对所述粒子系统S执行的分子动力学计算中的某一时间间隔t的计算中的至少1个计算。
3.如权利要求1或2所述的模拟方法，其特征在于，所述第2重正化因子γ为0。
4.如权利要求1 3中任一项所述的模拟方法，其特征在于，在将η设为正整数时，所述第1重正化因子为2η。
5.如权利要求1 4中任一项所述的模拟方法，其特征在于，所述维数d为3。
6.一种模拟方法，其特征在于，具有(a)工序，在考虑包含N个粒子、各粒子质量为m、且能够以表示对粒子间距离的依赖性的无因次函数f与相互作用系数ε之积来表示粒子间的相互作用势能的粒子系统S 时，使用配置有该粒子系统S的空间的维数d确定大于1的第1重正化因子α和0以上的第3重正化因子δ，并通过转换式N' =N/Cid求出重正化的粒子数N'，通过转换式m' = Hia s求出重正化的粒子质量m'，及(b)工序，对包含重正化的粒子数为N'个的粒子、各粒子具有重正化的粒子质量m'、 且由所述无因次函数f与所述相互作用系数ε之积表示粒子间的相互作用势能的粒子系统S'执行分子动力学计算。
7.如权利要求6所述的模拟方法，其特征在于，还具有(c)工序，所述(c)工序如下 在将所述工序(b)中用分子动力学计算求出的所述粒子系统S'中所含的某一粒子的位置矢量表示为q'，将所述工序(b)中用分子动力学计算求出的所述粒子系统S'中所含的某一粒子的运动量矢量表示为P'，将所述工序(b)中用分子动力学计算求出的所述粒子系统S'中所含的某一粒子的速度矢量表示为ν'时，使用在所述工序(a)中确定的所述第1 重正化因子α和所述第3重正化因子δ中的至少一方，执行通过转换式q = q' α求出所述粒子系统S中所含的某一粒子的位置矢量q的计算、通过转换式P = ρ' /α M求出所述粒子系统S中所含的某一粒子的运动量矢量ρ的计算、及通过转换式ν = ν' α 8/2求出所述粒子系统S中所含的某一粒子的速度矢量ν的计算中的至少1个计算。
8.一种模拟方法，其特征在于，具有(d)工序，在考虑包含N个粒子、各粒子质量为m、且能够以表示对粒子间距离的依赖性的无因次函数f与相互作用系数ε之积来表示粒子间的相互作用势能的粒子系统S、 并且在配置有所述粒子系统S的空间考虑XYZ正交坐标系时，分别使用大于1的X方向的重正化因子α x、Y方向的重正化因子α γ及Z方向的重正化因子α ζ，通过转换式N' =N/ (αχαγαζ)求出重正化的粒子数N'，并且确定0以上的重正化因子δ，通过转换式％'= ma /求出有关X方向的重正化的粒子质量％'，通过转换式m/ =ma/求出有关Y方向的重正化的粒子质量％'，通过转换式mz' = ma zs求出有关Z方向的重正化的粒子质量 mz'；及(e)工序，对包含重正化的粒子数为N'个的粒子、各粒子具有重正化的粒子质量m'、 且由所述无因次函数f与所述相互作用系数ε之积表示粒子间的相互作用势能的粒子系统S'，在X方向的运动方程式中使用有关X方向的重正化的粒子质量mx'来执行分子动力学计算，在Y方向的运动方程式中使用有关Y方向的重正化的粒子质量％'来执行分子动力学计算，在Z方向的运动方程式中使用有关Z方向的重正化的粒子质量mz'来执行分子动力学计算。
9.如权利要求8所述的模拟方法，其特征在于，还具有(f)工序，所述(f)工序如下 将在所述工序(e)中用分子动力学计算求出的所述粒子系统S'中所含的某一粒子的位置矢量表示为q'，将该位置矢量q'的X、Y及Z方向各自的成分表示为qx'、q/及Ck'，将在所述工序(e)中用分子动力学计算求出的所述粒子系统S'中所含的某一粒子的运动量矢量表示为P'，将该运动量矢量P'的X、Y及Z方向各自的成分表示为px'、ργ'及ρζ'， 将在所述工序(e)中用分子动力学计算求出的所述粒子系统S'中所含的某一粒子的速度矢量表示为ν'，将该速度矢量ν'的X、Y及Z方向各自的成分表示为Vx'、νγ'及νζ'时， 使用所述X方向的重正化因子a x、Y方向的重正化因子a Y、Z方向的重正化因子Ciz及在所述工序(d)中确定的重正化因子δ，执行分别通过转换式qx=qx' ax、qY = qY' ^及 qz = qz' az求出所述粒子系统S中所含的某一粒子的位置矢量q的X、Y及Z方向各自的成分％、％及％的计算、分别通过转换式ρχ = ρχ' /αχδ/2,ργ = ργ' / α//2及ρζ = ρζ' / αζδ/2求出所述粒子系统S中所含的某一粒子的运动量矢量P的Χ、Υ及Z方向各自的成分 px、Py及Pz的计算及分别通过转换式Vx = Vx' αχδ/2,νγ = vY' Cm αζδ/2 求出所述粒子系统S中所含的某一粒子的速度矢量ν的X、Y及Z方向各自的成分\、νγ及 νζ的计算中的至少1个计算。
10.如权利要求8或9所述的模拟方法，其特征在于，所述X方向的重正化因子a χ、γ方向的重正化因子α Z方向的重正化因子CIz中，至少2个互不相同。
11. 一种模拟程序，其特征在于，用于在计算机中执行权利要求1 10中任一项所述的模拟方法。
全文摘要
提供采用了分子动力学的新的模拟方法。(a)在考虑包含N个粒子、各粒子质量为m、且能够以表示对粒子间距离的依赖性的无因次函数f与相互作用系数ε之积εf来表示粒子间的相互作用势能的粒子系统S时，使用配置有粒子系统S的空间的维数d确定大于1的第1重正化因子α、0以上d以下的第2重正化因子γ及0以上的第3重正化因子δ，并通过转换式N′＝N/αd求出重正化的粒子数N′，通过转换式m′＝mαδ/αγ求出重正化的粒子质量m′，通过转换式ε′＝εαγ求出重正化的相互作用系数ε′。(b)对包含重正化的粒子数为N′个的粒子、各粒子具有重正化的粒子质量m′，且由所述无因次函数f和重正化的相互作用系数ε′之积ε′f表示粒子间的相互作用势能的粒子系统S′执行分子动力学计算。
文档编号G06F19/00GK102257500SQ200980150638
公开日2011年11月23日 申请日期2009年10月29日 优先权日2008年12月19日
发明者市岛大路 申请人:住友重机械工业株式会社