ReaxFF-弹性墙壁：添加外电场、弹性墙壁

LAMMPS中，外加电场时，只考虑电场对原子的“受力”的影响，不考虑对能量（进而对化学反应）的影响。在AMS中，电场的影响不仅仅体现在力上，也体现在能量上。

模拟体系存在外加电场的情况下，如果有原子穿越Cell的边界，就会出现EEM电荷（EEM：ReaxFF计算原子电荷的方法模型）unphysically large相关报错。为了避免穿越Cell边界，因此需要设置一个弹性墙壁，限制原子的活动范围。当然这就变成了一个非周期性的体系，为了避免分子扩散，同理也是需要这个弹性墙壁的。

下面以一个简单的例子，演示如何添加外电场，以及如何添加弹性墙壁。思路如下：首先创建一个常规的盒型体系，去掉周期性（Main中周期性设置为None），然后将其修改为球形，设置外加电场后，设置球形弹性墙壁，球半径略大于模型半径。

本教程适用于AMS2022及其以上版本，2019及以前的版本，请参考：在ReaxFF中设置阻碍原子穿越的弹性墙壁与外加电场

注意

无论是2019版还是新版，凡是添加外电场，都是一个非周期性的体系，区别是旧版是一个方形的，新版是一个球形的。旧版虽然是一个Cell，而且Periodicity设置为Bulk，但是不允许原子穿越Cell边界，实际上就是一个非周期体系。

建模

这里以均相体系H₂O和O₂混合物为例，由于弹性墙壁是球形（而非旧版的正交盒形），因此我们在创建均相体系之后，可以把模型修改为球形（这样能确保分子扩散后，密度尽量不发生大的变化）。当然非均相体系也是可以的，只是需要注意球半径的设置。

创建均相混合物

首先我们创建H₂O和O₂的混合物，AMSinput → Edit → Builder：

注意设置好分子种类与个数之后，窗口底部会显示密度，如果没有问题，点击Generate即可生成该体系。

将体系修改为球形

Main → Periodicity设置为None。菜单栏View → Axes显示坐标轴，选中尽量接近坐标原点的某个原子，然后Select → Select within Radius，输入一个半径例如12 Å（建议这个值略小于Cell直径的一半，确保球在Cell内），OK。这样一个半径为12 Å的一个球被选中了，此时球表面有的分子只有部分原子被选中，我们希望完整的分子都被选中，因此我们按Ctrl m键，则相关完整分子都被选中了（也意味着半径略有扩大）。注意，确保球心为坐标原点非常重要，否则有可能分子在弹性墙壁之外！

说明

上述操作的目的，是因为直接取消周期性，那么Cell边界的分子将会由于周期性的失去，而变成半截分子碎片，因此创建一个直径显著小于边长的球，保留下来的球形体系的所有分子都是完整的分子。

Select - Invert Selection反选，然后del键删除，这样留下来的，就是刚才半径约为12 Å的一个球：

那么球形弹性墙壁，大致可以设置为12 Å多一些（如果弹性墙壁半径没有包纳全部原子，运行将报错Atoms too close）。本例中，设为17埃才正常运行。

设置外加电场

Model → Electric Field，三个数值是电场沿x、y、z方向的分量，例如这里我们设置z方向电场强度为0.01V/Å（注意V/Å是一个非常大的单位，1Å的距离将会产生1V的电压差）

该设置，在生成的*.run文件中，体现在System关键词中

System
......省略
   ElectrostaticEmbedding
      ElectricField 0.0 0.0 0.01
   End
......省略
End

设置球形弹性墙壁

照常设置基本分子动力学参数，例如力场、系综、温度、步数、保存轨迹的频率等等，可以参考其他教程。设置完毕后保存作业，在生成的*.run文件中，我们添加如下字段：

EngineAddons
   WallPotential
      Enabled Yes
      Gradient 5.0
      Prefactor 1.0
      Radius 17
   End
End

这里Gradient是描述弹性墙壁的Potential强度变化梯度，单位是1/Å，这个值设置的越小，墙壁越陡峭，可以理解为墙壁的厚度为这个值。因此一般建议这个值设置的很小例如5。这里Prefactor是Potential的指前因子，Radius是以坐标原点为中心的球半径，单位为Å，这里我们设置为17 Å，足以将之前的原子全部包围。

这个字段可以添加在*.run文件任何字段之间，例如在MolecularDynamics字段后面（该例子是一个NVT系综）：

MolecularDynamics
    NSteps 50000
    InitialVelocities
        Temperature 298
    End
    Thermostat
        Type NHC
        Temperature 298
        Tau 5
    End
End
EngineAddons
   WallPotential
      Enabled Yes
      Gradient 5.0
      Prefactor 1.0
      Radius 17
   End
End
System
   ElectrostaticEmbedding
      ElectricField 0.0 0.0 0.01
   End
......省略
End

保存*.run文件，并在AMSJobs窗口运行该作业（或在命令行运行该*.run文件）。

结果检验

SCM → Movie：

可以看到经过2 ps以后，分子仍然被限定在一个非常狭小的范围，没有扩散到真空中。

补充

有的文献中，添加了正弦、方波电场，例如：利用分子动力学模拟研究微波辐射对焦炭气化机理的非热效应（Renewable Energy 2023，应该是使用老版本（2019.1及之前的版本）添加的外电场（教程）。

老版本的格式里面：

  cat > eregime.in <<eor
  #Electric field regimes
  #start #V direction Magnitude(V/Angstrom)
  0000 1 x  0.001
  1000 1 x -0.001
  2000 1 y  0.001
  3000 1 y -0.001
  4000 2 x -0.001  y -0.001
  5000 2 x  0.001  y  0.001
  eor

上述表示：

第0步开始，外加电场沿x方向，强度为0.001V/Angstrom
第1000步，电场变为沿x方向，强度为-0.001V/Angstrom
第2000步，电场变为沿y方向，强度为0.001V/Angstrom
第3000步，电场变为沿y方向，强度为-0.001V/Angstrom
第4000步，电场变为x方向分量-0.001V/Angstrom；y方向分量-0.001V/Angstrom
第5000步，电场变为x方向分量0.001V/Angstrom；y方向分量0.001V/Angstrom
其中第二列数字是表示“分量”数，例如xy两个分量，就是2，只有x一个分量，就是1，如果xyz三个分量，则是3

该格式允许用户设置每个step时，电场的大小。因此用户实际上可以通过其他工具生成需要的正弦函数，横坐标为步数，纵坐标为场强，然后替换到这部分内容就可以了。

目录