C++ 中的多线程分子模拟

Question

我正在用 C++ 开发一个分子动力学模拟代码，它基本上将原子位置和其他属性作为输入，并根据牛顿运动定律模拟它们的运动。核心算法使用所谓的 Velocity Verlet 方案，如下所示：

//  iterate through time (k=[1,#steps])
double Dt = 0.002; // time step
double Ttot = 1.0; // total time
double halfDt = Dt/2.0;

for (int k = 1; k*Dt <= Ttot; k++){
    for (int i = 0; i < number_particles; i++)
        vHalf[i] = p[i].velocity + F[i]*halfDt; // step 1

    for (int i = 0; i < number_particles; i++)
        p[i].position += vHalf[i]*Dt; // step 2

    for (int i = 0; i < number_particles; i++)
        F[i] = Force(p,i); // recalculate force on all particle i's

    for (int i = 0; i < number_particles; i++)
        p[i].velocity = vHalf[i] + F[i]*halfDt; // step 3
}

p 是一个类对象数组，用于存储粒子位置、速度、质量等信息，Force 是一个函数，它使用 Lennard-Jones 势等计算粒子上的净力。

我的问题是关于完成计算所需的时间；我所有的子程序都在处理数字方面进行了优化（例如，使用x*x*x 提高到三次方而不是pow(x,3)），但主要问题是时间循环通常会执行数百万次迭代，并且通常有接近一百万个粒子。有没有办法使用多线程来实现这个算法？根据我的理解，多线程本质上打开了另一个进出 CPU 内核的数据流，这将允许我同时运行两个不同的模拟；我想使用多线程来使其中一个模拟运行得更快

Answer 1

我建议使用OpenMP。

您的特定用例可以轻松并行化。 并行化应该很简单：

double Dt = 0.002; // time step
double Ttot = 1.0; // total time
double halfDt = Dt/2.0;

for (int k = 1; k*Dt <= Ttot; k++){

    #pragma omp parallel for
    for (int i = 0; i < number_particles; i++)
        vHalf[i] = p[i].velocity + F[i]*halfDt; // step 1
        p[i].position += vHalf[i]*Dt; // step 2

    #pragma omp parallel for
    for (int i = 0; i < number_particles; i++)
        F[i] = Force(p,i); // recalculate force on all particle i's    
        p[i].velocity = vHalf[i] + F[i]*halfDt; // step 3
}

大多数流行的编译器和平台都支持 OpenMP。