LAMMPS教程

LAMMPS教程王延颋中国科学院理论物理研究所中科院超算中心培训北京2012年9月17-18日121.分子动力学模拟基本概念3计算科学理论科学实验科学模型框架求解模型参数指导1.1.计算科学─理论与实验科学的桥梁求解解析理论无法求解的问题模拟多体问题，得到更贴近实际体系的结果模拟实验做起来困难或不可能做的条件计算用到大量近似和模拟理想化实验条件，所以还是更贴近理论计算物理的作用41.2.分子模拟─研究原子分子层面的物性计算物理包括数值求解和计算机模拟两大类对物质的模拟大致分为宏观、介观、微观三个层次分子模拟方法分为MonteCarlo(MC)和MolecularDynamics(MD)两大类分子模拟对Boltzmann分布的重要性采样：统计物理是理论基础.5解多体薛定谔方程：HΨ({r})=EΨ({r})方程中包含所有原子核和电子的哈密顿量和波函数原子间相互作用由主要由电子决定计算系统最低能量构型，原子光谱等计算化学反应势垒，速率等Gaussian软件包1.3.第一性计算─电子结构方法6解多体薛定谔方程：HΨ(ρ({r}))=EΨ(ρ({r}))用电子密度代替电子坐标，从而大大减少计算自由度适于进行较大体系的计算VASP软件包仍然难以做有限温度的模拟1.4.第一性计算─密度泛函理论7只计算势能，不用计算力基于Boltzmann分布，只能模拟平衡态体系模拟步长可以很大1.5.蒙特卡罗（MC）模拟8iijijjttFFriiittmFa1iiittttvva11iiittttrrv通过拟合实验数据或第一性计算获得的数据确定经验力场经验力场的精度直接影响模拟结果的质量并行计算环境下可模拟几百万个原子，时间尺度约几十纳秒ijF数值方法求解多体牛顿方程.1.6.全原子分子动力学模拟热耦调节速度以实现恒温模拟9分子动力学模拟软件AMBERCHARMM10目的：把整个原子作为一个质点进行模拟（全原子模拟，AtomisticSimulation或All-atomSimulation，也叫做ForceFieldMethod或MolecularMechanics)，以减少计算自由度，加大可计算的体系的空间和时间尺度，简化数据处理和分析。方法：一般的做法是根据原子间相互作用的物理特性，预先设定一个有待定参数的二体或多体的相互作用的经验势的函数形式，然后根据第一性计算的数据或实验结果拟合经验势的参数。误差：因为描述体系的自由度被大大减少，全原子模型不可能重建系统的所有性质。拟合参数时，往往选择一组最关心的物理性质进行拟合，以求误差尽量小，而放松对其它性质的要求。所以要根据待研究的物理问题适当选取全原子模型。1.7.分子建模11约化单位的换算1.7.1.约化单位（reducedunit）数值模拟时使用的内部单位，需要乘上常数才能对应于实际体系的真实物理单位（国际单位制，SI-units）。1）给定四个基本物理量的单位:长度L，质量M，时间t，电荷电量Q2）计算其它物理量：能量温度压力质量密度数量密度介电常数B/TEk31/nL3/ML2ANQLE22*/EMLt3/PEL其中kB是Boltzmann常数，NA是Avogadro常数。12最常用的描述原子间范德华力的经验势。最广泛使用的是12-6LJ：126126()4Vrrr126126ˆ()242rVrrrrFr惰性气体的原子间相互作用仅用LJ就基本可以完全描述。氩原子之间的相互作用，wikipedia截断距离（cutoffdistance）：对于短程作用，大于cutoff的贡献是常数。24crVrrdr三维空间中，以上积分收敛的为短程作用，发散的为长程作用。1.7.2.Lennard-Jones势13因为金属中的价电子可以自由运动，所以一般要用多体作用描述金属体系的力场。GlueModelEAM(EmbeddedAtomModel)12iijijijijVFrriijijjijVrUr其中rij是两个原子间的距离，是类型为和的原子之间的二体势，是类型为的原子j产生的电子电量密度在i处的值，F是一个嵌入函数，代表把类型为的原子i嵌入电子云中需要的能量。可以用于单一金属或者合金体系。被广泛应用于多种金属及其合金。只适用于单一金属。较好地平衡了表面和内部的结构和能量。1.7.3.金属体系的力场14成键作用（BondedInteractions）：Bonds，ValenceAngles，DihedralAngles(TorsionalAngles),ImproperDihedralAngles非成键作用（NonbondedInteractions）：范德华力和静电力PandeGroupinStanfordAMBERCHARMMOPLSGROMOS1.7.4.化学和生物体系的力场15BondDihedralAngleValenceAngleImproperDihedralAngleKaptonunitEMIm+16类似于从第一性原理层面到全原子层面，粗粒化方法力图从全原子层面进一步简化到粗粒化（coarse-graining，CG）层面，以期大大提高计算的时间和空间尺度。困难在于全原子层面上，原子间相互作用并不集中在局部。而在第一性层面上，电子及其相互作用基本局限在相应的原子核周围。不同的粗粒化方法着重于重建不同的物性，如结构或扩散特性等。一些粗粒化方法假定作用势的函数形式，然后用全原子模拟的结果定参数。另一类从结构函数（如径向分布函数）出发，反推出作用势。我们的方法从全原子作用势出发，通过数学变换较严格地得到粗粒化力场。1.7.5.粗粒化方法17事先做全原子模拟，得到全原子力场。假设中心二体力的粗粒化力场形式。最小化全原子力场和粗粒化力场之间的差值。*W.Noid,P.Liu,Y.Wangetal.J.Chem.Phys.128,244115(2008).从全原子力场出发严格建立粗粒化力场，不预先设定力场的函数形式。大大减少计算自由度。系统扩散过程加快。可以在粗粒化层面去除一些原子自由度（如水分子）保证较好地重建结构性质理论中没有考虑可移植性！粗粒化方法I：MultiscaleCoarse-Graining(MS-CG)*2,,113CGINNIAAICGICGIFFNN18粗粒化方法II：EffectiveForceCoarse-Graining(EF-CG)**Y.Wang,W.Noid,P.Liu,G.A.VothPhys.Chem.Chem.Phys.11,2002(2009).显式计算原子间作用力具有比MS-CG更好的可移植性11MNijijijrr1111ˆˆˆMNijDijRDDDMNijijijDijDrRRRRFRRRRRFrRRRRRrR19势能面（potentialenergysurface）：由不同构型形成的势能的集合。系综（ensemble）：系统在给定宏观条件下所有状态的集合。两个基本假设：等几率原理与各态历经。等几率原理（principleofequalweights）：一个热力学体系有相同的几率访问每一个微观态（注意：不是能量的等几率！一个能量一般会对应很多微观态）。由等几率原理推导得出Boltzmann分布：exp()/jjPEQexp()iBiQEkT其中配分函数（partitionfunction）1.8.平衡态统计物理基本概念各态历经（egodicity）：只要系统演化无穷长时间，总有几率历经势能面上的所有点。即在极限情况下，系综平均和时间平均是等价的。20exp/jjjAAEQ时间平均：分子动力学模拟（MolecularDynamics,MD）系综平均：蒙特卡罗模拟（MonteCarlo,MC）0111limMtitiAAtdtAttM微正则系综(MicrocanonicalEnsemble):NVE皆为常数。正则系综(CanonicalEnsemble):NVT皆为常数。巨正则系综(GrandcanonicalEnsemble):μVT皆为常数，粒子数不固定。等压-等温系综(Isobaric-IsothermalEnsemble):NPT皆为常数。等张力-等温系综(Isotension-IsothermalEnsemble):模拟盒子的形状可变。常用系综21pp1ENiiE=H=E+pV2k112NiiiEmv动能B1ijijijkTNp=VdVfr常用热力学量温度势能21B1NiiiTmdNkv其中d是空间维数焓可以理解为NPT下的有效总内能压强lnBS=kN,V,E熵其中Ω是系统的总微观状态数22BlnF=ETS=kTQHelmholtz自由能G=F+pVE-TS+pVGibbs自由能NPT下的自由能NVT下的自由能T,pT,VGF=NN化学势23空间的连续性：离散模型，如伊辛（Ising）模型，连续模型边界条件：自由、刚性、周期周期性边界条件(PeriodicBoundaryCondition,PBC)：模拟的盒子中的粒子与无穷多的镜像中的粒子有相互作用，从而可以用~103-106个粒子模拟~1023个粒子的体系。totij,1'EL2i,jEnrn1.9.模拟与采样1.9.1.基本概念24特征长度(characteristiclength)：某一特定物理量在空间的相关性的长度。原则上，模拟盒子的边长应该大于所关心的物理量的特征长度。具体操作上，可以通过变化模拟尺寸来了解有限尺度效应(finitesizeeffect)的影响。作用势的截断距离(cutoffdistance)：小于模拟盒子的边长的一半以避免与同一粒子的两个镜像同时作用。有简单截断、截断平移、最小镜像法三种处理方法。采样(sampling)：本质在于在有限时间内进行重要性采样(importancesampling)，即采样对系综平均贡献最大的瞬时量的子集。一般采用均匀时间间隔的采样。初始构型(InitialConfiguration)：尽量接近平衡态。一般需要一段初始的模拟过程以让初始构型达到平衡。在这段初始的模拟过程中不采样。需要某些参数来量化观察系统是否平衡（如液体的体积很容易平衡，势能其次，而扩散系数则较难）。样本的相关度(Correlation)：离得越近的采样样本相关度越大。相关的样本不影响平均值，但是影响误差范围。25Isokineticsthermostat:每一步对速度进行修正，得到需要的温度。Berendsenthermostat：Andersenthermostat:每一步从温度为T的Maxwell-Boltzmann分布中随机产生一个速度赋予被选中的粒子。Nosé-Hooverthermostat:增加额外的自由度，在扩展的微正则系综内实现实际体系的正则系综。1.9.2.热耦（ThermostatorHeatBath)2scaleB031,22iiiiiTNkTmvvvT12011TtTT26周期边界条件下长程力的计算方法。把点电荷用屏蔽电荷分为实空间的短程力和长程力，短程力在实空间做截断，长程力傅里