c - 使用 Open-MP 的多核并行计算。
问题描述
我正在努力弄清楚如何使用 OpenMP 并行化此代码,不胜感激。下面是基本代码和描述。
在模拟一组软粒子(例如流体中的蛋白质)时,当一对粒子重叠时,它们之间会产生排斥力。这项任务的目标是使用并行计算来加速这些排斥力的计算,使用带有 Open-MP 的多个内核。
在力斥力函数中,假设粒子具有单位半径。粒子在尺寸 L × L × L 的“模拟盒”中。尺寸 L 的选择使得粒子的体积分数为 φ = 0.3。模拟框具有周期性(环绕)边界条件,这解释了为什么我们需要使用余数函数来计算两个粒子之间的距离。如果粒子重叠,即两个粒子之间的距离 s 小于 2,则斥力与 k(2-s) 成正比,其中 k 是力常数。力沿着连接两个粒子的向量。
- 编写一个程序来测试代码的正确性。这可以通过计算正确的力并将它们与优化代码计算的力进行比较来完成。在您的报告中提供证据证明您的程序使用您的测试程序可以正常工作
- 与提供的基线代码相比,您的加速代码快了多少?包括不同问题规模的时间安排。请务必在报告中包含您的代码列表。
并行化的代码
#include <stdlib.h>
#include <stdio.h>
#include <math.h>
#include <sys/time.h>
double get_walltime() {
struct timeval tp;
gettimeofday(&tp, NULL);
return (double) (tp.tv_sec + tp.tv_usec*1e-6); }
void force_repulsion(int np, const double *pos, double L, double krepulsion, double *forces)
{
int i, j;
double posi [4]; double rvec [4];
double s2, s, f;
// initialize forces to zero
for (i=0; i<3*np; i++)
forces[i] = 0.;
// loop over all pairs
for (i=0; i<np; i++)
{
posi[0] = pos[3*i ];
posi[1] = pos[3*i+1]; posi[2] = pos[3*i+2];
for (j=i+1; j<np; j++)
{
// compute minimum image difference
rvec[0] = remainder(posi[0] - pos[3*j ], L);
rvec[1] = remainder(posi[1] - pos[3*j+1], L);
rvec[2] = remainder(posi[2] - pos[3*j+2], L);
s2 = rvec [0]* rvec [0] + rvec [1]* rvec [1] + rvec [2]* rvec [2];
if (s2 < 4)
{
s = sqrt(s2);
rvec[0] /= s; rvec[1] /= s;
rvec[2] /= s;
f = krepulsion*(2.-s);
forces[3*i ] += f*rvec[0];
forces[3*i+1] += f*rvec[1];
forces[3*i+2] += f*rvec[2];
forces[3*j ] += -f*rvec[0];
forces[3*j+1] += -f*rvec[1];
forces[3*j+2] += -f*rvec[2]; }
} }
}
int main(int argc, char *argv[]) {
int i;
int np = 100; // default number of particles
double phi = 0.3; // volume fraction
double krepulsion = 125.; // force constant
double *pos; double *forces;
double L, time0 , time1;
if (argc > 1)
np = atoi(argv[1]);
L = pow(4./3.*3.1415926536*np/phi, 1./3.);
// generate random particle positions inside simulation box
forces = (double *) malloc(3*np*sizeof(double));
pos = (double *) malloc(3*np*sizeof(double));
for (i=0; i<3*np; i++)
pos[i] = rand()/(double)RAND_MAX*L;
// measure execution time of this function
time0 = get_walltime ();
force_repulsion(np, pos, L, krepulsion, forces);
time1 = get_walltime ();
printf("number of particles: %d\n", np);
printf("elapsed time: %f\n", time1-time0);
free(forces);
free(pos);
return 0; }
解决方案
从理论上讲,它会像这样简单:
void force_repulsion(int np, const double *pos, double L, double krepulsion,
double *forces)
{
// initialize forces to zero
#pragma omp parallel for
for (int i = 0; i < 3 * np; i++)
forces[i] = 0.;
// loop over all pairs
#pragma omp parallel for
for (int i = 0; i < np; i++)
{
double posi[4];
double rvec[4];
double s2, s, f;
posi[0] = pos[3 * i];
//...
汇编:
g++ -fopenmp example.cc -o example
请注意,我没有检查正确性。确保您不会在并行内部有全局变量(因为我更新了您的代码..)
推荐阅读
- python - 如何在字符串中查找一个字符并将其替换为python中的以下字符
- react-native - React Native Redux - 限制 redux 状态更改时导致的组件重新渲染(针对数组子级)
- c# - 为什么 Math.Round() 没有像我在 C# 中所期望的那样产生结果
- signature - 如何解决 Zeek 签名中的规则定义两次错误?
- go - Go 中的依赖注入
- python - 在 Python 中使用不同长度的列表中的元素执行数学运算
- javascript - _lodash 按 itemId、颜色和大小对对象中的数组进行排序和分组并获得总和
- javascript - Javascript / Vue3 - Mixins - Return 'null' by default
- lua - 后备表是 lua 中的常见做法吗?
- python - 如何在 Pygame 中连续平滑地移动图像?