CN112733416B - 一种计算粒子间相互作用力的方法及系统 - Google Patents

Abstract

本发明公开了一种计算粒子间相互作用力的方法和系统，涉及分子动力学模拟系统技术领域，包括基于从傅里叶空间中随机抽取的部分傅里叶项、所述傅里叶空间的电荷分布、粒子的带电量，得到所述粒子的库仑长程力；其中，所述傅里叶空间与计算所述库仑长程力的函数相对应。采样方法包括重要性采样，随机抽取的部分傅里叶项的数量级为百或千。该方法还包括根据库仑定律得到库仑短程力；根据其他参数得到粒子的伦纳德琼斯力和成键相互作用力；最后得到相互作用力。该方法可用于分子动力学模拟中，减少计算的时间复杂度，减小数据通信量，从而提高计算性能和并行效率，降低通信方面的消耗。

Description

一种计算粒子间相互作用力的方法及系统

技术领域

本发明涉及分子动力学模拟技术领域，尤其涉及一种计算粒子间相互作用力的方法及系统。

背景技术

分子动力学模拟是描述微纳尺度体系演化的一种重要的计算机模拟方法，被广泛应用于生物分子、药物设计和能源材料等领域，是科学计算最为重要的方法之一。它通过计算粒子之间的相互作用力和牛顿运动方程得到所有粒子的速度和位置，从这些轨迹得到一系列的热力学和动力学性质，如密度分布、结合能、扩散系数、迁移率等。传统的分子动力学力场将势能函数分解为短程的键能相互作用、范德华相互作用和长程的库仑相互作用，并依次计算它们对势能的贡献。

在分子动力学模拟的力计算中，库仑长程相互作用的处理具有最高的时间复杂度O(N²)，其中N是体系中总的粒子数。在大规模生物体系中，N的取值可能会达到亿(10⁸)的数量级。因此，库仑长程相互作用的处理耗时巨大，往往要占据90％以上的计算资源。另一方面，库仑相互作用又是不可忽略的。不计算库仑力就没办法准确地描述系统中的多体关联，也无法准确模拟大部分软物质系统的功能和性质。例如DNA和RNA等生物大分子的结构和功能的调控，胶体颗粒的自组装，水分子的结构性质以及细胞膜上钠-钾泵离子通道的形态和功能等，都是由库仑相互作用来起主要调控作用的。因此，发展快速计算库仑力的方法是十分必要的，也是分子模拟领域近半个世纪以来最亟待解决的问题之一。

前人的研究中，提出了许多在当时的时代技术背景下有效的库仑相互作用处理方法。一类重要的方法是基于球谐函数或其他基函数对自由空间中的三维泊松方程的格林函数做展开的思想而提出的。这类方法的鼻祖是Greengard和Rokhlin于1987年提出的快速多级子(FMM，Fast Multipole Method)方法。该方法通过构造拉普拉斯方程基本解的多级展开和局部展开来快速计算库仑力，并在随后的三十年间推广至任意满足二阶椭圆偏微分方程的核函数及部分径向基函数，其优点是具有O(N)的最优线性复杂度。Ewald求和方法是另一类计算库仑力的方法，包括经典Ewald方法，其复杂度为O(N^3/2)。及后续使用格矢求和及快速傅里叶变换(FFT，Fast Fourier Transform)改进的Ewald方法，即PPPM埃瓦尔德方法(particle-particle particle mesh Ewald，粒子-粒子-粒子网格埃瓦尔德方法)，其复杂度为O(NlogN))。无论哪种方法，它们都是在傅里叶空间中利用全部傅里叶项来进行计算，其计算时间复杂度都是很高的。

此外，现有技术的计算数据通信量过高也成为了瓶颈。近年来随着硬件性能的提升，特别是算力更强的高性能服务器的出现，并行效率低也成为了现有技术的一个普遍存在且亟待解决的问题。当使用核数增加时，现有方法中数据通信时间的占比甚至会大大超过计算时间的占比。基于FMM的方法需要构造复杂的自适应分层树结构。例如，二维问题要构建四叉树，三维问题要构造八叉树。在实现这些方法时需要进行复杂的进程间通信处理，使得这些方法的并行效率下降严重。FMM方法的另一个缺点是构建复杂数据结构所需的处理时间过长。经典的Ewald求和方法每次计算需要传输的数据量是O(N^3/2)，这表示随着粒子数的增加，通信时间会快速增长，从而导致这些方法的并行效率降低。而作为基于FFT的方法，PPPM方法更是需要CPU与CPU之间广泛而密集的通信，并且复杂度为O(NlogN)，无法达到最优的线性复杂度。

因此，本领域的技术人员致力于开发一种用于分子动力学模拟中计算粒子间相互作用力的方法及系统，来减少计算的时间复杂度，减小计算数据通信量，从而提高计算性能和并行效率。

发明内容

有鉴于现有技术的上述缺陷，本发明所要解决的技术问题是如何从傅里叶空间中抽取部分傅里叶项来得到库仑长程力，一方面减少计算粒子库仑力所需的傅里叶项的数量，使所需的傅里叶项的数量变为一个与分子动力学模拟系统中的粒子数量无关的常数，保证计算高效，另一方面还要保证计算结果准确。

为实现上述目的，本发明提供了一种计算粒子间相互作用力的方法，包括：基于从傅里叶空间中随机抽取的部分傅里叶项、所述傅里叶空间的电荷分布、粒子的带电量，得到所述粒子的库仑长程力；其中，所述傅里叶空间与计算所述库仑长程力的函数相对应。

进一步地，通过重要性采样的方法从所述傅里叶空间中随机抽取所述部分傅里叶项。

进一步地，通过直接采样的方法从所述傅里叶空间中随机抽取所述部分傅里叶项。

进一步地，通过接收-拒绝采样的方法从所述傅里叶空间中随机抽取所述部分傅里叶项。

进一步地，随机抽取的所述部分傅里叶项的数量级为百或千。

进一步地，还包括：根据所述粒子的位置坐标，得到所述粒子的相对距离；基于所述相对距离，得到所述粒子的库仑短程力。

进一步地，还包括：结合所述库仑长程力与所述库仑短程力，得到所述粒子的库仑力。

进一步地，还包括：基于所述粒子的所述位置坐标、范德华半径以及范德华吸引势的特征能量，得到所述粒子的伦纳德琼斯力。

进一步地，还包括：基于所述粒子的键长的平衡位置距离、键角的平衡位置角度，得到所述粒子的成键相互作用力。

进一步地，还包括：结合所述粒子的所述库仑力、所述伦纳德琼斯力和所述成键相互作用力，得到所述粒子间相互作用力。

本发明还提供一种计算粒子间相互作用力的系统，包括：

库仑长程力模块，被配置为能够基于从傅里叶空间中随机抽取的部分傅里叶项、所述傅里叶空间的电荷分布、粒子的带电量，得到所述粒子的库仑长程力；其中，所述傅里叶空间与计算所述库仑长程力的函数相对应。

进一步地，还包括：采样模块，所述采样模块被配置为能够通过重要性采样的方法从所述傅里叶空间中随机抽取所述部分傅里叶项。

进一步地，还包括：采样模块，所述采样模块被配置为能够通过直接采样的方法从所述傅里叶空间中随机抽取所述部分傅里叶项。

进一步地，还包括：采样模块，所述采样模块被配置为能够通过接收-拒绝采样的方法从所述傅里叶空间中随机抽取所述部分傅里叶项。

进一步地，随机抽取的所述部分傅里叶项的数量级为百或千。

进一步地，还包括库仑短程力模块，所述库仑短程力模块被配置为根据所述粒子的位置坐标，得到所述粒子的相对距离；以及基于所述相对距离，得到所述粒子的库仑短程力。

进一步地，还包括伦纳德琼斯力模块，所述伦纳德琼斯力模块被配置为基于所述粒子的所述位置坐标、范德华半径以及范德华吸引势的特征能量，得到所述粒子的伦纳德琼斯力。

进一步地，还包括成键相互作用力模块，所述成键相互作用力模块被配置为基于所述粒子的键长的平衡位置距离、键角的平衡位置角度，得到所述粒子的成键相互作用力。

本发明还提供一种计算粒子间相互作用力的装置，包括存储器、处理器以及存储在所述存储器中并能够在所述处理器上运行的计算机程序，其特征在于，所述处理器被配置为能够在执行所述计算机程序时实现计算粒子间相互作用力的方法的步骤。

本发明还提供一种计算机可读存储介质，其上存储有计算机程序，所述计算机程序被处理器执行时能够实现计算粒子间相互作用力的方法的步骤。

与现有技术方案相比，本发明的有益技术效果在于：

(1)使得库仑长程力计算误差的期望为零，在提升效率的同时还保持了动力学的性质。

(2)每次需要计算的傅里叶项的数量变为一个常数，与空间中的粒子数无关，减小了数据通信量，从而大大降低了通信方面的消耗，使得算法具有极高的并行效率，并保持O(N)的最优线性复杂度。

(3)降低了全原子分子动力学模拟计算过程中库仑力计算的时间开销。

以下将结合附图对本发明的构思、具体结构及产生的技术效果作进一步说明，以充分地了解本发明的目的、特征和效果。

附图说明

图1是本发明的一个较佳实施例的流程图；

图2是本发明的一个较佳实施例与PPPM方法的均方位移对比图；

图3是本发明的一个较佳实施例与PPPM方法的径向分布函数对比图；

图4是本发明的一个较佳实施例与PPPM方法的单步计算时间对比图；

图5是本发明的一个较佳实施例与PPPM方法的并行效率对比图。

具体实施方式

以下参考说明书附图介绍本发明的多个优选实施例，使其技术内容更加清楚和便于理解。本发明可以通过许多不同形式的实施例来得以体现，本发明的保护范围并非仅限于文中提到的实施例。

在经典的分子动力学模拟中，库仑力是关于距离的倒数的函数，这种距离倒数

的函数计算是最耗时的部分。一个常用的方法是分解成两部分

和

其中，前一部分

是误差补函数除以距离r，在实空间中关于r呈指数衰减，即衰减率和O(exp(-r))是同一数量级，因此这一部分衰减很快，可以直接在实空间中做截断来实现快速计算，求得短程库仑力。而后一部分

是误差函数除以距离r，在实空间中的衰减速度是O(1/r)，因此衰减很慢，导致两个相距很远的原子之间也会产生库仑相互作用。由此可以看出，主要的工作量在计算

从而求得库仑长程力。

要想减小求得库仑长程力的误差函数除以距离r这一部分的

的计算工作量，通常的想法是将其从实空间变换到傅里叶空间中，进而通过在傅里叶空间中做截断来快速计算。经过傅里叶变换后，第j个粒子的库仑长程力在傅里叶空间中的求解表达式为

其中结构因子

式中，Im表示取虚部，V表示模拟系统的体积，q_j为粒子的电量，r_j表示粒子间的距离。其中k是三维向量，需要取遍

个k才能算出一个粒子相对准确的力，计算量非常庞大，对N个粒子直接计算的最优时间复杂度达到了

为了解决这一问题，就要想办法减少计算所需的三维向量k，也就是减少计算所需的傅里叶项。如前所述，库伦力的实空间计算拆分出来误差函数除以距离和误差补函数除以距离。其中的误差函数除以距离部分需要变换到傅里叶空间计算，因此就要在与计算粒子的库仑长程力对应的傅里叶空间中进行采样，抽取部分傅里叶项进行计算，即对库仑力的长程部分(即误差函数除以距离)做傅里叶变换后做随机采样。一方面，为了能够减少计算量，采样的数量要取得适宜，例如，采样的样本数量级可以从百到千，那么这样的计算量相比于现有技术动辄上亿的计算量大大减小。另一方面，还要保证计算的准确度，一般来说，采样的样本数量越多，相应的，计算结果准确度会越高。根据不同的计算精度要求，可以在傅里叶空间中做相应的截断，取不同数量级的样本。所以，要兼顾计算样本数量的大小和计算结果准确度。

为了确保在使用较少样本，较少计算量的前提下，获得准确的结果，那么从傅里叶空间中随机抽取部分傅里叶项而获得计算样本的采样方法就显得尤为重要，它要能保证计算的误差期望值尽量小，理想情况是为零，以确保计算结果准确。对此，发明人观察到在求解库仑长程力的表达式中有一项

可以被看作为高斯函数，其三维求和的结果收敛到一个常数。因此，可以将其看作一个三维概率分布。即对每个三维向量k对应的计算式

以

的概率对其进行重要性采样，抽样P个k，{k₁，k₂…k_p},每个k_i被抽样的概率不同，重要性越高的k_i，其采样概率就越大。而每个k_i的重要性是依据其对整体系统的贡献来决定的，也就是，根据第i个k值所计算出的高斯函数的值在所有k值计算的高斯函数之和中的占比，占比越高，对应的k_i贡献就越大。

接下来，利用采样得到的P个采样值来估计总体的结果。假设总的采样数为P，P的数量级一般为几百，也可以采样到上千，依据实际需要而定。采样数量越少，则相应的计算量越小。总的重要性为

即可得到第j个粒子的库仑长程力的一个近似值表达式为

由于式(1)和(2)的期望值不变，因此可以认为经过重要性采样后，将P个三维向量k求得的傅里叶空间库仑长程力的值代入求解牛顿方程，经过一定时间模拟后得到的平衡态分布和根据全部k值求得的傅里叶空间库仑长程力计算出的平衡态分布是相同的。

假设体系中有一组数量为N的源电荷

和一个目标粒子B，坐标分别为

如果每个源电荷A_l的带电量为q_l,则这组源电荷在带电量为q_B的粒子B处产生的电势为

其中

这部分为误差补函数除以距离，衰减很快，可以通过高效的截断方法直接计算，即不用算全部N个源电荷，只需要其中很少一部分的贡献就可以计算出库仑短程力。而

对应的静电力，即库仑长程力，需要变换到傅里叶空间中计算，其表达式为

其中

直接计算的最优时间复杂度为

在本发明的一个实施例中，每次抽取100个样本，则库仑长程力的快速计算表达式为

总的计算量为O(100)，是一个常数，大大减少了总的计算量。这意味着，时间复杂度变成了最优的常数O(P)，减小了计算时间复杂度。而k的数量决定了计算的数据通信量。在本实施例中，无需根据全部傅里叶频率来计算库仑长程力，而仅需要利用100个三维向量k的三维求和来计算库仑长程力，大大减小了数据通信量。

在本发明的其他实施例中，还可以应用其他采样方法，如直接采样，接收-拒绝采样等，在计算的傅里叶空间中抽取P个k值，根据前述公式

计算目标点B的库仑长程力。这样的计算方式因为仅需要P个样本来进行计算，不需要使用傅里叶空间中的全部傅里叶项来参与计算，因此计算工作量相比于现有技术也大大减小，减少了计算的时间复杂度。根据对计算准确度的不同要求，在使用直接采样，接收-拒绝采样这些采样方式时，可以使采样样本的数量级取百、千、万等。从保证计算的准确度角度来说，效果最好的仍然是采用重要性采样方法来抽取部分傅里叶项，如前所述，因为该采样方法是对重要性越大的项赋予越高的概率，因而能保证计算结果误差的期望值为零，所以计算结果准确，而且计算效率高。

在当前的主流分子动力学模拟软件(例如Gromacs和Lammps)的数据结构中，如果总的计算资源包含L个CPU线程，粒子的信息会被分为L份分别存储在每个线程的局部存储空间中，每个线程约存储有N/L个粒子的信息。如果在力的计算过程中，某一部分的计算需要全部粒子的信息，那么就必须通过至少一次全局通信来将所有粒子的信息传递到每个线程的局部存储空间上。例如在Ewald方法中，对结构因子

的计算便需要全部N个粒子的信息。因此对总数量为

的全部k，计算ρ(k)都需要一次全局数据传输，总的数据通信量为

而在本发明提出的方法中，通过使用重要性采样的方法，利用抽出的P个k对应的结构因子逼近总粒子的计算结果，使得总的计算数据通信量为一个常数O(P)，大大减小了计算时的数据通信量，从而减少了通讯消耗。

相比于现有技术，该方法还减小了计算的时间复杂度。在经典分子动力学模拟中，特别是对生物系统中的大规模体系的模拟，总的粒子数N可能会达到亿的数量级。因此，这里所说的时间复杂度是指在计算N个源电荷对N个目标点的库仑相互作用时，总的计算量同总粒子数N的数量关系。本发明中提出的方法对

和

分别使用空间截断和重要性采样的方法，使得计算

的时间复杂度为O(CN)，计算

的时间复杂度为O(PN)，其中C是截断区域里的粒子数，P是傅里叶空间的抽样数量。对于固定的截断半径和傅里叶空间抽样数，C和P是和N无关的两个常数，因此该方法总的时间复杂度是O(N)，达到了线性最优。

在本发明的另一个实施例中，如图1所示，显示了本发明所提计算粒子间相互作用力的方法在分子动力学模拟中具体应用时的流程图。

首先，设置一系列计算参数，其中主要包括四类参数，第一类参数是计算原子库仑力要用到的，包括原子的电荷，库仑短程力的近场截断距离，库仑长程力的高斯函数中的误差控制常数α以及每步的采样数量P。第二类参数是计算原子的伦纳德琼斯力所需要用到的，包括原子的范德华半径以及范德华吸引势的特征能量。第三类是计算原子的成键相互作用力所需要的用到的，包括键长的平衡位置距离，键的力学常数，键角的平衡位置角度及其键角的力学常数。第四类参数为其他参数，包括动力学积分步长、温度、总积分步数、动力学轨迹，力的输出频率。

接着，从模拟系统的起始结构读取初始状态下的原子坐标，并根据设定的温度随机对每个原子分配初始速度。后续可以在每一步积分过程中更新原子的坐标，以便计算原子的相对距离。

接下来，判断是否达到了模拟停止条件。在一些实施例中，模拟停止条件可以为是否已经完成了设定的总积分步数。如果已经完成，则停止。如果未完成，则执行下列操作：

利用重要性采样的方法，对可以看作是高斯分布的

进行采样，采样得到P个k；其中，α为计算库仑长程力的高斯函数中的误差控制常数，P为采样数量，k为傅里叶空间中的傅里叶项，即傅里叶空间频率，是一个三维向量。

利用重要性采样得到的P个k，结合公式

得到每个粒子在傅里叶空间的库仑长程力；其中，总的重要性

结构因子ρ(k_l)是傅里叶空间的电荷分布，

q_j为所述粒子的带电量，V为模拟系统的体积。其中，采样数量P的数量级可以为百或千。

接下来，可以运用库仑定律直接求得粒子的库仑短程力。

将粒子的库仑长程力与库仑短程力相加，可以得到粒子的库仑力。

根据粒子坐标，粒子的范德华半径以及范德华吸引势的特征能量，可以得到粒子伦纳德琼斯力。

根据键长的平衡位置距离，键的力学常数，键角的平衡位置角度及其键角的力学常数，可以得到粒子成键相互作用力。

原子的库仑力、伦纳德琼斯力和成键相互作用力全部加起来，可以得到原子所受到的合力。根据牛顿第二定律，原子的合力除以其质量，可以得到原子运动的加速度，再结合原子当前时刻的速度和位置坐标，做一步牛顿积分，便可以得到下一个时刻粒子的速度和位置坐标。然后再根据新时刻的原子的位置，重新计算新时刻原子的受力，循环反复整个过程，直到跑完了程序设定的总积分步数。

最后根据动力学轨迹、力的输出频率这些程序输出参数来控制程序输出结果。

以上就是本发明所提计算粒子间相互作用力的方法在分子动力学模拟中的完整应用过程。

在本发明的其他实施例中，还可以在上述流程中每次做更新计算时，采用不同的采样方法，例如直接采样、接收-拒绝采样方法抽取P个k值，再结合公式计算，以充分发挥各种采样方法的优点。

在本发明的另一些实施例中，还可以在上述流程中每次做更新计算时，取不同的采样样本数量，根据实际情况，取百或千数量级的样本，而无需固定一个样本数值，以适应不同计算情况的需求。

在实现方面，申请人设计了耦合Ewald和随机分批思想的长程相互作用快速随机方法的高性能并行化程序。主要考虑两方面的性能：一是尽可能的减少运算量，特别是减少在并行模拟中效率低下的全局变量的运算；二是尽可能的减少集群中不同处理器间的通讯时间，如通过传输数据和运算的交替进行缩短处理器等待数据的时间，并使用紧致而高效的数据结构来更有效地利用计算资源。计算程序提供多核并行版本，使用MPI(跨语言的通讯协议)处理在分布式内存体系结构的不同计算节点之间的数据通信，使用OpenMP(共享内存并行系统的多线程程序设计方案)来协调每个计算节点的计算负载，降低代码阅读复杂性的同时最大限度提升计算资源的使用率，同时充分利用高性能集群大量计算单元的优势，让属于同一台机器的多个CPU分别承担运算和处理器间通讯的不同任务。申请人已经将实现本方法的计算机代码成功地整合入了全原子分子动力学模拟软件lammps中，并验证了其在全原子模拟中的正确性和高并行效率。

在本发明的另一个实施例中，将实现本方法的计算机代码成功地整合入了全原子分子动力学模拟软件lammps中，并且已经在包含有1000万水分子的纯水的全原子分子动力学模拟中，与传统PPPM方法进行了准确性和计算效率的对比验证。图中RBE(random batchEwald随机分批埃瓦尔德)代表本发明所提方法。比较了水分子的动力学信息，即均方位移，如图2所示。还比较了结构信息，即氧原子-氧原子之间的径向分布函数，如图3所示。在这个水分子的实例中，本发明所提方法中的采样数量P为100，误差控制常数α为0.01，水分子模型以及相关的势能参数采用SPCE模型，伦纳德琼斯力和库仑相互作用力的短程截止距离设置为1.2纳米，PPPM与本发明所提RBE(random batch Ewald，随机分批埃瓦尔德)这两种方法各自总共模拟了20万步，积分步长为1飞秒。对比发现本方法与传统的PPPM方法得到的结果完全一致。同时本方法的计算效率，即动力学模拟跑一步所需要的时间，在一万个CPU核心上比PPPM方法要低10倍左右，如图4所示。这意味着，本发明的时间复杂度更低。

此外，从图4中可以注意到本方法在CPU核心比较少的时候，与经典的PPPM方法的计算效率相近；而在CPU核心越多的时候，PPPM方法的计算效率达到了瓶颈，其每一步计算的所需时间没有明显改进。PPPM方法的计算效率达到瓶颈实际上是因为，PPPM使用了快速傅里叶变换做加速，使得该方法的通讯时间消耗会随着使用CPU核数的增长而增长。最终，通讯时间的增长率会超过计算时间的下降率，使得计算效率达到瓶颈。而本发明所提方法，由于数据通信量小，所以通讯消耗很低，即使CPU核数很高的时候也不会有明显的上升。因此，计算效率仍可以继续提升。所以，本发明所提方法更加适用于针对包含大量原子的全原子系统，在更大规模的多核CPU的超级计算机上进行分子动力学模拟。

图5是本方法与PPPM方法的并行效率对比图。从图中可以看出，本方法获得了更高的计算并行效率，这是因为本方法减少了CPU核与核之间的交换。

综上所述，本发明基于随机分批思想改进了经典Ewald求和中傅里叶空间的计算方法，将重要性采样方法、直接采样或者接收-拒绝采样等用于傅里叶空间的采样，从高斯分布中抽取少部分样本来进行计算，从而大大减少了计算量，并且保证了计算结果的准确性。同现有技术相比，在避免了使用通讯密集型方法FFT的同时，使得方法的并行效率大大提高。在实现方面，通过MPI(跨语言的通讯协议)和OpenMP(共享内存并行系统的多线程程序设计方案)混合编程实现了本发明所提的方法，降低了全原子分子动力学模拟过程中计算库仑力的时间开销，大大节约了计算所需的时间，使得本发明同现有技术相比更适用于未来在高性能集群上运行的大规模分子动力学模拟程序。在计算机辅助药物设计过程中，例如蛋白质的动力学模拟，蛋白质与小分子的对接，材料的分子动力学模拟等，都可以采用本发明以节约计算时间。

本申请提供的技术方案可以是系统、方法、装置和/或计算机程序产品。计算机程序产品可以包括计算机可读存储介质，其上载有用于使处理器实现本申请的各个方面的计算机可读程序指令。

在一些实施例中，本申请还提供一种计算机装置、设备或终端。该计算机装置、设备或终端包括通过系统总线连接的处理器、存储器、网络接口、显示屏和输入装置。其中，处理器用于提供计算和控制能力，存储器包括非易失性存储介质、内存储器。非易失性存储介质存储有操作系统和计算机程序。内存储器为非易失性存储介质中的操作系统和计算机程序的运行提供环境。网络接口用于与外部的终端通过网络连接通信。计算机程序被处理器执行时以实现本申请公开的各种方法、流程、步骤，或者处理器执行计算机程序时实现本申请公开的实施例中各个模块或单元的功能。显示屏可以是液晶显示屏或者电子墨水显示屏，输入装置可以是显示屏上覆盖的触摸层，也可以是外壳上设置的按键、轨迹球或触控板，还可以是外接的键盘、触控板或鼠标等。

示例性的，计算机程序可以被分割成一个或多个模块或单元，这些模块或单元被存储在存储器中，并可由处理器执行，以实现本申请的技术方案。这些模块或单元可以是能够完成特定功能的一系列计算机程序指令段，该指令段用于描述计算机程序在装置、设备或终端中的执行过程。

上述的装置、设备或终端可以是桌上型计算机、笔记本、移动电子设备、掌上电脑及云端服务器等计算设备。本领域技术人员应当理解，图中所示出的结构，仅仅是与本申请方案相关的部分结构的框图，并不构成对本申请方案所应用于其上的装置、设备或终端的限定，具体的装置、设备或终端可以包括比图中所示更多或更少的部件，或者组合某些部件，或者具有不同的部件布置。

处理器可以是中央处理单元(Central Processing Unit，CPU)，也可以是其他通用或专用的处理器、微处理器、数字信号处理器(Digital Signal Processor，DSP)、专用集成电路(Application Specific Integrated Circuit，ASIC)、可编程门阵列(Field-Programmable Gate Array，FPGA)或者其他可编程逻辑器件、分立门或者晶体管逻辑器件、分立硬件组件等。处理器是上述的装置、设备或终端的控制中心，利用各种接口和线路连接装置、设备或终端的各个部分。

存储器可用于存储计算机程序、模块和数据，处理器通过运行或执行存储在存储器内的计算机程序和/或模块，以及调用存储在存储器内的数据，实现装置、设备或终端的各种功能。存储器可主要包括程序存储区和数据存储区，其中，程序存储区可存储操作系统、至少一个功能所需的应用程序(比如声音播放功能、图像播放功能等)等；数据存储区可存储根据应用所创建的各类数据(比如多媒体数据、文档、操作历史记录等)等。此外，存储器可以包括高速随机存取存储器，还可以包括非易失性存储器，例如硬盘、内存、插接式硬盘、智能存储卡(Smart Media Card，SMC)、安全数字(Secure Digital，SD)卡、闪存卡(Flash Card)、磁盘存储器件、闪存器件、或其他易失性固态存储器件。

本申请还提供一种计算机可读存储介质，其上存储有计算机程序，计算机程序被处理器执行时实现上述方法的步骤。本领域普通技术人员可以理解实现上述实施例方法中的全部或部分流程，是可以通过计算机程序来指令相关的硬件来完成，所述的计算机程序可存储于一非易失性计算机可读取存储介质中，该计算机程序在执行时，可包括如上述各方法的实施例的流程。其中，本申请所提供的各实施例中所使用的对存储器、存储、数据库或其它介质的任何引用，均可包括非易失性和/或易失性存储器。非易失性存储器可包括只读存储器(ROM)、可编程ROM(PROM)、电可编程ROM(EPROM)、电可擦除可编程ROM(EEPROM)或闪存。易失性存储器可包括随机存取存储器(RAM)或者外部高速缓冲存储器。作为说明而非局限，RAM以多种形式可得，诸如静态RAM(SRAM)、动态RAM(DRAM)、同步DRAM(SDRAM)、双数据率SDRAM(DDRSDRAM)、增强型SDRAM(ESDRAM)、同步链路(Synchlink)DRAM(SLDRAM)、存储器总线(Rambus)直接RAM(RDRAM)、直接存储器总线动态RAM(DRDRAM)、以及存储器总线动态RAM(RDRAM)等。

上述的装置或终端设备集成的模块和单元，如果以软件功能单元的形式实现并作为独立的产品销售或使用时，可以存储在计算机可读存储介质中。基于这样的理解，本申请实现所公开的各种方法中的全部或部分流程，也可以通过计算机程序来指令相关的硬件来完成，计算机程序可存储于计算机可读存储介质中，该计算机程序在被处理器执行时，可实现上述各个方法的步骤。其中，计算机程序包括计算机程序代码，计算机程序代码可以为源代码形式、对象代码形式、可执行文件或某些中间形式等。计算机可读介质可以包括：能够携带计算机程序代码的任何实体或装置、记录介质、U盘、移动硬盘、磁碟、光盘、计算机存储器、只读存储器(ROM，Read-Only Memory)、随机存取存储器(RAM，Random Access Memory)、电载波信号、电信信号以及软件分发介质等。需要说明的是，计算机可读介质包含的内容可以根据司法管辖区内立法和专利实践的要求进行适当的增减。

在一些实施例中，本申请公开的各种方法、流程、模块、装置、设备或系统可以在一个或多个处理装置(例如，数字处理器、模拟处理器、被设计成用于处理信息的数字电路、被设计成用于处理信息的模拟电路、状态机、计算设备、计算机和/或用于以电子方式处理信息的其他机构)中被实现或执行。该一个或多个处理装置可以包括响应于以电子方式存储在电子存储介质上的指令来执行方法的一些或所有操作的一个或多个装置。该一个或多个处理装置可以包括通过硬件、固件和/或软件被配置而专门设计成用于执行方法的一项或多项操作的一个或多个装置。以上所述，仅为本申请较佳的具体实施方式，但本申请的保护范围并不局限于此，任何熟悉本技术领域的技术人员在本申请揭露的技术范围内，根据本申请的技术方案及其发明构思加以等同替换或改变，都应涵盖在本申请的保护范围之内。

本申请的实施方式可以在硬件、固件、软件或其各种组合中进行，还可以作为存储在机器可读介质上的且可以使用一个或多个处理装置读取和执行的指令来实现。在一些实施方式中，机器可读介质可以包括用于存储和/或传输呈机器(例如，计算装置)可读形式的信息的各种机构。例如，机器可读存储介质可以包括只读存储器、随机存取存储器、磁盘存储介质、光存储介质、快闪存储器装置以及用于存储信息的其他介质，并且机器可读传输介质可以包括多种形式的传播信号(包括载波、红外信号、数字信号)以及用于传输信息的其他介质。虽然在执行某些动作的特定示例性方面和实施方式的角度可以在以上公开内容中描述固件、软件、例程或指令，但将明显的是，这类描述仅出于方便目的并且这类动作实际上由机器设备、计算装置、处理装置、处理器、控制器、或执行固件、软件、例程或指令的其他装置或机器产生。

在本申请的权利要求书和说明书中，用来执行指定功能的模块或者使用功能性特征描述的模块，意在涵盖能够执行该功能的任何方式，例如：执行该功能的电路元件的组合，用来执行或实现该功能的软件、硬件以及软件和硬件的组合，或者任何形式的软件、固件、代码及其与适当电路或其他装置的组合。由各种模块提供的功能被以权利要求书所主张的方式组合在一起，由此应当认为，是可以提供这些功能的任何模块、部件、元件都等价或等效于权利要求书中限定的模块。根据电路的等效变换的原理，也可以对本申请中一些实施例的电路结构进行变更、修改，例如：将电流源变换为电压源、串联结构变换为并联结构等，从而获得更多样化的实施例，但这些变更和修改均属于本申请公开的范围。

以上详细描述了本发明的较佳具体实施例。应当理解，本领域的普通技术无需创造性劳动就可以根据本发明的构思作出诸多修改和变化。因此，凡本技术领域中技术人员依本发明的构思在现有技术的基础上通过逻辑分析、推理或者有限的实验可以得到的技术方案，皆应在由权利要求书所确定的保护范围内。

Claims (14)

1.一种计算粒子间相互作用力的方法，其特征在于，包括：

通过重要性采样的方法从傅里叶空间中随机抽取部分傅里叶项，基于从所述傅里叶空间中随机抽取的部分傅里叶项、所述傅里叶空间的电荷分布、粒子的带电量，得到所述粒子的库仑长程力；其中，所述傅里叶空间与计算所述库仑长程力的函数相对应。

2.如权利要求1所述的计算粒子间相互作用力的方法，其特征在于，随机抽取的所述部分傅里叶项的数量级为百或千。

3.如权利要求1所述的计算粒子间相互作用力的方法，其特征在于，还包括：

根据所述粒子的位置坐标，得到所述粒子的相对距离；

基于所述相对距离，得到所述粒子的库仑短程力。

4.如权利要求3所述的计算粒子间相互作用力的方法，其特征在于，还包括：

结合所述库仑长程力与所述库仑短程力，得到所述粒子的库仑力。

5.如权利要求4所述的计算粒子间相互作用力的方法，其特征在于，还包括：

基于所述粒子的所述位置坐标、范德华半径以及范德华吸引势的特征能量，得到所述粒子的伦纳德琼斯力。

6.如权利要求5所述的计算粒子间相互作用力的方法，其特征在于，还包括：

基于所述粒子的键长的平衡位置距离、键角的平衡位置角度，得到所述粒子的成键相互作用力。

7.如权利要求6所述的计算粒子间相互作用力的方法，其特征在于，还包括：

结合所述粒子的所述库仑力、所述伦纳德琼斯力和所述成键相互作用力，得到所述粒子间相互作用力。

8.一种计算粒子间相互作用力的系统，其特征在于，包括：

库仑长程力模块，被配置为能够基于从傅里叶空间中随机抽取的部分傅里叶项、所述傅里叶空间的电荷分布、粒子的带电量，得到所述粒子的库仑长程力；其中，所述傅里叶空间与计算所述库仑长程力的函数相对应；

还包括采样模块，所述采样模块被配置为能够通过重要性采样的方法从所述傅里叶空间中随机抽取所述部分傅里叶项。

9.如权利要求8所述的计算粒子间相互作用力的系统，其特征在于，随机抽取的所述部分傅里叶项的数量级为百或千。

10.如权利要求8所述的计算粒子间相互作用力的系统，其特征在于，还包括库仑短程力模块，所述库仑短程力模块被配置为根据所述粒子的位置坐标，得到所述粒子的相对距离；以及基于所述相对距离，得到所述粒子的库仑短程力。

11.如权利要求8所述的计算粒子间相互作用力的系统，其特征在于，还包括伦纳德琼斯力模块，所述伦纳德琼斯力模块被配置为基于所述粒子的位置坐标、范德华半径以及范德华吸引势的特征能量，得到所述粒子的伦纳德琼斯力。

12.如权利要求8所述的计算粒子间相互作用力的系统，其特征在于，还包括成键相互作用力模块，所述成键相互作用力模块被配置为基于所述粒子的键长的平衡位置距离、键角的平衡位置角度，得到所述粒子的成键相互作用力。

13.一种计算粒子间相互作用力的装置，包括存储器、处理器以及存储在所述存储器中并能够在所述处理器上运行的计算机程序，其特征在于，所述处理器被配置为能够在执行所述计算机程序时实现根据权利要求1-7中任一项所述的计算粒子间相互作用力的方法的步骤。

14.一种计算机可读存储介质，其上存储有计算机程序，其特征在于，所述计算机程序被处理器执行时能够实现根据权利要求1-7中任一项所述的计算粒子间相互作用力的方法的步骤。

Images