在预处理一些大型、柔性模型时，会遇到表现为”bad contact”、”inf(infinity) energy / force”、”can not be settled”、”step XXX.pdb 刷屏”这些情况，往往是由结构不合理或力场不合适，导致瞬间的力过高、体系崩溃。如果是一些小意外，比如盒子留的边再大1 A就好了、初始建模时自己粗心了，那倒还好。
然而，有时候，在做蛋白、聚合物、某些大型体系时，我们拿到一些结构复杂的低质量模型，其中可能有很多处不合理的接触，这时候可以考虑一些特别的技巧，包括利用软势。
这里以某个球蛋白为例，这个蛋白拿到手就是一个建模精度很差的结构，pdb文件里不含氢、也没有经过优化或者复查。它有约1000个氨基酸，加上氢之后有约3.3万个原子。它的二、三级结构没大问题，但细节上精度很不好。
(, 下载次数 Times of downloads: 13)
当使用gmx pdb2gmx和一些MD相关的可视化软件、工具包加氢后，运行正常的能量最小化，很快报错，提示13639号原子具有“无限大”的受力崩溃。
打开VMD，选取"same residue as (within 3 of index 13639)"以显示周围的情况。果然，这里的两个苯丙氨酸残基有问题。注意画红圈的部位，很明显H原子重叠了，而且上方的N原子也离苯环太近，VMD甚至默认以为C-N之间成键了，而这两个残基并不是直接相连的，蛋白质里面也不可能存在这种乱连的情况。这固然一部分原因是gromacs并没有考虑这种情况，只是简单的按照一般键长、键角加了氢，但本质原因还是这个结构本身精度不好。
(, 下载次数 Times of downloads: 19)
如果是小分子，也许用gaussview还是chemcraft，稍微拉开它们，再做优化也就罢了。然而这个蛋白质有3万多个原子，里面这些“打架”的情况不可能人工处理。如果编程处理，比如将间距<1 A的原子在附近重新指定一个还算合理的位置，也并非不可，但要做到以后也能用，设计算法、调试，恐怕一两天是不够的。
不过，即使VMD以为他们胡乱成键了，gmx pdb2gmx给出的top文件是严格按照序列来的。键、键角、二面角都严格按照残基指定，因此有了拓扑信息，对于MD程序来说，两个苯丙氨酸内部的键、键角、二面角互不干扰，只是LJ和静电部分因为接触太近而给出了过高的值，使得瞬间力和速度过大。
这种时候不妨将“LJ和静电”换掉，也就是换一个pair_style。在LAMMPS中，可以使用pair_style soft。
(, 下载次数 Times of downloads: 11)
LJ势和用于计算静电的Columb势都是距离越近，受力急剧升高的函数。如果原子间距离仅有sigma的1/2，那么它的LJ相互作用强度接近平衡距离时强度的4000倍（2^12-2^6）。可以利用软势（soft）暂且代替它，“温和”地推开距离太近的原子。
具体的做法是：
1、将已经生成的gro文件和top文件输入Intermol，这是一个转换不同MD软件间格式的软件。生成Lammps适用的in文件和data文件，Intermol将它们命名为*_converted.input和*_convert.lmp。我们将其重命名为soft_min.in和soft_min.data。in文件包含模拟设置；data文件包含原子坐标和连接性，力场参数可以写在in文件里，也可以写在data文件里。
2、将in文件里以下内容删除：
- special_bonds lj 0.0 0.0 0.5 coul 0.0 0.0 0.8333
- pair_style lj/cut/coul/long 9.0 9.0
- pair_modify tail yes
- kspace_style pppm 1e-8
- pair_coeff 1 1 0.1700000 3.2500000
- pair_coeff 2 2 0.0157000 1.0690800
- pair_style soft 1.0
- pair_coeff * * 1.0
- dump 1 all atom soft_min.lammpstrj
- minimize 1e-6 1e-6 10000 10000
最终的文件：
- units real
- atom_style full
- boundary p p p
- bond_style hybrid harmonic
- angle_style hybrid harmonic
- dihedral_style hybrid multi/harmonic
- read_data soft_min.data
- pair_style soft 1.0
- pair_coeff * * 1.0
- thermo 50
- thermo_style custom ebond eangle edihed eimp epair pe
- dump 1 all atom 200 soft_min.lammpstrj
- min_style cg
- minimize 1e-6 1e-6 10000 10000
- write_data soft_min.data
单核执行lammps，不到一分钟即完成。这个蛋白本身的键长、键角等是没有多少阻力的，而体系的总势能却直接下降了3个数量级，可见初始结构的原子重叠问题非常严重。
(, 下载次数 Times of downloads: 11)
最后，原先位置的两个苯丙氨酸残基的结构变成了这样：
(, 下载次数 Times of downloads: 14)
而体系总体的二、三级结构不变，其余部分也得到了修正。最后，再利用VMD将LAMMPS的轨迹转为pdb文件，即可以为gromacs正常使用。
(, 下载次数 Times of downloads: 15)
几个可能会关注的问题： 1、如果只是因为LJ和Coulumb在近距离过高，在gmx里将它们的参数（sigma、epsilon、电荷）降低、或直接关闭是否可行？
不可行，这不会改变函数在近距离直线上升的形式。另外，如果直接关闭，那么像残基外围的氢之间就完全没有了斥力，仍然会意外重叠。
2、用软势推开重叠原子之后，还需要换回经典LJ和Coulumb能量最小化吗？
需要，之后就正常做模拟，这一步只相当于代替我们手动精修初始结构。最终还是在真正要用的分子力场中跑。
3、gmx也有soft-core函数（本意用来做自由能计算的），是否也能用？
实测对于重叠比较严重、偏离比较明显的结构，在修正初始结构的情景下，其形式和效果都不如LAMMPS包含的这个cosine函数好用。
4、能用软势做正常模拟吗？它还有什么作用
它没有明确的物理意义，不能够用来做正常模拟。还可以用来优化缠绕的聚合物、将物质“蓬松地”填到空间里（巧妙设定强度和cutoff）、用来操纵结构制作一些特定的示意图等。
5、lammps和gromacs之间结构文件互相转化还有哪些办法常用？
lammpstrj文件本质上就和xyz文件差不多，只含有类型、坐标、电荷，很多软件都能直接读写。如果要包含连接性信息，则需要读写data文件。VMD可以通过topo readlammpsdata和topo writelammpsdata指令读写；OVITO可以读取其他MD程序的结构，并制作data文件；lammps本身也可以输出pdb和xtc等类型的文件。