第五讲分子动力学模拟的Lammps实现课件.ppt

B、使用L-J势模拟裂纹的扩展裂纹lowerupperleftupperleftlowerxyLammps计算输入文件•# 2d LJ crack simulation(问题的基本初始化)•dimension 2•#2维分子动力学模拟•boundary s s p•#x,y方向不存在周期边界条件，z方向周期边界条件•atom_style atomic•#原子类型（金属）•Neighbor 0.3 bin•#建近邻表参数（rc+0.3）bin表示为近邻表类型•neigh_modify delay 5•#间隔多少载荷步重新形成近邻表Lammps计算输入文件•# create geometry创建初始几何构形•Lattice hex 0.93•#指定晶格类型（二维hex）和晶格常数•Region box block 0 100 0 40 -0.25 0.25•#定义一个区域•create_box5 box•#在指定区域建立一个simulation box,5表示原子类型的种类数•create_atoms 1 box•#在simulation box中创建类型为1的原子（原子位置初始化）•Mass 1 1.0•Mass 2 1.0•Mass 3 1.0•Mass 4 1.0•Mass 5 1.0•#指定单个原子的质量Lammps计算输入文件•# LJ potentials(指定原子作用势)•pair_stylelj/cut 2.5•#指定lj势，截断半径为2.5•pair_coeff * * 1.0 1.0 2.5•#指定lj势参数•# define groups（便于加载）•Region 1 block INF INF INF 1.25 INF INF•Group lower region 1•#定义lower组（便于施加外加速度）•Region 2 block INF INF 38.75 INF INF INF•Group upper region 2•#定义upper组（便于施加外加速度）•Group boundary union lower upper•#定义总边界组•Group mobile subtract all boundary•#定义可动原子组（便于统计温度）Lammps计算输入文件•regionleftupper block INF 20 20 INF INF INF•regionleftlower block INF 20 INF 20 INF INF•groupleftupper region leftupper•groupleftlower region leftlower•#定义左上、左下原子组（便于指定裂纹的存在）•setgroup leftupper type 2•setgroup leftlower type 3•setgroup lower type 4•setgroup upper type 5•#指定原子类型（便于指定裂纹的存在）Lammps计算输入文件•# initial velocities初始化速度•computenew mobile temp•#定义温度的计算（可动区域内统计平均）•compute new2 mobile stress/atom •#定义原子应力的计算（整个区域）•Velocity mobile create 0.01 887723 temp new•#按指定的温度（0.01）计算方法，初始化原子的速度•Velocity upper set 0.0 0.3 0.0•#upper原子组y方向的速度为0.3•Velocity mobile ramp vy 0.0 0.3 y 1.25 38.75 sum yes•#mobile原子的速初始度从0到0.3线性变化•# fixes施加约束•fix1 all nve•#nve系综的积分算法•fix2 boundary setforce NULL 0.0 0.0 •#边界boundary上力条件，钢化原子，便于加载！！Lammps计算输入文件•# run运行计算•timestep0.003•#时间间隔步•Thermo 200•#每200步输出热动力学统计量•thermo_modify temp new•#计算温度通过new指示的方法计算•neigh_modify exclude type 2 3•#原子2，3之间作用取消（也就是通过不使他们在近邻表中出现实现）•Dump 1 all atom 500 dump.crack•#每隔500步将原子信息写入文件dump.crack•Dump 2 mobile custom 500 dump2.crack tag x y z c_new2[2]•run5000•#进行5000步的模拟C.使用EAM势模拟Ni的循环行为uuuzxy常温30K条件下金属材料模拟中Lammps的单位•# 3d metal 拉伸模拟•#模拟条件的初始化•Units metal•#指定模拟中的单位类型•boundary s s s•#指定模拟的边界条件•atom_style atomic•#指定原子类型，原子的属性•Lattice fcc 3.52•#指定材料的晶格类型和晶格常数•region box block -5 5 -5 5 -15 15•#xlo,xhi,ylo,yhi,zlo,zhi =区域box的上下限•create_box 1 box•#建立只有1中原子类型的simulation box•create_atoms 1 box•#在simulation box中创建类型为1的原子（原子位置初始化）Lammps计算输入文件•pair_styleeam•#定义原子之间相互作用势•pair_coeff* * Ni_u3.eam•#定义作用势的参数，通过文件输入•Neighbor 2.0 bin•#定义近邻表参数(r+rc),rc=2.0埃•neigh_modify delay 5•#定义更新近邻表的加载步间隔Lammps计算输入文件•regionlower block INF INF INF INF INF -13•#生成lower区域，便于加载•regionupper block INF INF INF INF 13 INF•#生成upper区域，便于加载•grouplower region lower•#定义lower组•groupupper region upper•#定义upper组•groupboundary union lower upper•#定义boundary组•groupmobile subtract all boundary•#定义mobile组，便于计算统计温度（调温）Lammps计算输入文件•#定义计算和边界钢化•compute new2 all centro/atom•#定义计算原子的适配度（针对所有的原子）•compute new3 mobile temp•定义计算统计温度（针对可动原子）•Fix 1 all nve•#时间积分算法（微正则系综）•fix 2 lower setforce NULL NULL 0.0•#使得lower group之间没有z方向的力，z方向钢化 •fix 3 upper setforce NULL NULL 0.0•#使得upper group之间没有z方向的力，z方向钢化Lammps计算输入文件•#定义控温和原子速度初始化•fix 4 mobile temp/rescale 10 30.0 30.0 5.0 1.0•#使用直接调温法调整mobile原子的速度，每10步一调•#使得目标温度为30+-1K，温度在30+-5rescale之外调•fix_modify 4 temp new3•#调温时，计算温度采用news计算•Velocity mobile create 30 887723 temp new3•#初始化mobile group中原子的初速度•#通过温度30K初始化，温度通过new3计算，887723为随机数Lammps计算输入文件•dump 2 all custom 200 d_t.txt tag x y z c_new2•#定义输出，每两百步输出一次•#tag为原子标号•#x，y，z为原子坐标•#c_new2为new2所计算的原子适配度•#这一参数可以判断原子是否是理想晶格•#或者在位错和堆垛层错上•timestep 0.001•#时间步长间隔（ps=10e-12秒）•#时间步长在飞秒量级10e-15秒•#小于原子自振周期的十分之一（0.1ps/10）Lammps计算输入文件Lammps计算输入文件•thermo_style custom step temp pe etotal press vol•#屏幕输出梁的控制，step加载步•#temp温度，pe总势能，etotal总能量•#press为压强,vol为simulation box的体积•thermo100•# 每100步在屏幕上打印统计热动力学量•thermo_modify temp new3•#输出的温度为new3的计算值•#弛豫和循环加载•run 2000•#不加任何载荷计算2000步•velocityupper set 0 0 0.4•velocitylower set 0 0 -0.4•run 10000•velocityupper set 0 0 -0.4•velocitylower set 0 0 0.4•run 20000•velocityupper set 0 0 0.4•velocitylower set 0 0 -0.4•run 10000•#循环载荷模拟Lammps计算输入文件D.使用EAM势模拟Ni的剪切行为vupperlowermobilezxy<100><011><011>Ni FCC Ni FCC 晶格晶格D.使用EAM势模拟Ni的剪切行为Z‘xY’<100><010><001><011><011>yz•# 3d metal shear simulation 3维金属剪切模拟•Units metal•#采用金属材料单位•boundary s s p•#施加z方向的周期边界条件•atom_style atomic•#定义原子之间的相互作用（没有键键作用）•Lattice fcc 3.52•#指定晶格为fcc，晶格常数3.52埃•Region box block 0 16.0 0 10.0 0 2.828427•#定义区域box （x上下限，y上下限，z上下限）•create_box 3 box•#在box区域内创建包含3种原子类型的simulation boxLammps计算输入文件Lammps计算输入文件•Lattice fcc 3.52 orient x 1 0 0 orient y 0 1 1 orient z & 0 -1 1 origin 0.5 0 0•#指定晶格类型、晶格常数、x,y,z坐标的晶向•# &续行符、晶胞的起始点（0.5表示半个晶格常数）•create_atoms 1 box•#在box区域内建立类型为1的原子•pair_styleeam•#定义原子之间作用势为嵌入原子势•pair_coeff * * Ni_u3.eam•#指定嵌入势参数文件Lammps计算输入文件•neighbor0.3 bin•#指定原子近邻列表参数rc+r1中的r1•neigh_modify delay 5•#指定近邻表更新频率•regionlower block INF INF INF 0.9 INF INF•regionupper block INF INF 6.1 INF INF INF•grouplower region lower•groupupper region upper•groupboundary union lower upper•#指定边界组，用于施加边界条件•groupmobile subtract all boundary•#指定可动部分组，用于计算温度Lammps计算输入文件•setgroup lower type 2•setgroup upper type 3•#设置原子类型lower为2，upper为3•Compute new3d mobile temp•#设定温度计算new3d•Compute new2d mobile temp/partial 0 1 1•#设定温度计算new2d，忽略x方向速度•compute new1d all stress/atom •#设定原子应力计算•Velocity mobile create 300.0 5812775 temp new3d•#初始化mobile区域原子的速度•#使new3d方法计算的温度为300k？？ •Fix 1 all nve•#使系综为nve系综，同时确定了方程的积分方法•Fix 2 boundary setforce 0.0 0.0 0.0•#使边界上原子x，y，z方向的受力均为0•Fix 3 mobile temp/rescale 10 300.0 300.0 10.0 1.0•#使用直接调温法，每10步准备调一次•#如果在300+-10k内不调•#如调的话要求其范围在300+-1k之内•fix_modify 3 temp new3d•#fix 3 调温时采用new3d的计算方案Lammps计算输入文件Lammps计算输入文件•#以下为弛豫过程•thermo25•#每25步在屏幕上打印热动力学统计量•thermo_modify temp new3d•#输出的温度按new3d方案计算•timestep0.001•#时间步长为0.001ps•run100•#弛豫100个载荷步Lammps计算输入文件•#以下为剪切加载•velocityupper set 1.0 0 0•#上边界原子施加x方向的速度•velocitymobile ramp vx 0.0 1.0 y 1.4 8.6 sum yes•#在y[1.4,8.6]区域定义0，1之间的vx插值速度•#该速度累加到弛豫的速度上•Unfix 3•#释放弛豫时候的温度约束3•fix 3 mobile temp/rescale 10 300.0 300.0 10.0 1.0•#重新施加温度调控约束3•fix_modify 3 temp new2d•#温度的计算采用new2d计算方案Lammps计算输入文件•#以下为剪切计算•Dump 1 all custom 100 dump.shear tag type & •x y z c_new2[4]•#每100步输出，c_new1d[4]为new1d计算的原子应力•thermo100•thermo_modifytemp new2d•#每次输出屏幕上的温度按new2d方案计算•reset_timestep0•#重新设当前为第0步•run3000•#剪切计算3000步。