首先，在使用正确的集成方法之前，甚至不要考虑“优化”。处理更多的计算机可能听起来像是解决问题的最简单方法，但实际上它比您想象的要困难得多（由于稍后解释的隐式 BLAS 多线程）。如果你使用显式欧拉，你可以通过改变方法做得更好。适当的自适应时间步进（可能是隐含的，取决于刚度）算法可以比显式 Euler 更有效几个数量级。甚至不要考虑并行化和保持显式欧拉，除非你有很好的理由（即双曲 PDE，所以你担心 SSP，但在这种情况下，只需使用 SSP 优化的积分器......你不会去有充分的理由）。

接下来，您的导数计算函数是迄今为止计算中最昂贵的部分。首先，您应该确保它已优化。语言在这里并不重要，因为所有时间都将在 BLAS 内核中进行矩阵乘法，但请确保它不是通过就地操作进行分配，尽可能使用 BLAS 或 SIMD 循环等。

一旦一切顺利，请考虑是否仍需要并行化。使用正确的 ODE 方法和性能良好的导数实现，10^5 可以相当快。如果你最终使用了隐式方法，你会想要做一些事情，比如将线性求解器换成一些迭代求解器，这会给你另一个很好的加速。

但由于这是一篇关于并行性的文章，我会继续下去，假设你想做一些极端的事情。

如何并行化：要知道的细节

由于和是矩阵，如果您使用具有标准 BLAS 设置的语言，那么您已经通过多线程进行并行化。只需检查您的 CPU 使用情况（例如：在 Linux 上），您就会看到在 Python、R、Julia、MATLAB 等中使用多线程进行足够大的矩阵乘法运算。因此，获得更多性能并不像将其并行化那么简单，因为如果您使用的是单台计算机，则很可能您已经以最佳方式使用了所有内核。$A$$D$htop

像 parareal 这样的并行时间算法只有在 ODE 便宜时才有意义，因为没有其他可用的并行性。时间并行算法效率不高，但它们试图通过允许并行性来弥补这一点。这些算法的开销截止可能相当高，我读过 32 个内核是一个很好的低估计。这些绝对不适合您描述问题的方式。

但是，还有其他两个级别的并行性。如果您的求解成本很高，那么同时在多个系统上并行求解将几乎是最佳加速。即在 16 核机器中一次解决 16 个系统的并行 for 循环或映射，由于与完全集成相比，数据传输时间几乎为 0，因此速度将提高约 16 倍。例如，如果您使用的是 Julia，您可以很容易地自己编写它，或者使用DifferentialEquations.jl 中内置的内容，或为任何其他语言找到合适的资源。请记住，BLAS 已经是多线程的（并且这些矩阵乘法很可能是循环中最昂贵的操作），因此只有当您可以将其外包给其他计算机/节点时，这才有意义。Julia 的原生并行性可以做到这一点，Python 的多处理可以做到这一点，您可以在 C/Fortran 中使用 MPI 等。如果您有一个大集群并且希望通过少量工作获得很多收益，这是解决您问题的一种方法。

最后，如果您有一个优化的导数函数，大部分时间都在矩阵运算中，但没有额外的节点挂起（或者如果有，您仍然可以应用它），您可以利用的另一个并行级别是方法内. 您已经通过多线程 BLAS 操作进行级别内并行，但您可以做得更好。GPU 本质上是并行线性代数机，所以如果你有一个闲置的，那么你可以只执行矩阵乘法，那么这是一个好主意。设置它并不难，因为像 CuBLAS 这样的东西只是 GPU 上的 BLAS，所以你只需要调用它。但是，我认为如果你天真地这样做，你会遇到记忆问题，所以这可能很困难。

TL;博士

停止使用显式欧拉：您应该做的第一件事是真正优化您的代码和方法。如果这还不够，请同时解决不同节点上的系统。如果这还不够，请使用 GPU/Xeon Phis/等加速线性代数。因为限制步骤是矩阵运算。

作为参考，这里是解释优化步骤的笔记本，这里是解释方法内 GPU 并行性的博客文章。