variable step equal 1000 # 500步输出一次 variable cell equal 1 # 定义变量 cell,扩展mxene范围,其值为 0.5 # ==================== 基础设置 ==================== units metal # 皮秒 atom_style charge # 包含电荷的原子类型 dimension 3 # 三维模拟 boundary p p p # 全周期边界 neighbor 2.0 bin neigh_modify every 5 delay 0 check yes # ==================== 数据文件读取 ==================== read_data Ti3C2OH2-5376-OVT.data # 读取mxene数据文件 Reading data file ... triclinic box = (0 0 0) to (150 25 100) with tilt (0 0 0) 1 by 1 by 1 MPI processor grid reading atoms ... 5376 atoms read_data CPU = 0.014 seconds # ==================== 分组 ==================== group mxene type 1 2 3 4 # 定义MXENE原子组 5376 atoms in group mxene # 获取 MXENE 的坐标边界 variable xmax equal bound(mxene,xmax)+${cell} variable xmax equal bound(mxene,xmax)+1 variable xmin equal bound(mxene,xmin)-${cell} variable xmin equal bound(mxene,xmin)-1 variable ymax equal bound(mxene,ymax)+${cell} variable ymax equal bound(mxene,ymax)+1 variable ymin equal bound(mxene,ymin)-${cell} variable ymin equal bound(mxene,ymin)-1 variable zmax equal bound(mxene,zmax)+${cell} variable zmax equal bound(mxene,zmax)+1 variable zmin equal bound(mxene,zmin)-${cell} variable zmin equal bound(mxene,zmin)-1 # 使用计算出的边界值定义mxene区域 region mxene_region block ${xmin} ${xmax} ${ymin} ${ymax} ${zmin} ${zmax} units box region mxene_region block 31.763324 ${xmax} ${ymin} ${ymax} ${zmin} ${zmax} units box region mxene_region block 31.763324 117.998985 ${ymin} ${ymax} ${zmin} ${zmax} units box region mxene_region block 31.763324 117.998985 -0.165619999999999 ${ymax} ${zmin} ${zmax} units box region mxene_region block 31.763324 117.998985 -0.165619999999999 25.44575 ${zmin} ${zmax} units box region mxene_region block 31.763324 117.998985 -0.165619999999999 25.44575 0.342348 ${zmax} units box region mxene_region block 31.763324 117.998985 -0.165619999999999 25.44575 0.342348 33.812089 units box # 计算mxene区域的体积 variable mxene_vol equal (v_xmax-v_xmin)*(v_ymax-v_ymin)*(v_zmax-v_zmin) # 定义 z 坐标范围的两个 region(单位为 box 单位)与对应的组 region 2DMXene-layer block INF INF INF INF 23.00 33.00 units box # 剥离层 MXene region 3DMXene-layer block INF INF INF INF 1.0 22.0 units box # 固定层 MXene group 2DMXenelayer region 2DMXene-layer 1792 atoms in group 2DMXenelayer group 3DMXenelayer region 3DMXene-layer 3584 atoms in group 3DMXenelayer # 在 mxene 中再筛选 group 2DMXene intersect mxene 2DMXenelayer 1792 atoms in group 2DMXene group 3DMXene intersect mxene 3DMXenelayer 3584 atoms in group 3DMXene # ==================== 创建水分子 ==================== # 定义水分子区域 region box_l block INF ${xmin} INF INF INF INF units box # mxene左侧区域 region box_l block INF 31.763324 INF INF INF INF units box region box_r block ${xmax} INF INF INF INF INF units box # mxene右侧区域 region box_r block 117.998985 INF INF INF INF INF units box region box_h block INF INF INF INF ${zmax} INF units box # mxene上侧区域 region box_h block INF INF INF INF 33.812089 INF units box region water_region union 3 box_h box_l box_r # 并集区域 molecule water_molecule water.txt # 读取水分子数据文件 Read molecule template water_molecule: # Water molecule. TIP3P geometry 1 molecules 0 fragments 3 atoms with max type 2 0 bonds with max type 0 0 angles with max type 0 0 dihedrals with max type 0 0 impropers with max type 0 # 0表示分子模板中的type种类数加0. create_atoms 0 random 2000 67895 water_region mol water_molecule 12345 overlap 3 maxtry 1000 Created 6000 atoms using lattice units in triclinic box = (0 0 0) to (150 25 100) with tilt (0 0 0) create_atoms CPU = 0.939 seconds group water region water_region 5999 atoms in group water # 计算 mxene 组的原子数目 variable mxene_count equal count(mxene) # 输出 mxene 组的原子数目到日志文件 print "Number of atoms in mxene group: ${mxene_count}" Number of atoms in mxene group: 5376 # 定义一个初速度区域和组,用于mxene初期的剥离 region vlayer block 2.5 4.5 INF INF INF INF units box group forceLayerBlock region vlayer 88 atoms in group forceLayerBlock group forcelayer intersect 2DMXene forceLayerBlock 0 atoms in group forcelayer # ———— 定义 water 组密度 ———— # 1) 先统计 water 组原子数 variable cnt equal count(water) # 2) 计算水分子密度(n/V) variable waterden equal v_cnt/(vol-v_mxene_vol) group mobile union water 2DMXene 7791 atoms in group mobile # ==================== 力场初始化 ==================== # ReaxFF力场配置(动态键处理) pair_style reaxff NULL safezone 50 mincap 500 minhbonds 200 # GPU加速计算 pair_coeff * * am0c17536_si_002.txt C H O Ti # 原子类型顺序必须与data文件严格一致 WARNING: Van der Waals parameters for element LI indicate inner wall+shielding, but earlier atoms indicate a different van der Waals method. This may cause division-by-zero errors. Keeping van der Waals setting for earlier atoms. (src/REAXFF/reaxff_ffield.cpp:255) WARNING: Changed valency_val to valency_boc for F (src/REAXFF/reaxff_ffield.cpp:300) WARNING: Changed valency_val to valency_boc for X (src/REAXFF/reaxff_ffield.cpp:300) fix qe all qeq/reaxff 1 0.0 10.0 1.0e-6 reaxff # 电荷平衡 write_data mxene-h2o.data # 写入数据文件 System init for write_data ... Last active /omp style is pair_style reaxff/omp Neighbor list info ... update: every = 5 steps, delay = 0 steps, check = yes max neighbors/atom: 2000, page size: 100000 master list distance cutoff = 12 ghost atom cutoff = 12 binsize = 6, bins = 25 5 17 2 neighbor lists, perpetual/occasional/extra = 2 0 0 (1) pair reaxff/omp, perpetual attributes: half, newton off, ghost, omp pair build: half/bin/newtoff/ghost/omp stencil: full/ghost/bin/3d bin: standard (2) fix qeq/reaxff/omp, perpetual, copy from (1) attributes: half, newton off pair build: copy stencil: none bin: none # ==================== 能量最小化 ==================== min_style cg # 共轭梯度法 min_modify dmax 0.01 # 限制最大位移步长(防止原子飞散) minimize 1e-4 1e-6 5000 10000 CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE Your simulation uses code contributions which should be cited: - pair reaxff command: doi:10.1016/j.parco.2011.08.005 @Article{Aktulga12, author = {H. M. Aktulga and J. C. Fogarty and S. A. Pandit and A. Y. Grama}, title = {Parallel Reactive Molecular Dynamics: {N}umerical Methods and Algorithmic Techniques}, journal = {Parallel Computing}, year = 2012, volume = 38, number = {4--5}, pages = {245--259} } - pair reaxff/omp and fix qeq/reaxff/omp command: doi:10.1177/1094342017746221 @Article{Aktulga17, author = {H. M. Aktulga and C. Knight and P. Coffman and K. A. O'Hearn and T. R. Shan and W. Jiang}, title = {Optimizing the Performance of Reactive Molecular Dynamics Simulations for Multi-Core Architectures}, journal = {International Journal of High Performance Computing Applications}, year = 2018 } CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE Switching to 'neigh_modify every 1 delay 0 check yes' setting during minimization Last active /omp style is pair_style reaxff/omp Per MPI rank memory allocation (min/avg/max) = 2.044e+04 | 2.044e+04 | 2.044e+04 Mbytes Step Temp E_pair E_mol TotEng Press 0 0 -1155140.9 0 -1155140.9 521947.43 31 0 -1198370.1 0 -1198370.1 -93247.228 Loop time of 4.45901 on 24 procs for 31 steps with 11376 atoms 1185.8% CPU use with 1 MPI tasks x 24 OpenMP threads Minimization stats: Stopping criterion = energy tolerance Energy initial, next-to-last, final = -1155140.88209652 -1198259.04599548 -1198370.10695099 Force two-norm initial, final = 10209.811 366.22545 Force max component initial, final = 133.91499 12.847462 Final line search alpha, max atom move = 0.00072554811 0.0093214517 Iterations, force evaluations = 31 31 MPI task timing breakdown: Section | min time | avg time | max time |%varavg| %total --------------------------------------------------------------- Pair | 3.0494 | 3.0494 | 3.0494 | 0.0 | 68.39 Neigh | 0 | 0 | 0 | 0.0 | 0.00 Comm | 0.001765 | 0.001765 | 0.001765 | 0.0 | 0.04 Output | 0 | 0 | 0 | 0.0 | 0.00 Modify | 1.4027 | 1.4027 | 1.4027 | 0.0 | 31.46 Other | | 0.005066 | | | 0.11 Nlocal: 11376 ave 11376 max 11376 min Histogram: 1 0 0 0 0 0 0 0 0 0 Nghost: 19772 ave 19772 max 19772 min Histogram: 1 0 0 0 0 0 0 0 0 0 Neighs: 2.09983e+06 ave 2.09983e+06 max 2.09983e+06 min Histogram: 1 0 0 0 0 0 0 0 0 0 Total # of neighbors = 2099832 Ave neighs/atom = 184.58439 Neighbor list builds = 0 Dangerous builds = 0 # 收敛标准:能量变化<1e-4 kcal/mol, 力<1e-6 kcal/(mol·Å) # 最大迭代步数:5000次,最大力评估次数:10000次 # ==================== 结构优化 ==================== write_data minimized.data # 优化后结构 System init for write_data ... Last active /omp style is pair_style reaxff/omp write_restart minimized.restart System init for write_restart ... Last active /omp style is pair_style reaxff/omp # ==================== 热平衡 ==================== # --- 第一阶段:MXene 热平衡 --- # 短时 NVT 让系统温度分布合理 reset_timestep 0 thermo ${step} thermo 1000 thermo_style custom step temp pe etotal press vol # 对 Group mxene 施加 NVT velocity mxene create 300.0 123456 rot yes dist gaussian fix mynvt mxene nvt temp 300 300 0.1 # 输出 MXene 热平衡轨迹 dump balance all custom ${step} mxene_balance.lammpstrj id type x y z dump balance all custom 1000 mxene_balance.lammpstrj id type x y z timestep 0.0005 run 2000 Last active /omp style is pair_style reaxff/omp Per MPI rank memory allocation (min/avg/max) = 2.044e+04 | 2.044e+04 | 2.044e+04 Mbytes Step Temp PotEng TotEng Press Volume 0 141.75824 -1198370.1 -1198161.7 -92655.38 375000 1000 nan -466881.63 nan nan 375000 2000 nan -466881.63 nan nan 375000 Loop time of 92.0266 on 24 procs for 2000 steps with 11376 atoms Performance: 0.939 ns/day, 25.563 hours/ns, 21.733 timesteps/s, 247.233 katom-step/s 898.3% CPU use with 1 MPI tasks x 24 OpenMP threads MPI task timing breakdown: Section | min time | avg time | max time |%varavg| %total --------------------------------------------------------------- Pair | 57.892 | 57.892 | 57.892 | 0.0 | 62.91 Neigh | 0.013463 | 0.013463 | 0.013463 | 0.0 | 0.01 Comm | 0.12457 | 0.12457 | 0.12457 | 0.0 | 0.14 Output | 0.012196 | 0.012196 | 0.012196 | 0.0 | 0.01 Modify | 33.972 | 33.972 | 33.972 | 0.0 | 36.92 Other | | 0.01219 | | | 0.01 Nlocal: 11376 ave 11376 max 11376 min Histogram: 1 0 0 0 0 0 0 0 0 0 Nghost: 19775 ave 19775 max 19775 min Histogram: 1 0 0 0 0 0 0 0 0 0 Neighs: 2.13295e+06 ave 2.13295e+06 max 2.13295e+06 min Histogram: 1 0 0 0 0 0 0 0 0 0 Total # of neighbors = 2132947 Ave neighs/atom = 187.49534 Neighbor list builds = 2 Dangerous builds = 0 unfix mynvt undump balance # --- 第二阶段:对水分子做速度初始化并保存平衡结构 --- reset_timestep 0 thermo ${step} thermo 1000 thermo_style custom step temp pe etotal press vol v_waterden # 为 water 组生成初始速度:300 K,高斯分布,去除质心动量 velocity water create 300.0 123456 dist gaussian mom yes fix mynve1 water nve fix mynve2 mxene nve # 计算 water 组的瞬时温度,并将其存储在计算变量 tpw 中。 compute tpw0 water temp # 对 water 组的原子使用 Berendsen 温度耦合器进行温度控制。300 300:初始和目标温度都为 300K。 fix wbr0 water temp/berendsen 300 300 0.1 # 修改 wbr 固定操作,使其使用计算出的温度 tpw。即温度控制会基于 compute tpw 计算的水分子的瞬时温度,而不是全局温度 fix_modify wbr0 temp tpw0 timestep 0.001 run 3000 Last active /omp style is pair_style reaxff/omp ERROR on proc 0: Non-numeric atom coords - simulation unstable (src/OPENMP/domain_omp.cpp:58) Last command: run 3000