分子动力学模拟的Lammps实现

第五讲分子动力学模拟的Lammps实现A.常用分子动力学模拟程序•NAMD：免费，主要针对与生物和化学软材料体系，程序设计水平高，计算效率高。有很好的分析辅助软件VMD。•AMBER主要针对生物体系，也适当兼容一般化学分子。有很好的内置势能模型，自定义新模型和新分子很方便，有很完善的维护网站。计算效率不高运算速度慢。•CHARMM主要针对生物体系，也包含部分化学体系。势能模型更新很快自定义新模型比较方便。计算效率低。•GROMACS免费主要针对生物体系，也适当照顾一般化学体系。算法好，计算效率高。界面友好，维护服务好。•TINKER免费一般性分子动力学软件，对生物体系略有偏重。优点支持多种模型。仍在开发中，某些方面还不完善。•DL-POLY一般性分子模拟软件，界面友好，计算效率高。维护服务很好.A.常用分子动力学模拟程序•LAMMPS免费一般性分子模拟软件。兼容当前大多数的势能模型，编程水平高，计算效率高。可以模拟软材料和固体物理系统。•MaterialsExplorer立足于Windows平台的多功能分子动力学软件。拥有强大的分子动力学计算及MonteCarlo软件包，是结合应用领域来研究材料工程的有力工具。MaterialsExplorer可以用来研究有机物、高聚物、生物大分子、金属、陶瓷材料、半导体等晶体、非晶体、溶液，流体，液体和气体相变、膨胀、压缩系数、抗张强度、缺陷等。MaterialsExplorer软件中包含2Body，3Body，EAM，AMBER等63个力场可供用户选择。MaterialsExplorer软件拥有完美的图形界面，方便使用者操作。A.常用分子动力学模拟程序B、使用L-J势模拟裂纹的扩展•L-J无量纲化计算中的单位B、使用L-J势模拟裂纹的扩展裂纹lowerupperleftupperleftlowerxyLammps计算输入文件•#2dLJcracksimulation(问题的基本初始化)•dimension2•#2维分子动力学模拟•boundaryssp•#x,y方向不存在周期边界条件，z方向周期边界条件。•atom_styleatomic•#原子类型（金属）•Neighbor0.3bin•#建近邻表参数（rc+0.3）bin表示为近邻表类型。•neigh_modifydelay5•#间隔多少载荷步重新形成近邻表Lammps计算输入文件•#creategeometry创建初始几何构形•Latticehex0.93•#指定晶格类型（二维hex）和晶格常数•Regionboxblock0100040-0.250.25•#定义一个区域•create_box5box•#在指定区域建立一个simulationbox,5表示原子类型的种类数•create_atoms1box•#在simulationbox中创建类型为1的原子（原子位置初始化）•Mass11.0•Mass21.0•Mass31.0•Mass41.0•Mass51.0•#指定单个原子的质量Lammps计算输入文件•#LJpotentials(指定原子作用势)•pair_stylelj/cut2.5•#指定lj势，截断半径为2.5•pair_coeff**1.01.02.5•#指定lj势参数•#definegroups（便于加载）•Region1blockINFINFINF1.25INFINF•Grouplowerregion1•#定义lower组（便于施加外加速度）•Region2blockINFINF38.75INFINFINF•Groupupperregion2•#定义upper组（便于施加外加速度）•Groupboundaryunionlowerupper•#定义总边界组•Groupmobilesubtractallboundary•#定义可动原子组（便于统计温度）Lammps计算输入文件•regionleftupperblockINF2020INFINFINF•regionleftlowerblockINF20INF20INFINF•groupleftupperregionleftupper•groupleftlowerregionleftlower•#定义左上、左下原子组（便于指定裂纹的存在）•setgroupleftuppertype2•setgroupleftlowertype3•setgrouplowertype4•setgroupuppertype5•#指定原子类型（便于指定裂纹的存在）Lammps计算输入文件•#initialvelocities初始化速度•computenewmobiletemp•#定义温度的计算（可动区域内统计平均）•computenew2mobilestress/atom•#定义原子应力的计算（整个区域）•Velocitymobilecreate0.01887723tempnew•#按指定的温度（0.01）计算方法，初始化原子的速度•Velocityupperset0.00.30.0•#upper原子组y方向的速度为0.3•Velocitymobilerampvy0.00.3y1.2538.75sumyes•#mobile原子的速初始度从0到0.3线性变化•#fixes施加约束•fix1allnve•#nve系综的积分算法•fix2boundarysetforceNULL0.00.0•#边界boundary上力条件，钢化原子，便于加载！！Lammps计算输入文件•#run运行计算•timestep0.003•#时间间隔步•Thermo200•#每200步输出热动力学统计量•thermo_modifytempnew•#计算温度通过new指示的方法计算•neigh_modifyexcludetype23•#原子2，3之间作用取消（也就是通过不使他们在近邻表中出现实现）•Dump1allatom500dump.crack•#每隔500步将原子信息写入文件dump.crack•Dump2mobilecustom500dump2.cracktagxyzc_new2[2]•run5000•#进行5000步的模拟C.使用EAM势模拟Ni的循环行为uuuzxy常温30K条件下金属材料模拟中Lammps的单位•#3dmetal拉伸模拟•#模拟条件的初始化•Unitsmetal•#指定模拟中的单位类型•boundarysss•#指定模拟的边界条件•atom_styleatomic•#指定原子类型，原子的属性•Latticefcc3.52•#指定材料的晶格类型和晶格常数•regionboxblock-55-55-1515•#xlo,xhi,ylo,yhi,zlo,zhi=区域box的上下限•create_box1box•#建立只有1中原子类型的simulationbox•create_atoms1box•#在simulationbox中创建类型为1的原子（原子位置初始化）Lammps计算输入文件•pair_styleeam•#定义原子之间相互作用势•pair_coeff**Ni_u3.eam•#定义作用势的参数，通过文件输入•Neighbor2.0bin•#定义近邻表参数(r+rc),rc=2.0埃•neigh_modifydelay5•#定义更新近邻表的加载步间隔Lammps计算输入文件•regionlowerblockINFINFINFINFINF-13•#生成lower区域，便于加载•regionupperblockINFINFINFINF13INF•#生成upper区域，便于加载•grouplowerregionlower•#定义lower组•groupupperregionupper•#定义upper组•groupboundaryunionlowerupper•#定义boundary组•groupmobilesubtractallboundary•#定义mobile组，便于计算统计温度（调温）Lammps计算输入文件•#定义计算和边界钢化•computenew2allcentro/atom•#定义计算原子的适配度（针对所有的原子）•computenew3mobiletemp•定义计算统计温度（针对可动原子）•Fix1allnve•#时间积分算法（微正则系综）•fix2lowersetforceNULLNULL0.0•#使得lowergroup之间没有z方向的力，z方向钢化•fix3uppersetforceNULLNULL0.0•#使得uppergroup之间没有z方向的力，z方向钢化Lammps计算输入文件•#定义控温和原子速度初始化•fix4mobiletemp/rescale1030.030.05.01.0•#使用直接调温法调整mobile原子的速度，每10步一调•#使得目标温度为30+-1K，温度在30+-5rescale之外调•fix_modify4tempnew3•#调温时，计算温度采用news计算•Velocitymobilecreate30887723tempnew3•#初始化mobilegroup中原子的初速度•#通过温度30K初始化，温度通过new3计算，887723为随机数Lammps计算输入文件•dump2allcustom200d_t.txttagxyzc_new2•#定义输出，每两百步输出一次•#tag为原子标号•#x，y，z为原子坐标•#c_new2为new2所计算的原子适配度•#这一参数可以判断原子是否是理想晶格•#或者在位错和堆垛层错上•timestep0.001•#时间步长间隔（ps=10e-12秒）•#时间步长在飞秒量级10e-15秒•#小于原子自振周期的十分之一（0.1ps/10）Lammps计算输入文件Lammps计算输入文件•thermo_stylecustomsteptemppeetotalpressvol•#屏幕输出梁的控制，step加载步•#temp温度，pe总势能，etotal总能量•#press为压强,vol为simulationbox的体积•thermo100•#每100步在屏幕上打印统计热动力学量•thermo_modifytempnew3•#输出的温度为new3的计算值•#弛豫和循环加载•run2000•#不加任何载荷计算2000步•velocityupperset000.4•velocitylowerset00-0.4•run10000•velocityupperset00-0.4•velocitylowerset000.4•run20000•velocityupperset000.4•velocitylowerset00-0.4•run10000•#循环载荷模拟Lammps计算输入文件D.使用EAM势模拟Ni的剪切行为vupperlowermobilezxy100011011NiFCC晶格D.使用EAM势模拟Ni的剪切行为Z‘xY’100010001011011yz•#3dmetalshearsimulation3维金属剪切模拟•Unitsmetal•#采用金属材料单位•boundaryssp•#施加z方向的周期边界条件•atom_styleatomic•#定义原子之间的相互作用（没有键键作用）•Latticefcc3.52•#指定晶格为fcc，晶格常数3.52埃•Regionboxblock016.0010.002.828427•#定义区域box（x上下限，y上下限，z上下限）•crea