用户工具

站点工具


atk:分子动力学基础

差别

这里会显示出您选择的修订版和当前版本之间的差别。

到此差别页面的链接

两侧同时换到之前的修订记录前一修订版
后一修订版
前一修订版
atk:分子动力学基础 [2016/09/03 17:09] – [NVT模拟] dong.dongatk:分子动力学基础 [2019/11/16 22:50] (当前版本) – [介绍] dong.dong
行 2: 行 2:
  
 =====介绍===== =====介绍=====
-Virtual NanoLab 和 ATK 可以使分子动力学模拟变得非常简单:只需要向构型(configuration)添加所需的计算器(calculator)和分子动力学模块(MolecularDynamics block),就可以开始计算了。对于不同类型的模拟,什么样的参数和设置是最合适的,在这里给出了一些指导方针。本文简要介绍 VNL-ATK 中分子动力学(MD)功能,并且一步步的解释了如何正确地设置模拟过程以得到想要的结果。 
  
 所谓分子动力学(MD)模拟,是在预先确定条件下(比如温度,压力,应力,外力等等)模拟原子和分子运动的一种方法。分子动力学模拟可以用来研究纳米尺度下的动力学过程,还可以用来计算相图、扩散系数和各种响应函数等诸多性质。 所谓分子动力学(MD)模拟,是在预先确定条件下(比如温度,压力,应力,外力等等)模拟原子和分子运动的一种方法。分子动力学模拟可以用来研究纳米尺度下的动力学过程,还可以用来计算相图、扩散系数和各种响应函数等诸多性质。
-=====方法=====+ 
 +QuantumATK(和 NanoLab)可以使分子动力学模拟变得非常简单:只需要向构型(configuration)添加所需的计算器(calculator)和分子动力学模块(MolecularDynamics block),就可以开始计算了。对于不同类型的模拟,什么样的参数和设置是最合适的,在这里给出了一些指导方针。本文简要介绍 QuantumATK 中分子动力学(MD)功能,并且一步步的解释了如何正确地设置模拟过程以得到想要的结果。 
 + 
 +<WRAP center info 100%> 
 +=== 提示 === 
 +**本教程使用特定版本的QuantumATK创建,因此涉及的截图和脚本参数可能与您实际使用的版本略有区别,请在学习时务必注意。** 
 +</WRAP> 
 +===== 方法原理 =====
  
 分子动力学模拟本质上是对给定初始构型中原子运动的牛顿方程求数值解。这通常是通过对有限时间步长进行积分计算得到的。原子之间的相互作用(也就是原子间力)可以通过不同的方法(可以是密度泛函理论(DFT)或经典力场)来进行计算。这些力决定了原子的加速度并且使得原子的位置和速度传递到下一个时间步长。多次重复这个过程会产生一系列的构型快照,这些快照描述了系统在相空间中的运动轨迹,从而可以进一步从中分析提取出想要的性质。 分子动力学模拟本质上是对给定初始构型中原子运动的牛顿方程求数值解。这通常是通过对有限时间步长进行积分计算得到的。原子之间的相互作用(也就是原子间力)可以通过不同的方法(可以是密度泛函理论(DFT)或经典力场)来进行计算。这些力决定了原子的加速度并且使得原子的位置和速度传递到下一个时间步长。多次重复这个过程会产生一系列的构型快照,这些快照描述了系统在相空间中的运动轨迹,从而可以进一步从中分析提取出想要的性质。
行 11: 行 17:
 在设置模拟之前你需要选择你所感兴趣的计算类型。要考虑得因素有:总能量是否守恒(孤立体系中的情况)?模拟体系与热浴(heat bath)之间的耦合时是否要保持温度不变?系统是否受到外加压力或者应力?基于这些考虑我们应该选择一套合适的模拟参数:时间步长大小, 积分步长的数量(模拟持续的时间),积分算法,初始温度,约束条件等等。 在设置模拟之前你需要选择你所感兴趣的计算类型。要考虑得因素有:总能量是否守恒(孤立体系中的情况)?模拟体系与热浴(heat bath)之间的耦合时是否要保持温度不变?系统是否受到外加压力或者应力?基于这些考虑我们应该选择一套合适的模拟参数:时间步长大小, 积分步长的数量(模拟持续的时间),积分算法,初始温度,约束条件等等。
  
-类似于真实的实验,进行分子动力学模拟需要一些实证上的知识和经验。在本例中,你将会了解分子动力学模拟的基本组成部分,以及如何使用 ATK 的分子动力学模块来进行分子动力学模拟。+类似于真实的实验,进行分子动力学模拟需要一些实证上的知识和经验。在本例中,你将会了解分子动力学模拟的基本组成部分,以及如何使用 QuantumATK 的分子动力学模块来进行分子动力学模拟。
  
 将这些技巧应用于一个简单、著名的测试体系——块体硅,你将逐步学习到如何选择时间步长大小,如何控制温度和压力,以及如何通过约束条件来固定某些原子。 将这些技巧应用于一个简单、著名的测试体系——块体硅,你将逐步学习到如何选择时间步长大小,如何控制温度和压力,以及如何通过约束条件来固定某些原子。
行 68: 行 74:
 {{ :atk:mdbasical05.png?600 |}} {{ :atk:mdbasical05.png?600 |}}
  
-总体上,系统的能量保持在一个相对于平均值小于10-5量级上的小涨落内。重要的是,平均值是固定的,没有发现漂移。这些发现表明所选分子动力学(MD)的设置较为恰当。+总体上,系统的能量保持在一个相对于平均值小于 $10^{-5}$ 量级上的小涨落内。重要的是,平均值是固定的,没有发现漂移。这些发现表明所选分子动力学(MD)的设置较为恰当。
  
 我们看看如果选择了一个很大的时间步长大小会发生什么。修改初始温度为 5000 K 和时间步长大小为 5 fs,重复上述模拟。你可以再一次使用MolecularDynamics 模块,也可以将 ''md-nve-1fs.py'' 直接拖拽到 **Editor** 按钮。记得改一下轨迹(trajectory)文件名。保存这个新修改的脚本为 ''md-nve-5000K-5fs.py'' 并运行。 我们看看如果选择了一个很大的时间步长大小会发生什么。修改初始温度为 5000 K 和时间步长大小为 5 fs,重复上述模拟。你可以再一次使用MolecularDynamics 模块,也可以将 ''md-nve-1fs.py'' 直接拖拽到 **Editor** 按钮。记得改一下轨迹(trajectory)文件名。保存这个新修改的脚本为 ''md-nve-5000K-5fs.py'' 并运行。
行 76: 行 82:
 {{ :atk:mdbasical06.png?600 |}} {{ :atk:mdbasical06.png?600 |}}
  
-核心要点是,大的时间步长看上去会提高模拟效率,但实际上会增运动方程数值积分方案的误。原子力在一个积分步长内近似守恒的基本假设失效。本质上,所选时间步长应该足够小以解析原子的最高振动频率(也就是说,其应该远小于最小的振动周期)。所以如果你有轻原子(比如氢),那么通常你需要选择比只有重原子(比如金)时更小的时间步长。如果在你的计算中有不同的元素,且温度较高或者原子远离其平衡构型(即粒子受到很大的力)时,也需要选择小的时间步长。如果你不知道使用多大的时间步长,对于多数系统你可以选择 1 fs 作为一个比较安全的开始。通过监测感兴趣的条件下的 NVE 模拟总能是否守恒,你可以对大一些的时间步长值进行评估。+核心要点是,大的时间步长看上去会提高模拟效率,但实际上会增运动方程数值积分方案的误。原子力在一个积分步长内近似守恒的基本假设失效。本质上,所选时间步长应该足够小以解析原子的最高振动频率(也就是说,其应该远小于最小的振动周期)。所以如果你有轻原子(比如氢),那么通常你需要选择比只有重原子(比如金)时更小的时间步长。如果在你的计算中有不同的元素,且温度较高或者原子远离其平衡构型(即粒子受到很大的力)时,也需要选择小的时间步长。如果你不知道使用多大的时间步长,对于多数系统你可以选择 1 fs 作为一个比较安全的开始。通过监测感兴趣的条件下的 NVE 模拟总能是否守恒,你可以对大一些的时间步长值进行评估。
  
 另一个需要知道的分子动力学模拟的重要方面是,一个系统需要一定的时间来与所设置的外界温度和压力达到平衡。在本例中,系统温度和总能总的涨落的量级达到稳态用了大约 400 fs。这点很重要,因为任何可观察量的测量需要在系统达到平衡之后才可以实现(除非你要研究非平衡现象)。在本例中你需要等待至少 400 fs。总的来说,分子动力学模拟的平衡时间最大可以直到几百纳秒(比如在生物分子或者聚合物系统),并且依赖于系统中出现的最长的弛豫时间。所以你应该仔细检查,所感兴趣的可观察量是否在你的模拟中达到稳态。 另一个需要知道的分子动力学模拟的重要方面是,一个系统需要一定的时间来与所设置的外界温度和压力达到平衡。在本例中,系统温度和总能总的涨落的量级达到稳态用了大约 400 fs。这点很重要,因为任何可观察量的测量需要在系统达到平衡之后才可以实现(除非你要研究非平衡现象)。在本例中你需要等待至少 400 fs。总的来说,分子动力学模拟的平衡时间最大可以直到几百纳秒(比如在生物分子或者聚合物系统),并且依赖于系统中出现的最长的弛豫时间。所以你应该仔细检查,所感兴趣的可观察量是否在你的模拟中达到稳态。
行 85: 行 91:
  
 ====设置脚本==== ====设置脚本====
-你将以室温下结构的平衡态为开始。用上述同样的方法和同样的模块(NewCalculator和Optimize ‣ MolecularDynamics) 来准备脚本。实际上,如果你先前的脚本还存储在**ScriptGenerator**里,你可以重新使用它。+你将以室温下结构的平衡态为开始。用上述同样的方法和同样的模块(NewCalculator 和 Optimize ‣ MolecularDynamics) 来准备脚本。实际上,如果你先前的脚本还存储在 **ScriptGenerator** 里,你可以重新使用它。
  
-双击MolecularDynamics模块打开设置。由于我们要运行NVT系综下的模拟,我们这回需要选择一个不同的 Type。对于NVT模拟,ATK提供了四种选项:+双击 MolecularDynamics 模块打开设置。由于我们要运行 NVT 系综下的模拟,我们这回需要选择一个不同的 Type。对于 NVT 模拟,QuantumATK 提供了四种选项:
   - NVT Berendsen((J. C. Berendsen, J. P. M. Postma, W. F. van Gunsteren, A. DiNola, and J. R. Haak: Molecular-Dynamics with Coupling to an External Bath. J. Chem. Phys. 81, 3684, 1982))    - NVT Berendsen((J. C. Berendsen, J. P. M. Postma, W. F. van Gunsteren, A. DiNola, and J. R. Haak: Molecular-Dynamics with Coupling to an External Bath. J. Chem. Phys. 81, 3684, 1982)) 
   - NVT Nose Hoover((W. G. Hoover: Canonical dynamics: Equilibrium phase-space distributions. Phys. Rev. A 31, 1695, 1985))   - NVT Nose Hoover((W. G. Hoover: Canonical dynamics: Equilibrium phase-space distributions. Phys. Rev. A 31, 1695, 1985))
行 93: 行 99:
   - NVT Nose Hoover Chain((G. J. Martyna, M. L. Klein, M. Tuckerman: Nosé–Hoover chains: The canonical ensemble via continuous dynamics . J. Chem. Phys. 97, 2635, 1992))   - NVT Nose Hoover Chain((G. J. Martyna, M. L. Klein, M. Tuckerman: Nosé–Hoover chains: The canonical ensemble via continuous dynamics . J. Chem. Phys. 97, 2635, 1992))
  
-在**MolecularDynamics**界面工具栏,选择NVT Nose Hoover chain类型。设置步数为5000和log interval为20。保存轨迹(trajectory)到文件traj-nvt-nosehooverchain-300K.nc选择Maxwell-Boltzmann初始速度为300K。时间步长默认1fs,Reservoir temperature (虚拟热浴温度)选择为300K,Thermostat timescale 选择100fs。后面的值决定了系统温度多快达到虚拟热浴温度。我们将会在本部分后面仔细研究这个参数的影响。Final temperature值总是和虚拟热浴温度同步,除非你在这里明确的指定一个不同的值。这个参数可被用来在模拟过程中线性增大或线性减小温度(参见参考手册NVTNoseHooverChain)。Thermostat chain length这个参数决定了调用多少个并发的调温器(thermostat)来控制温度。默认值三可以在多数情况下够用。只有当你面临一个强的持久的温度振荡时,可以增加该值以获得更稳定的模拟。+在 **MolecularDynamics** 界面工具栏,选择 NVT Nose Hoover chain 类型。设置步数为 5000 和 log interval 为 20。保存轨迹(trajectory)到文件 ''traj-nvt-nosehooverchain-300K.nc'' 选择 Maxwell-Boltzmann 初始速度为 300 K。时间步长默认 1 fs,Reservoir temperature (虚拟热浴温度)选择为 300 K,Thermostat timescale 选择 100 fs。后面的值决定了系统温度多快达到虚拟热浴温度。我们将会在本部分后面仔细研究这个参数的影响。Final temperature值总是和虚拟热浴温度同步,除非你在这里明确的指定一个不同的值。这个参数可被用来在模拟过程中线性增大或线性减小温度(参见参考手册 NVTNoseHooverChain)。Thermostat chain length 这个参数决定了调用多少个并发的调温器(thermostat)来控制温度。默认值三可以在多数情况下够用。只有当你面临一个强的持久的温度振荡时,可以增加该值以获得更稳定的模拟。
  
 {{ :atk:mdbasical07.png?600 |}} {{ :atk:mdbasical07.png?600 |}}
  
-点击 OK关闭分子动力学设置,发送脚本到**JobManager**并运行。+点击 OK 关闭分子动力学设置,发送脚本到 **JobManager** 并运行。
  
-模拟结束,轨迹文件将会出现在**Lab floor**。使用**Movie Tool**打开。+模拟结束,轨迹文件将会出现在 **Lab floor**。使用 **Movie Tool** 打开。
  
 {{ :atk:mdbasical08.png?600 |}} {{ :atk:mdbasical08.png?600 |}}
  
-观察一下温度曲线。你会发现在一个短暂的初始平衡时期之后,曲线在所选的模拟热浴温度300K附近轻微涨落。实际上,这与之前部分NVE模拟的结果看上去没有太大的不同。与NVE模拟不同的是,系统最终平衡温度不依赖于通过Maxwell-Boltzmann分布设定的初始温度。现在让我们更细致的观察。+观察一下温度曲线。你会发现在一个短暂的初始平衡时期之后,曲线在所选的模拟热浴温度 300 K 附近轻微涨落。实际上,这与之前部分 NVE 模拟的结果看上去没有太大的不同。与 NVE 模拟不同的是,系统最终平衡温度不依赖于通过 Maxwell-Boltzmann 分布设定的初始温度。现在让我们更细致的观察。
  
-打开上一个模拟的**Script Generator**。双击**MolecularDynamics**模块改变虚拟热浴(Reservoir)温度为1000K。Maxwell-Boltzmann温度为默认300K。最后,将文件名改为 traj-nvt-nosehooverchain-1000K.nc。关闭设置,运行脚本。+打开上一个模拟的 **Script Generator**。双击 **MolecularDynamics** 模块改变虚拟热浴(Reservoir)温度为 1000 K。Maxwell-Boltzmann 温度为默认 300 K。最后,将文件名改为 ''traj-nvt-nosehooverchain-1000K.nc''。关闭设置,运行脚本。
  
 结果应该与下图类似。 结果应该与下图类似。
行 111: 行 117:
 {{ :atk:mdbasical09.png?600 |}} {{ :atk:mdbasical09.png?600 |}}
  
-温度曲线起始于初始值300K(通过Maxwell-Boltzmann分布选取)。接着,由于系统与热浴接触,温度迅速提升到虚拟热浴(reservoir)温度1000K+温度曲线起始于初始值 300 K(通过 Maxwell-Boltzmann 分布选取)。接着,由于系统与热浴接触,温度迅速提升到虚拟热浴(reservoir)温度 1000 K
  
-你也会注意到总能的增加。与NVE模拟相比,这并不表示参数选取不合适,而是反映了系统在与虚拟热浴交换能量这一的事实。+你也会注意到总能的增加。与 NVE 模拟相比,这并不表示参数选取不合适,而是反映了系统在与虚拟热浴交换能量这一的事实。
  
-与上一个模拟相比,在模拟过程中原子的动作增大。如果你将虚拟热浴温度提升到更高的值(比如3000K),你将会发现晶体开始融化和部分非晶化。对晶体高温退火实际上是一个产生非晶结构材料的通常方法。你将会在[[atk:生成无定形结构]]这个例子中学到更多关于这个方法。+与上一个模拟相比,在模拟过程中原子的动作增大。如果你将虚拟热浴温度提升到更高的值(比如 3000 K),你将会发现晶体开始融化和部分非晶化。对晶体高温退火实际上是一个产生非晶结构材料的通常方法。你将会在 [[atk:生成无定形结构]] 这个例子中学到更多关于这个方法。
  
-上面说过,调温器(thermostat)时间标度参数可以控制向目标温度接近过程的粒子速度。现在你将使用一个不同的调温器时间标度值来重复同样的模拟。打开上一个模拟的脚本,打开MolecularDynamics设置。保持其他设置不变,改变Thermostat timescale 值为500fs。为了完全获取这个持续很久的平衡相,增加步数到10000,改变文件名为 traj-nvt-nosehooverchain-1000K-500fs.nc。关闭设置,开始模拟。+上面说过,调温器(thermostat)时间标度参数可以控制向目标温度接近过程的粒子速度。现在你将使用一个不同的调温器时间标度值来重复同样的模拟。打开上一个模拟的脚本,打开 MolecularDynamics 设置。保持其他设置不变,改变 Thermostat timescale 值为 500 fs。为了完全获取这个持续很久的平衡相,增加步数到10000,改变文件名为 ''traj-nvt-nosehooverchain-1000K-500fs.nc''。关闭设置,开始模拟。
  
 模拟结果应该与下图类似。 模拟结果应该与下图类似。
行 125: 行 131:
 定性地,结果与上一个模拟非常类似。然而,你会注意到温度达到虚拟热浴温度过程现在变得很慢。 定性地,结果与上一个模拟非常类似。然而,你会注意到温度达到虚拟热浴温度过程现在变得很慢。
  
-总体上讲,一个小的调温器(thermostat)时间标度参数导致一个与虚拟热浴的紧密耦合,和一个虚拟热浴温度的系统温度。然而,如此紧密的调温器耦合通常伴随着对粒子自然的动力学显著的干扰。所以,当你需要对扩散或者振动等动力学性质有一个准确的测量时,你需要使用一个大的热耦合系数或者考虑在NVE系综上进行模拟以避免调温器(thermostat)对系统不符合自然规律的作用。 +总体上讲,一个小的调温器(thermostat)时间标度参数使得与虚拟热浴的耦合更加紧密从而使系统温度更加密的保持在虚拟热浴温度。然而,如此紧密的调温器耦合通常伴随着对粒子自然的动力学显著的干扰。所以,当你需要对扩散或者振动等动力学性质进行准确的测量时,你需要使用比较大的热耦合系数或者考虑在 NVE 系综上进行模拟以避免调温器(thermostat)对系统不符合自然规律的作用。 
-=====NVT Berendsen和Langevin模拟===== +=====NVT Berendsen 和 Langevin 模拟===== 
-除了Nose-Hoover chain 温度器(thermostat)之外,你还可以选择其他算法,比如 **NVT Berendsen**, **Langevin**, 或者**NVT Nose-Hoover**。最后的这个方法是Nose-Hoover-Chain方法更早和更简单的版本,它的呈现更多的是出于历史原因。如果可能,你应该优先选取Nose-Hoover-Chain方法,因为它给出许多更好的结果。 +除了 Nose-Hoover chain 温度器(thermostat)之外,你还可以选择其他算法,比如 **NVT Berendsen**, **Langevin**, 或者**NVT Nose-Hoover**。最后的这个方法是 Nose-Hoover-Chain 方法更早和更简单的版本,它的呈现更多的是出于历史原因。如果可能,你应该优先选取 Nose-Hoover-Chain 方法,因为它给出许多更好的结果。 
-Berendsen 调温器(thermostat)提供一种算法有效的抑制温度振荡,使得在一些情况下系统更稳定。但是,相比于Nose-Hoover-Chain,它不能精确复制正则系综,尽管偏差非常小。由于这个原因,Berendsen方法应该被用在平衡态模拟,主要是在需要强健温度控制的情况中。 +Berendsen 调温器(thermostat)提供一种算法有效的抑制温度振荡,使得在一些情况下系统更稳定。但是,相比于 Nose-Hoover-Chain,它不能精确得到正则系综,尽管偏差非常小。由于这个原因,Berendsen 方法主要用在平衡态模拟,主要是在需要强健温度控制的情况中。 
-你可以为NVT模拟选择的第三个调温器(thermostat)类型是**Langevin**算法。这个算法求解Langevin方程,明确地在牛顿运动方程中包含摩擦以及随机碰撞力以模拟粒子与虚拟热浴之间的相互作用。 +你可以为 NVT 模拟选择的第三个调温器(thermostat)类型是 **Langevin** 算法。这个算法求解 Langevin 方程,明确地在牛顿运动方程中包含摩擦以及随机碰撞力以模拟粒子与虚拟热浴之间的相互作用。 
-不同于其他温度器(thermostat)类型,Langevin温器(thermostat)将每个粒子分别与虚拟热浴耦合。这产生一个很紧密的耦合,好像系统沉浸在一个虚拟的粘性介质中。但是也以一种显著的方式抑制了自然的动力学。出于此原因,Langevin调温器(thermostat)应主要被用来产生结构或者系综取样,而不是计算动力学性质。 +不同于其他温度器(thermostat)类型,Langevin 温器(thermostat)将每个粒子分别与虚拟热浴耦合。这产生一个很紧密的耦合,好像系统沉浸在一个虚拟的粘性介质中。但是也以一种显著的方式抑制了自然的动力学。出于此原因,Langevin 调温器(thermostat)应主要被用来产生结构或者系综取样,而不是计算动力学性质。 
-在Langevin NVT模拟中,耦合程度可以用 Friction参数来设置。它被定义为量纲为时间的倒数,意味着它与Berendsen或者Nose Hoover温度器(thermostat)的耦合参数的工作方式相反。增加friction值会导致与虚拟热浴更紧密的耦合,和对动力学更严格的修正。+在 Langevin NVT 模拟中,耦合程度可以用 Friction 参数来设置。它被定义为量纲为时间的倒数,意味着它与 Berendsen 或者 Nose Hoover 温度器(thermostat)的耦合参数的工作方式相反。增加 friction 值会导致与虚拟热浴更紧密的耦合,和对动力学更严格的修正。
 =====NPT模拟===== =====NPT模拟=====
 在这个部分你将学到如何进行恒定压力下的模拟。 在这个部分你将学到如何进行恒定压力下的模拟。
  
 ====设置计算==== ====设置计算====
-设置一个与NVT模拟包含相同模块的脚本(BulkConfiguration,NewCalculator, Optimize ‣ MolecularDynamics)。+设置一个与 NVT 模拟包含相同模块的脚本(BulkConfiguration, NewCalculator, Optimize ‣ MolecularDynamics)。
  
-打开**MolecularDynamics**设置。通过查看 Type下面的选项,你会找到两个运行NPT模拟的选项:**NPT Berendsen**和**NPT Melchionna**算法。+打开 **MolecularDynamics** 设置。通过查看 Type 下面的选项,你会找到两个运行 NPT 模拟的选项:**NPT Berendsen** 和 **NPT Melchionna** 算法。
  
-**NPT Berendsen**方法的工作原理和NVT积分器相同,所以不能精确产生正确的物理系综,尽管偏差通常也很小。而且,ATK中Berendsen的实现只考虑压力标量$P=1/3(P_{xx}+P_{yy}+P_{zz})$,(其中$P_{ij}=−σ_{ij}$为从应力张量σ得到的压力张量的分量),视系统为各项同性。如果模拟像液体或者立方晶体这样的各项同性系统,采用标量压力是一个合适的选择。对于各向异性系统,这个方法无法说明不同方向的不同模量。所以,NPT Berndsen方法应只能被用于各项同性材料的平衡模拟。+**NPT Berendsen** 方法的工作原理和 NVT 积分器相同,所以不能精确产生正确的物理系综,尽管偏差通常也很小。而且,QuantumATK 中 Berendsen 的实现只考虑压力标量 $P=1/3(P_{xx}+P_{yy}+P_{zz})$,(其中 $P_{ij}=−σ_{ij}$ 为从应力张量 σ 得到的压力张量的分量),视系统为各项同性。如果模拟像液体或者立方晶体这样的各项同性系统,采用标量压力是一个合适的选择。对于各向异性系统,这个方法无法说明不同方向的不同模量。所以,NPT Berendsen 方法应只能被用于各项同性材料的平衡模拟。
  
-为了将压力分量$P_{ij}$与固定的压力分别耦合,你更应该选择使用NPT Melchionna方法((S. Melchionna, G. Ciccotti, B.L. Holian: Hoover NPT dynamics for systems varying in shape and size. Mol. Phys. 78, 533, 1993))。平衡定温定压1bar下的硅超胞的合适的设置如下图所示。+为了将压力分量 $P_{ij}$ 与固定的压力分别耦合,你更应该选择使用 NPT Melchionna 方法((S. Melchionna, G. Ciccotti, B.L. Holian: Hoover NPT dynamics for systems varying in shape and size. Mol. Phys. 78, 533, 1993))。平衡定温定压 1 bar 下的硅超胞的合适的设置如下图所示。
  
 {{ :atk:mdbasical11.png?600 |}} {{ :atk:mdbasical11.png?600 |}}
  
-分子动力学(Molecular Dynamics)编组框包含了一般的设置,这些我们现在应该很熟悉。10000个积分步长足以展示恒压器(barostat)的效果。Log interval可设置为50步以析涨落。注意这将产生相当大的文件,或许你可以考虑在评估之后删除轨迹文件。+分子动力学(Molecular Dynamics)设置框包含了一般的设置,这些我们现在应该很熟悉。10000 个积分步长足以展示恒压器(barostat)的效果。Log interval 可设置为 50 步以析涨落。注意这将产生相当大的文件,或许你可以考虑在评估之后删除轨迹文件。
  
-另一个编组框包括恒温器(thermostat)和恒压器(barostat)的设置。恒温器(thermostat)的设置我们已经在之前的部分讨论过。设置External stress 为1bar。设置Barostat time scale 为100fs。改变Bulk modulus为98GPa(硅的文献值)。后面的参数不是特别重要,它只是告诉恒温器(thermostat)体积如何随着应力的变化而改变和确保压力平衡的实际时间标度与所设时间标度的一致性。+另一个设置框包括恒温器(thermostat)和恒压器(barostat)的设置。恒温器(thermostat)的设置我们已经在之前的部分讨论过。设置External stress 为 1 bar。设置 Barostat time scale 为 100 fs。改变 Bulk modulus 为 98 GPa(硅的文献值)。后面的参数不是特别重要,它只是告诉恒温器(thermostat)体积如何随着应力的变化而改变和确保压力平衡的实际时间标度与所设时间标度的一致性。
  
 最后,应力耦合编组框允许你选择恒压器(barostat)与哪个压力/应力张量分量进行耦合,这意味着系统只与静水压力/应力耦合,而忽视剪切应力分量。这些设置对当前的例子十分有效。 最后,应力耦合编组框允许你选择恒压器(barostat)与哪个压力/应力张量分量进行耦合,这意味着系统只与静水压力/应力耦合,而忽视剪切应力分量。这些设置对当前的例子十分有效。
  
-要知道, ATK NPT Melchionna方法总是将所有压力分量当做互相独立的。并且,不仅是单个外界压力,你也可以为每一个分量设定独立的参考应力值来模拟力学剪切或者蠕变实验。ATK参考手册[[http://www.quantumwise.com/documents/manuals/latest/ReferenceManual/index.html/ref.nptmelchionna.html|Reference Manual]]和例子[[http://quantumwise.com/publications/tutorials/item/576-simulating-a-creep-experiment-of-polycrystalline-copper|Simulating a creep experiment of polycrystalline copper]]为此功能提供了更详细的解释。+要知道, QuantumATK NPT Melchionna 方法总是将所有压力分量当做互相独立的。并且,不仅是单个外界压力,你也可以为每一个分量设定独立的参考应力值来模拟力学剪切或者蠕变实验。QuantumATK参考手册[[http://www.quantumwise.com/documents/manuals/latest/ReferenceManual/index.html/ref.nptmelchionna.html|Reference Manual]]和例子[[http://quantumwise.com/publications/tutorials/item/576-simulating-a-creep-experiment-of-polycrystalline-copper|Simulating a creep experiment of polycrystalline copper]] 为此功能提供了更详细的解释。
 ====编辑脚本==== ====编辑脚本====
-将脚本送到Editor来增加一些额外的功能。为了可视化体积和压力随时间变化的涨落,将这几行代码附在脚本的末尾:+将脚本送到 Editor 来增加一些额外的功能。为了可视化体积和压力随时间变化的涨落,将这几行代码附在脚本的末尾:
 <code python> <code python>
 # Get the number of snapshots in the trajectory # Get the number of snapshots in the trajectory
行 185: 行 191:
  
 ====模拟结果==== ====模拟结果====
-将脚本送到Job Manager来进行模拟。这将会花费几分钟来完成。使用Movie Tool中的trajectory打开轨迹,你就可以监测在模拟过程中晶胞的尺寸和形状的变化。温度在所选的模拟热浴温度300K周围涨落。+将脚本送到 **Job Manager** 来进行模拟。这将会花费几分钟来完成。使用 **Movie Tool** 中的 trajectory 打开轨迹,你就可以监测在模拟过程中晶胞的尺寸和形状的变化。温度在所选的模拟热浴温度300K周围涨落。
  
-你可以打开Fluctuations_NPT_tau100.png来查看体积的涨落。它应该与下图类似。+你可以打开 ''Fluctuations_NPT_tau100.png'' 来查看体积的涨落。它应该与下图类似。
  
 {{ :atk:MDbasical12.png?500 |}} {{ :atk:MDbasical12.png?500 |}}
  
-上面的曲线显示了晶胞体积显著的振荡。这个振荡是Nose-Hoover类的温度和压力控制下的典型结果。只要他们的量级保持在合理的范围内,并且在模拟过程中不会急剧增加,那么就没有必要担心。振荡的周期可以通过设置恒压器(Barostat)时间尺度来改变。这个值选的越小,振荡频率越高。振荡的振幅随着恒压器(barostat)时间尺度的减小而增加。+上面的曲线显示了晶胞体积显著的振荡。这个振荡是 Nose-Hoover 类的温度和压力控制下的典型结果。只要他们的量级保持在合理的范围内,并且在模拟过程中不会急剧增加,那么就没有必要担心。振荡的周期可以通过设置恒压器(Barostat)时间尺度来改变。这个值选的越小,振荡频率越高。振荡的振幅随着恒压器(barostat)时间尺度的减小而增加。
  
-下面的曲线是模拟过程中的压力,显示了类似的行为。它在目标压力1bar附近振荡,但是振幅较大。即时压力值如此大的涨落其实在分子动力学模拟中很正常,因为压力是基于应力来计算的。而应力对于体积变化非常敏感。此外,压力是一个宏观性质,对于如此小的系统考虑相应的平均值比即时值更合理。+下面的曲线是模拟过程中的压力,显示了类似的行为。它在目标压力 1 bar 附近振荡,但是振幅较大。即时压力值如此大的涨落其实在分子动力学模拟中很正常,因为压力是基于应力来计算的。而应力对于体积变化非常敏感。此外,压力是一个宏观性质,对于如此小的系统考虑相应的平均值比即时值更合理。
  
-NPT模拟一个典型的应用是模拟相变,比如对石英在不同温度下进行一系列模拟并测量平均晶胞体积随温度变化的函数,来研究α-石英到β-石英相变过程,如文献((M. H. Müser, and K. Binder: Molecular dynamics study of the α-β transition in quartz: elastic properties, finite size effects, and hysteresis in the local structure. Phys. Chem. Minerals 28, 746, 2001))所报道。+NPT 模拟一个典型的应用是模拟相变,比如对石英在不同温度下进行一系列模拟并测量平均晶胞体积随温度变化的函数,来研究α-石英到β-石英相变过程,如文献((M. H. Müser, and K. Binder: Molecular dynamics study of the α-β transition in quartz: elastic properties, finite size effects, and hysteresis in the local structure. Phys. Chem. Minerals 28, 746, 2001))所报道。
 =====有约束条件的分子动力学模拟===== =====有约束条件的分子动力学模拟=====
-在进行分子动力学模拟过程中你可以固定选中原子的位置。设置这种约束最简单的方式是打开**Molecular Dynamics**工具然后点击 Add constraints 打开Constraints editor 工具。+在进行分子动力学模拟过程中你可以固定选中原子的位置。设置这种约束最简单的方式是打开 **Molecular Dynamics** 工具然后点击 Add constraints 打开 Constraints editor 工具。
  
 {{ :atk:MDbasical13.png?400 }} {{ :atk:MDbasical13.png?400 }}
  
-在工具栏里你会找到一个表格列出了所有被标签标记的原子组。你可以通过在展示窗口中(用鼠标拖一个矩形或者点住 Ctrl键选择多个原子)选择所需原子并点击Add tag from selection来添加新的原子组。可以在Constraint栏中对每一个原子组指定约束条件。第一种将对应原子组里的所有原子固定在它们的初始位置,而第二种将原子组内的原子整体视为一个刚性体。这意味着这些原子根据它们质心所受的力进行协同的移动。在当前的实行中只可能平移运动,而不能进行刚性体旋转。+在工具栏里你会找到一个表格列出了所有被标签标记的原子组。你可以通过在展示窗口中(用鼠标拖一个矩形或者点住 Ctrl 键选择多个原子)选择所需原子并点击 Add tag from selection 来添加新的原子组。可以在 Constraint 栏中对每一个原子组指定约束条件。第一种将对应原子组里的所有原子固定在它们的初始位置,而第二种将原子组内的原子整体视为一个刚性体。这意味着这些原子根据它们质心所受的力进行协同的移动。在当前的实行中只可能平移运动,而不能进行刚性体旋转。
  
-当使用约束时,你必须注意至少让一个原子保持无约束的。像NVT NoseHoover和NPT Melchionna这样的恒温器(thermostat)需要至少两个自由原子,因为它们还固定了系统的质心。恒温器(thermostat)会自动适应由约束所导致的自由度数的减少。+当使用约束时,你必须注意至少让一个原子保持无约束的。像 NVT NoseHoover 和 NPT Melchionna 这样的恒温器(thermostat)需要至少两个自由原子,因为它们还固定了系统的质心。恒温器(thermostat)会自动适应由约束所导致的自由度数的减少。
  
-然而,你需要知道Movie Tool或者MD Analyzer中一些功能不能与约束共同使用。比如,当分子动力学模拟中有约束时,Movie Tool中显示的温度可能不会反应系统自由原子实际的温度,而会显示一个较低温度。你可以通过在脚本里添加动能(kinetic energy)代码来监测温度,如下例所示:+然而,你需要知道 **Movie Tool** 或者 **MD Analyzer** 中一些功能不能与约束共同使用。比如,当分子动力学模拟中有约束时,Movie Tool 中显示的温度可能不会反应系统自由原子实际的温度,而会显示一个较低温度。你可以通过在脚本里添加动能(kinetic energy)代码来监测温度,如下例所示:
  
 <code python> <code python>
行 218: 行 224:
 </code> </code>
 <WRAP center box 100%> <WRAP center box 100%>
-注意+=== 注意 ===
  
-通常不建议在有约束原子的情况下实行定压力/应力(比如NPT)模拟,因为系统约束部分内产生的应力分布没有从全局应力张量中移除。+通常不建议在有约束原子的情况下实行定压力/应力(比如 NPT)模拟,因为系统约束部分内产生的应力分布没有从全局应力张量中移除。
 </WRAP> </WRAP>
  
 =====器件构型===== =====器件构型=====
-在**DeviceConfiguration**中实现NVE/NVT分子动力学模拟本质上是和在**BulkConfiguration**中一样的方式。拓展到中间区域的电极中的原子自动被约束不动以维持与电极的相容性。虽然技术上可以在**DeviceConfiguration**中实现定压力/应力(NPT)模拟,由于约束的存在,不建议这种计算。对于**DeviceConfiguration**晶胞的优化有多种方法供选择(参见例子[[http://quantumwise.com/support/tutorials/item/878-geometry-optimization-of-interfaces|Geometry optimization of interfaces]])。+在 **DeviceConfiguration** 中实现 NVE/NVT 分子动力学模拟本质上是和在 **BulkConfiguration** 中一样的方式。拓展到中间区域的电极中的原子自动被约束不动以维持与电极的相容性。虽然技术上可以在 **DeviceConfiguration** 中实现定压力/应力(NPT)模拟,由于约束的存在,不建议这种计算。对于 **DeviceConfiguration** 晶胞的优化有多种方法供选择(参见例子[[http://quantumwise.com/support/tutorials/item/878-geometry-optimization-of-interfaces|Geometry optimization of interfaces]])。 
 + 
 +===== 参考 =====
  
 +  * 原文:[[http://docs.quantumwise.com/tutorials/md_basics.html]]
   * 本文翻译:王吉章   * 本文翻译:王吉章
atk/分子动力学基础.1472893773.txt.gz · 最后更改: 2016/09/03 17:09 由 dong.dong

© 2014-2022 费米科技(京ICP备14023855号