|
LAMMPS讲解35-ReaxFF力场之Qeq电荷计算参数验证ReaxFF已经存在很多的力场参数。很多时候能够找到包含所需元素的力场参数。但是力场参数往往是根据特定体系而开发的，是否适用于自己的体系则需要验证力场参数。验证力场参数的基本流程就是对比自己体系中的分子各种能量（如键，键角，二面角）与DFT计算的能量。如果吻合较好那么找到的力场参数就可以用，如果不好就要基于找到的力场参数进行进一步优化。本节介绍力场参数的验证过程。严格的验证力场需要对比各项能量和参数，包括键结、非键结和电荷平衡参数。本例介绍如何进行电荷平衡参数验证 这个可以通过计算一个完整分子和分子片段内原子的电荷，对比DFT和ReaxFF的计算结果。以合成酯为例，合成酯中只包含CHO三种元素。可以使用这篇文章（Chowdhury Ashraf and Adri C.T. van Duin* Extension of the ReaxFF Combustion Force Field toward Syngas Combustion and Initial Oxidation Kinetics J. Phys. Chem. A 121, 5, 1051-1068）中所开发出的力场参数。把力场参数下载下来，命名为ffield2016.reax.cho。目前在LAMMPS中只用Qeq方法进行电荷动态平衡。该方法要求有三个参数分别是 chi eta gamma 对于LAMMPS来说，可以直接采用ReaxFF中力场文件中的参数，也可以自己单独写一个文件定义这些参数。通常我们会直接采用力场文件中的参数。这里要验证的由lammps的例子中给出的电荷平衡参数。对于合成酯采用以下几个分子进行验证。
为什么要选择这几种。这就需要对研究对象有充分的了解。这里要研究合成酯的热解问题，那么首先合成酯的完整分子自然要考虑，一半合成酯也可以考虑一下。然后根据前人研究合成酯会发生以下热解过程，那么热解过程中的反应产物也要考虑。而不含O的自由基和分子就不用考虑了，因为待验证的力场在开发时已经充分验证了。
这里采用LAMMPS计算ReaxFF所得到的电荷并与采用量子化学的电荷进行对比。要计算首先要有输入文件的坐标。可以采用Avogadro进行坐标生成。这里没有采用MS是因为MS太大了，这么个小事不需要这么大的软件。使用Avogadro分别生成上述三个分子的pdb文件和xyz文件。在Avogadro中画完分子之后一定要做几下结构优化只有再导出文件。Qeq方法只需要分子几何构型和原子的实验性质就可以预测分子内部的每个原子的电荷。先将导出的pdb文件中转换为data文件（一个分子就够了），然后在lammps中计算每个原子的部分电荷。计算in文件如下 units real atom_style full neighbor 2.0 bin neigh_modify every 1 delay 0 check no read_data half.data pair_style reaxff lmp_control pair_coeff * * ffield2016.reax.cho C H O timestep 0.25 velocity all create 300.0 4928459 rot yes dist gaussian thermo 100 thermo_style custom step temp pe ke thermo_modify flush yes fix 1 all nve fix 2 all temp/berendsen 300.0 300.0 50.0 fix 3 all qeq/reaxff 1 0.0 10.0 1e-6 param.qeq run 100000 dump 1 all custom 100 dump.lammpstrj id mol type q x y z dump_modify 1 sort id run 100000 这样就能得到由ReaxFF计算得到的电荷。然后使用导出的xyz文件，调用Multiwfn软件和Gaussian程序使用量子化学计算得到的电荷。这里我们采用1.2CM5方法计算电荷。计算的采用的量子化学信息为B3LYP-D3(BJ)/def2-SVP。由于1.2*CM5电荷就是在Truhlar等人提出的原始的CM5电荷基础上乘上1.2，这相当于增大了原子电荷的数量级，等效地体现出了溶剂环境对溶质的极化作用。这里的模拟没有溶剂化作用，所以需要对得到的电荷除以1.2得到CM5电荷。关于如何联用Multiwfn软件和Gaussian程序计算CM5电荷，参见网址计算适用于OPLS-AA力场做模拟的1.2*CM5原子电荷的懒人脚本 - 思想家公社的门口：量子化学·分子模拟·二次元 (sobereva.com)。这里要求Gaussian程序是Linux版下的。如需获取Linux版的Gaussian可加微信baolu_yao。 以下是上述三个分子中ReaxFF计算得到的原子电荷和量子化学得到的原子电荷对比。从三个图中可以看到这里所采用的Qeq参数以很高的精度再现了量子化学计算得到的电荷，表明采用Qeq参数是合理的。具体的Qeq参数见下表。
感谢鲍路瑶老师的分享，内容来自于鲍老师分享出来的资料 如有需要添加微信：lmp_zhushou 进入微信群，帮助他人，共建社区 获取完整版lammps讲义可以加微信lmp_zhushou或加入QQ群994359511 |
