ReaxFF的模拟健壮性、性能调优与一点经验
文/石墨 首发/计算化学公社
1、前言
自从前帖《利用Lammps ReaxFF研究反应动力学一例》、《LAMMPS中密度不均一体系域分解问题,和性能调优的一点经验》、《Lammps reaxff物种统计分析工具》发布以来,受到大家的广泛关注,不胜荣幸。在模拟实践中,大家也常常遇到许多问题,有些具有共性,在这里愿共同探讨、抛砖引玉。
ReaxFF本身是效率与精度间的妥协、是分子模拟的一个主要矛盾的冲突点和平衡点。分子模拟处处是矛盾。我们无法消除固有矛盾,只能谦卑地调和矛盾。有一天(可能用不了几年)ReaxFF会被机器学习势、半经验AIMD等等取代,进入矛盾的下一个阶段,不过我们现在仍然有讨论ReaxFF的必要。
2、LAMMPS ReaxFF计算的方式和流程
分子动力学模拟无非是原子小球受周围环境影响进行运动。根据运动方程和相互作用函数(力场),对运动进行演绎。这样,模拟大致的流程如下:
(1)确定原子周围可能有影响的其他原子
(2)基于周围环境,为原子计算受力、加速度,更新速度和位置
(3)每过一段时间,计算和输出一些我们感兴趣的量
其中步骤(1)称为邻居搜索,它的范围我们姑且称作”感知范围”,这一步往往很耗时。步骤(2)也就是我们熟知的受力计算和动力学走步,我们希望性能全都在这里释放。步骤(3)也会造成一定耗时。
3、健壮性
模拟的健壮性主要受3部分影响:
(1)结构的合理性、力场适配性、外界附加条件(温度、外力等)的强度。
(2)动力学走步的步长,它决定模拟的精细度,步长越短,运动也就越平缓、模拟越健壮。
(3)“感知范围”,由两部分组成:(1)键级计算范围,也就是原子考虑与它的邻居发生成键、断键、干涉其他键的范围,这个在control文件中,对有机物为主的体系,默认的4.5 A是个还不错的距离,对于一些无机物,可以考虑加大。键级计算范围对模拟的精度和结果影响很大。(2)缓冲范围,也就是在键级计算范围以外预留的缓冲区的大小,即neigh设定的距离。缓冲越保守(缓冲范围大、缓冲重新计算频率快),越不容易因为瞬时受力、瞬时高速使模拟出错。
4、性能
影响性能的因素有4点:
(1)步长,显而易见。
(2)模拟的密度,密度越大,原子感知范围内的其他原子越多,性能开销越大。
(3)感知范围越大,性能开销显著增大。
(4)额外计算和文件输出的开销。
5、健壮性与性能的妥协
步长直接影响模拟精度,需优先确定,对于温和缓慢的反应,步长可以设为0.25-0.5 fs,高温、大外力、不平衡的体系需要低至0.1 fs。密度也是我们无法决定的。因此,性能调优的关键在于感知范围的调整。
感知范围越小,那么原子的运动越“盲目”,自然模拟不健壮。感知范围太大了,凭空增加开销——一般来说6-8 A以外的原子对当前原子的化学性质影响很小,而且也不会毫无征兆地“侵入当前原子的领地”造成异常。
6、调优策略
(1)确定当前体系适用的步长
(2)根据体系的性质(结构本身是否平衡、温度、外力的大小)确定control文件中键级计算范围,和in文件中neigh的半径和频率。这里我们姑且将短半径、低频率称为“激进”(向健壮性作了妥协,换取了性能),反之则为“保守”。
(3)确定额外计算和文件输出的频率。由于邻居搜索和ReaxFF的物种输出(reax/c/species)、键级输出(reax/c/bonds)都是CPU、内存为瓶颈的计算,而且reax/c/species等有时还会覆盖neigh频率,尽量使后者为邻居搜索频率的整数倍,以使得任务周期同步。例如neigh频率为每50步计算一次,reax/c/species可设为100 5 500(100步计算1次,5次取平均,500步输出一次),轨迹输出为1000步输出一次。
7、上中下三策
(1)激进策略:温度不高/本身反应性不强/模拟已经大差不差开始堆时间,步长0.25-0.5,neigh频率50-100,半径1.5-2.0。
(2)温和策略:步长0.2-0.25,neigh频率20-50,半径2.0-2.5。
(3)保守策略:温度很高/反应性强/不平衡性强/初始结构不太好稍微跑会儿放松下,步长0.1,neigh频率10、半径>=3.0。
8、特殊情况
如果保守策略也跑不了,那90%概率是初始结构不好、力场不适配、条件偏极端等,一般来说3000 K或10 km/s MSST冲击使含能材料爆炸,保守策略也是能正确跑的。此时应该:
(1)检查力场的适用性,力场有的为高温/低温开发,有为稳定聚合物/不稳定的含能材料开发等,假如温和的反应用了含能材料的力场,那必然招致意外的分解;
(2)结构的合理性,如果结构还能救,可以先能量最小化,再用非常保守的策略一点点放松,如0.005 fs跑5000步、0.02 fs跑5000步...渐渐增至0.1 fs。
(3)外力和边界条件施加的合理性,例如外界条件先不要立刻施加,平稳后再一点点施加。另外像unit metal下,输入速度的单位是A / ps,力的单位是eV / A,有时会弄错。周期性重复基底的边界处理、边界有反射壁或LJ壁的也要小心。
9、GPU计算时的设置
由于目前版本LAMMPS的限制,截止2024年初,一般性的解决方案是用单张如Tesla V100一类的,具有双精度浮点的显卡进行计算,CPU使用单进程。LAMMPS需要正确安装KOKKOS模块。此时的运行指令为lmp -in [in文件] -k on g 1 -sf kk -pk kokkos newton on neigh half。1.0密度下,激进策略通常能达到100万原子步/秒。
要试用Tesla V100,可在autodl平台上租几个小时。二手显卡也很便宜。
|