刘建川・GROMACS 实战 | 拉伸、剪切、压缩实操教程
一、模拟前期体系弛豫
正式开展形变模拟前，必须先完成体系的平衡弛豫：首先搭建有机分子或聚合物的初始模型，采用 NPT 系综，在目标温度与压强下运行足够长时间，让体系密度、结构达到充分平衡状态，这是后续形变模拟结果可靠的前提。
二、形变模拟 mdp 参数配置
通过在 mdp 文件中配置各项异性压力耦合与 deform 参数，即可灵活实现拉伸、压缩、剪切三类形变。控压算法可选用 Parrinello-Rahman 或 Berendsen，压力耦合类型设置为各项异性（anisotropic）；压力张量共 6 个分量，依次对应 xx、yy、zz、xy/yx、xz/zx、yz/zy，因此参考压力 ref_p、压缩率 compressibility均需对应填写 6 个数值。
核心规则：发生形变的方向，压缩率必须设置为 0；形变速率与形变方向由 deform参数的 6 个分量控制，单位为 nm/ps。
2.1 单轴拉伸模拟
以 z 轴方向单轴拉伸为例，形变速率设置为 0.001 nm/ps，对应 z 方向压缩率设为 0，完整参数如下：
pcoupl = Parrinello-Rahman
pcoupltype = anisotropic
tau_p = 1.0
ref_p = 1.0 1.0 1.0 1.0 1.0 1.0
compressibility = 4.5e-5 4.5e-5 0.0 4.5e-5 4.5e-5 4.5e-5
deform = 0.0 0.0 0.001 0.0 0.0 0.0
2.2 单轴压缩模拟
实现压缩仅需将对应方向的 deform 数值改为负值即可。以 z 轴压缩为例，将 z 方向形变速率设为 -0.001 nm/ps，参数调整如下：
pcoupl = Parrinello-Rahman
pcoupltype = anisotropic
tau_p = 1.0
ref_p = 1.0 1.0 1.0 1.0 1.0 1.0
compressibility = 4.5e-5 4.5e-5 0.0 4.5e-5 4.5e-5 4.5e-5
deform = 0.0 0.0 -0.001 0.0 0.0 0.0
同理，修改 deform 对应分量的数值，即可在 x、y 等其他方向实现拉伸或压缩。
2.3 剪切形变模拟
若要实现 x-y 平面内的剪切形变，将 deform 的 xy 分量设置为对应速率，同时将对应剪切方向的压缩率保留正常值，非形变的正方向压缩率设为 0，参数示例如下：
pcoupl = Parrinello-Rahman
pcoupltype = anisotropic
tau_p = 1.0
ref_p = 1.0 1.0 1.0 1.0 1.0 1.0
compressibility = 0.0 0.0 4.5e-5 4.5e-5 4.5e-5 4.5e-5
deform = 0.001 0.001 0.0 0.0 0.0 0.0
形变速率没有固定标准，需要根据自身体系的尺寸、材料特性进行测试调整，保证模拟过程稳定且结果具备物理意义。下图是本案例对某一有机聚合物进行单轴拉伸的结构图：
三、应力数据提取与结果分析
模拟完成后，可使用 gmx energy命令提取各方向应力随时间的变化数据，结合设定的形变速率，即可绘制应力-应变曲线。
需要注意的是，原生 gmx energy计算的应力精度有限，若需更严谨的应力计算结果，可参考相关文献中的计算方法，也可使用定制修改版 GROMACS 工具，如 GROMACS-LS 与 MDStress library，相关工具与教程可参考：https://vanegaslab.org/software
掌握以上参数配置与分析方法，就可以在 GROMACS 中灵活开展有机分子、聚合物体系的拉伸、剪切、压缩力学性能模拟，为材料微观力学机制研究提供可靠的模拟数据支撑。
