使用粗粒化珠簧模型对均聚物进行模拟，但体系在平衡处理后存在大量粒子重叠不报错

jrfjrf123

想请教论坛的大家，

我目前在做一个单分散的，50条链，每条链N=200个珠子的珠簧模型（Bead-spring model）体系的验证实验

力场采用的lammps lj/expand 力场形式：
                        
bond采用harmonic键形式：


angle采用cosine势形式：



其中ε=1， σ=1， Δ=0，K_bond=200, r0=1, K_theta=1.0, 均采用LJ单位，

我借鉴了[size=9.978pt]Rolf Auhl的方法（Auhl, Rolf, et al.The Journal of chemical physics 119.24 (2003): 12718-12728.[size=9.978pt]），先设置rc=2^(1/6)σ的全排斥势，然后进行体系的压缩使得[size=9.978pt]密度达到温度T=1.0状态下的密度和各个[size=9.978pt]分子链的Flory特征比（C）[size=9.978pt]和均方末端距（Ree）接近相应链刚度下的无扰链模型，也就是如下图所示的C随K_theta变化图所示，[size=9.978pt]可见在K_theta=1.0附近其C值应该在2.0附近，通过[size=13.304px]Auhl文章中的计算方法我得到我这个体系在这一步后的C=2.16，满足无扰链的状态。（具体过程见in.balance）
           

但之后在in.annealing模拟文件里按文章里介绍的步骤进行push-off，也就是使用一个soft势函数推开分子链创建排除体系并切换回LJ势，rc设置由1.12变为2.5，平衡一段时间后，体系出现了非常多的重叠情况，而且都是1-3粒子，也就是相邻一个粒子的粒子对间存在非常严重的重叠，下图中红蓝两色的粒子就是高度重叠的间位粒子，图中粒子显示半径为0.3



想请问一下为什么出现了这样高度重叠的情况lammps没有报错而是继续跑完了全程，而且输出的log文件中E_pair 也是一个很小的数值

个人猜测是这些间位原子没有计算LJ势，但不知道是为什么造成了这个原因

具体的计算文件在附件中，计算顺序是先用in.balance计算，再用in.annealing计算前者输出的结构文件
[size=13.304px]

[size=13.304px]

jrfjrf123

问题已经解决了，确实是1-3原子的lj势能没有计算，

由于special_bond 的默认设置是special_bond lj 0 0 0 ， 存在键角拓扑的情况下默认不计算lj势能

将其设置为special_bond lj 0 1 1 就行了