我需要一些关于如何找到适当算法用于我的模拟的建议。

Question

我正在进行一项蒙特卡罗模拟来模拟一些粒子。我的代码中有几个瓶颈，但主要的问题在于，在某些尝试中，我需要更新全部粒子的属性。代码是用c++编写的，目前我有几个循环来实现：
1. 一个循环来存储所有粒子的旧属性并更新新属性。
2. 一个二维循环的相互作用。
3. 另一个二维循环的相互作用（我无法将其与第一个合并）。
4. 一个循环来存储接受步骤/拒绝步骤。

我希望通过交换来删除步骤4，但我找不到方法。所有粒子都是一个类，其中包含名为properties和nextProperties或oldProperties的几个元素。您会如何处理这个问题？

- Yotam

---

没有看到代码的情况下回答这个问题很困难。你能否提供一个简单的例子来展示你遇到的问题？ - Björn Pollex

你可以尝试计时特定的部分，看看瓶颈在哪里。一旦找到了瓶颈，发布慢速代码将有助于SO社区的好心人帮助你找到解决方案。 - Carl Winder

Answer 1

听起来你可能需要使用双缓冲。基本上，您将维护对两个粒子对象数组的指针 - 称它们为accepted和trial。在试验开始时，将accepted数组中粒子的属性复制到trial数组中，并进行任何所需的修改。如果试验成功，则只需交换指针，以使原来的trial数组成为accepted，反之亦然。

此外，您说只有一些试验涉及昂贵的更新。如果是这样，您可能会对快速变量拖动或集合更新等技术感兴趣。

- Ilmari Karonen

---

我不知道你的C++技能如何，你知道我是否可以持有一个迭代器指针吗？ - Yotam

如果你有 vector<Foo> v1 和 vector<Foo> v2，那么 v1.swap(v2) 应该是很便宜的。 - clstrfsck

Answer 2

步骤#2和#3将是O(N^2)，而#1和#4是O(N)，所以你应该专注于这些。

如果您绝对需要在每个粒子之间进行计算，那么您无能为力。但很可能，相距一定距离的粒子对彼此几乎没有影响，或者您只需要处理最近的k个粒子（对于固定的k）。在这种情况下，八叉树（或kd树，aabb树或类似物）是减少成对计算数量的最佳选择。

特别地，您可能希望了解Barnes-Hut方法，该方法用于减少重力计算的复杂性。

- celion

---

第三步是静电，它在许多方面类似于重力。在第四步中，我正在进行大量的复制，这可能会耗费很多时间，但我会看看你的方法。 - Yotam

Answer 3

我认为主要瓶颈在于你的二维交互。如果你正在进行N^2交互，你需要考虑是否可以使用平均场近似进行工作，基本上将计算域分割成单元格，为每个迭代的每个单元格计算一个平均场，然后仅计算与同一单元格中的粒子的相互作用，加上周围单元格的平均场修正。 在这里你可以阅读更多关于此优化以及何时适用的详细信息。

- lurscher

---

前N^2步是逐一处理单元格。第二个步骤不是（我非常确定这是真正的瓶颈）。但是，我需要检查这对我的目的是否可接受。 - Yotam

这取决于你进行第二次N^2遍历的原因以及该遍历所执行的操作，然后我们可以帮助你。您能详细解释一下这一步骤吗？ - lurscher

第三步是对系统中所有带电粒子求和。当有1000个这样的粒子时，代码运行非常缓慢。我还使用了一种称为Lekner求和的方法来处理系统的周期性属性。这种方法迫使我使用一个映射。 - Yotam

你为什么不使用单元格来计算第三步呢？如果你想要一个错误估计器，你可以增加每个粒子周围的最小单元格阈值，使均场近似开始生效，并使用差异来进行估计。 - lurscher

我不使用单元格方法的主要原因是我需要问我的顾问是否可以这样做。然而，这只能在星期天发生。 - Yotam