[发明专利]一种随机分批埃瓦德算法的近场加速方法

说明书

本发明公开了一种随机分批埃瓦德算法的近场加速方法，涉及分子动力学模拟技术领域。在分子动力学模拟埃瓦德算法中，对于长程部分，在傅里叶空间中对频率随机采样，计算对应的受力；对于短程部分，采用双层近邻列表方法。本发明利用随机方法近似受力，在对体系模拟的数值积分过程中，随机力波动的期望为零且方差可控，使得最终分布接近实际分布，因此具有良好的收敛性，计算复杂度为线性，比现有方法更快；只需要在CPU核之间通讯有限个频率样本，不需要密集通讯，并行效率高；对小截断半径内的力的作用准确计算，而对球壳内粒子随机选择并计算受力，减少CPU内存占用和计算量。

技术领域

本发明涉及分子动力学模拟技术领域，尤其涉及一种随机分批埃瓦德算法的近场加速方法。

背景技术

在分子动力学模拟中，应用的场景一般是整体电中性的带电离子体系，外加周期性边界条件，一般使用埃瓦德算法(Ewald算法)来处理这个具有长程相互作用的问题。Ewald算法是由德国物理学家保罗·彼得·埃瓦尔德(Paul Peter Ewald)发明的一种计算周期性系统静电相互作用势能的方法。其核心想法是把电荷间的库仑相互作用拆分成长程和短程两个部分，其中长程部分通过傅里叶变换快速计算；而短程部分和伦纳德琼斯(Lennard-Jones)势能都是快速衰减的项，通常使用截断方法进行处理；另外还有原子间共价键、键角等相互作用。伦纳德琼斯势能是由数学家约翰·莱纳德·琼斯(John Lennard-Jones)提出的用来描述惰性气体分子间的相互作用。具体的算法如下：

模拟体系包含N个离子，整体满足电中性条件：

外加周期性边界条件，可以得到体系整体的相互作用势能：

其中前三项都是粒子间的键能，分别代表键长、键角、以及多个键角形成的二面角所对应的能量，这个部分的能量只与相邻粒子有关。后两项分别表示库仑相互作用以及伦纳德琼斯势能。其中n∈Z³表示三维整数向量，∑_n'表示排除i＝j并且n＝0的情形。对于库仑相互作用，直接使用截断是不准确的，经典的Ewald方法是将拆分成长程和短程两个部分：和其中erfc(x)为误差补函数，衰减迅速，因此可以使用直接截断方法处理；erf(x)为误差函数，衰减缓慢，因此将通过傅里叶变换转变到傅里叶空间中快速衰减的级数，再对频率进行截断来计算。因此对粒子i得到的整体受力表达式为：

其中第二项为伦纳德琼斯势能对应的受力，结构因子Im表示取虚部，r_ij＝r_j-r_i。在Ewald方法中，选择合适的α以及截断长度k_c和r_s，使得长程和短程相互作用计算复杂度均为O(N^3/2)。

现在的方法通过对最后一项傅里叶空间中的受力应用快速傅里叶变换以及插值等加速手段，达到降低整体计算复杂度的效果，最终可以达到O(NlogN)。对于短程受力而言，目前的方法是对每个粒子建立近邻列表，计算近邻列表中粒子产生的整体受力。此外，还应用了元胞列表法，使得整体计算复杂度下降至线性；以及考虑比实际需要略大的近邻列表，使得降低近邻列表更新频率。

现有技术计算复杂度较高：对于库仑相互作用长程部分，现有的方法复杂度较高；在处理短程受力时，截断半径内相互作用粒子数量仍然较多(几百个)，计算复杂度大。现有的方法在处理库伦相互作用长程部分时，使用多次快速傅里叶变换，在多核并行模拟时消耗大量时间用于信息通讯，降低了并行效率。每个粒子都需要储存较大近邻粒子列表，其占用的CPU内存较高，在模拟非常大的粒子体系时，会受到内存瓶颈限制。

因此，本领域的技术人员致力于开发一种随机分批埃瓦德算法的近场加速方法(IRBE)。计算复杂度为线性，比现有方法更快；只需要在CPU核之间通讯有限个频率样本，不需要密集通讯，并行效率高；构建双层近邻列表，只需要存储内层近邻的信息，减少CPU内存占用和计算量。