C++分子动力学模拟程序

C++ Molecular Dynamics Simulation Code

本文关键字:模拟程序 动力学 C++      更新时间:2023-10-16

我正在进行一个模拟气体硬球模型的项目。(类似于理想气体模型。)

我已经写了我的整个项目,它正在发挥作用。为了让你了解我所做的事情,有一个循环可以执行以下操作:(伪代码)

Get_Next_Collision(); // Figure out when the next collision will occur
Step_Time_Forwards(); // Step to time of collision
Process_Collision(); // Process collision between 2 particles
(Repeat)

对于大量粒子(比如N个粒子),必须进行O(N*N)检查,以确定下一次碰撞何时发生。遵循上述程序显然效率低下,因为在绝大多数情况下,粒子对之间的碰撞不受其他地方碰撞处理的影响。因此,希望具有某种形式的优先级队列,该队列存储每个粒子的下一个事件。(事实上,由于碰撞涉及2个粒子,因此只会存储一半数量的事件,因为如果a与B碰撞,那么B也会与a碰撞,并且同时碰撞。)

我发现很难编写这样一个事件/冲突优先级队列。

我想知道是否有任何已经编写的分子动力学模拟器,我可以去查看源代码,以了解如何实现这样的优先级队列。

在谷歌上搜索后,我很清楚有很多MD程序已经编写完成,但其中许多程序要么过于复杂,要么不合适。

这可能是因为它们具有巨大的功能,包括产生可视化的能力或计算粒子之间相互作用力的模拟的能力,等等。

有些模拟器不适合,因为它们为不同的模型进行计算,即:与具有弹性碰撞的节能硬球模型不同的模型。例如,与电势相互作用的粒子或非球形粒子。

我试着查看LAMMPS的源代码,但它非常庞大,我很难理解它

我希望这是足够的信息,我正在努力做什么。如果没有,我可能会添加更多的信息。

位置感知系统的基本版本可能如下所示:

  1. 将宇宙划分为一个立方体网格(每个立方体有a面,体积为a^3),其中每个立方体足够大,但比系统的总体积小得多。每个网格立方体被进一步划分为4个子立方体,理论上可以将其粒子提供给相邻的立方体(并用于计算)
  2. 每个网格立方体都会注册包含在其中的粒子,并知道其相邻网格立方体包含的粒子
  3. 将粒子的可观测宇宙定义为半径为(网格尺寸/2)。定义timestep=(griddim/2) / max_speed。这假设来自最多四个相邻网格立方体的粒子理论上可以在该时间段内相互作用
  4. 对于每个网格立方体中的每个粒子,运行传统的碰撞检测算法(使用mini_timestep < timestep,检查每个粒子是否与其可观测宇宙中的其他粒子发生碰撞。将碰撞存储到按时间排序的任何结构中,甚至只是一个按碰撞时间排序的数组中
  5. 在mini_timestep中发生的第一次碰撞将您的宇宙(和宇宙时钟)重置为(last_time+time_to_collide),其中time_to_collide<mini_timestep。我想这与你目前的算法没有什么不同。重要提示:粒子的绝对坐标会更新,但它们所属的网格立方体和子立方体不会更新
  6. 重复步骤5,直到大timestep通过。按每个栅格正方形更新粒子的所有权

该系统的优点是,对于每个时间窗口,我们都有(假设粒子均匀分布)O(universe_particles*grid_size)而不是O(universe_particles*universe_size)碰撞检查。在良好的条件下(取决于宇宙大小、速度和粒子密度),你可以将计算效率提高几个数量级。

我不明白"优先级队列"方法是如何工作的,但我有一种替代方法可以帮助您。这就是我认为@Boyko Perfanov"利用当地"的意思。

可以将粒子分类到"桶"中,这样就不必将每个粒子相互对照(O(n²))。这利用了这样一个事实,即粒子只有在彼此非常接近的情况下才能碰撞。创建表示小面积/体积的桶,并填充当前位于桶的面积/体积中的所有粒子(O(n)最坏情况)。然后将桶内的所有颗粒与桶内的其他颗粒进行比较(O(m*(n/m)²)平均情况下,m=桶数)。桶需要重叠才能工作,否则也可以检查相邻桶中的粒子。

更新:如果粒子可以移动比桶大小更长的距离,一个明显的"解决方案"是减少时间步长。然而,这将再次增加算法的运行时间,并且只有在有最大速度的情况下才能工作。

即使没有最高速度,另一种适用的解决方案是创建一个额外的"高速"铲斗。由于速度分布通常是高斯曲线,因此不需要将太多粒子放入桶中,因此"桶方法"仍然比O(n²)更有效。