分子动力学资料-moldy和GULP结合

实验4:分子动力学资料

在本次实验中，我们将采用MOLDY以及GULP作为我们分子动力学(MD)的软件．MOLDY使用手册见下面连接．

/~keithr/mold y-manual/moldy.html

这是非常好的一本手册，它是一个全面学习分子动力学的非常好的资源．

因为GULP首先是一个静态的能量代码，所以它的MD特性并是那么好，在GULP使用手册以及实验１的GULP资料中都有介绍．

另外的一些推荐的连接：

http://www.fisica.uniud.it/~ercolessi/md/md/

http://www.fisica.uniud.it/~ercolessi/md/f90/

第一个连接是介绍MD,是Furio Ercolessi写的．它非常棒并且在阐述分子动力学及原子间作用势上简单易懂．强烈推荐你读这个链接．

第二个连接是一些关于分子动力学代码例子，用Fortran90写的。你看到这些代码样本时，会觉得代码本身不是很复杂（不包括复杂的能量方法）．

还有很多其它有用的分子动力学开源代码，下面网页中例举了很多．

/resources/MDPackages.html

包括CHARMM,，MOSCITO,，MARVIN，DYNAMO等等许多．

第1页

MD简要摘要

分子动力学遵守经典力学,特别地，对所有原子求解牛顿运动方程:

这里i表示原子.所以原则上,对于任何初始参数(位置,速度),将来任意时刻的位

置,速度,加速度都是确定的。原子所受的力,可以由势能对原子位置求导获得:

.

这里U是势能.当U的形式是解析的(如势能的形式),计算作用在原子上力是相对

容易的。因此,我们能求出加速度,因为原子质量是已知的.理论上我们能积分加速度得到速度,再积分速度得到位置.

在实践中,所有积分由数值计算给出。控制积分细度(fineness)的参数称为时间步长,δt。知道t时刻的位置,速度和加速度,我们就能得到时刻t+δt的位置,速度以及加速度。

在课堂上,你学了Verlet算法,Verlet是非常通用且简单的得到时间积分的算法。Mold y采用Beeman算法,属于一种预测-矫正算法。大致介绍如下:

x表示坐标,和分别表示"预测"和"矫正"的速度。GULP求解积分采用

leapfrog Verlet,velocity Verlet或者是Gear Predictor-Corrector算法。在起始参数下,选定时间步长,我们能计算t+δt时刻的参数。注意起始速度是随机选取的,但遵从Maxwell-Boltzmann分布。不幸的是,我们开始的速度可能会很不正确,所以在采样(sample)热力学量（例如势能和动能）前我们需要系统先运行几个时间步长。

第2页

这称为达到平衡时间,图示如下:

温度由换算原子的运动速度控制。怎样做？根据均分理论。

动能=3/2kT

在经典力学中

动能=1/2mv2

因此,温度由换算原子的运动速度来改变。注意,当这个发生时,我们不再严格遵循牛顿运动方程。所以我们在平衡时间段内开启温度换算(temperature scaling),而在模拟采样时关闭温度换算。

为了进行MD模拟,我们要有一个被模拟的系统。在我们事例中,我们建立Kr 的超晶胞。为什么必须是超晶胞?回忆实验1中的周期性边界条件。

我们需要用超晶胞看原子运动中的长波波动。

概括MD模拟中需要的

A.能量的方法(能量的导数用来得到加速度)

B.积分的方法(Mold y用Beeman算法,GULP中自己选)

C.积分中的时间步长(你需要找到它)

D.计算的系统(在超晶胞下,原子位置和原子质量)

第3页

调查不同的积分方法

实验4中所有脚本和文件被复制到一个本地目录：

username$cp/home/nedward/LAB4/scripts/*

我们利用GULP比较Verlet,velocity Verlet and5次Gear预测矫正积分算法，计算EAM Ni。

用于这次实验任务的文件放在/home/nedward/LAB4。

Gulp.in.example

为简便起见，我只讲与实验1描述不一样的输入文件中的行。我也会略掉关于超晶胞原子位置描述的几行。

第1行:不使用opti和单一的关键词,我们在MD中将采用md和conv关键词。关键词md表示这是个分子动力学模拟,conv意味着我们做等体积计算。

第8行:时间步长增量，单位为皮秒。你将改变这行。

第9行:温度，单位是开氏温度。在本次练习中你不需要改变它,但当你选择做附加提高题时，你需要编辑这行。

第4页

第10-11行:达到平衡的时间和平衡后的时间（Equilibration and production times）如果数字是个小数,单位为皮秒。如果不是小数,那么是时间步长的整数倍。你可能想改变达到平衡的时间和平衡后的时间，核实你得到平衡结构的平均值。

第12行:声明计算中所用的积分方法。可能的输入为leapfrog verlet,velocity verlet,以及gear。你将改变这行。

在此练习中你将采用的脚本是run_gulp.j,分两部分说明:首先是脚本的开头允许你使用队列,然后是构造一个Gulp输入文件。这两部分在你的机器中的一个文件(run_gulp.j)里.

开头(the header):

第1行:调用shell(不要改变它).

第3-7行:对于队列体系,你可以相应的改变(3)和(5).

第16行:转到你所在目录.

第5页