使用OpenMP并行化两个for循环。

huangapple
go评论93阅读模式
英文:

Parallelise 2 for loops with OpenMP

问题

以下是代码部分的中文翻译：

  1. 所以我有这个函数，我必须使用OpenMP静态调度来并行化，以适应n个线程
  3. void computeAccelerations(){
  4.     int i,j;
  5.     for(i=0;i<bodies;i++){
  6.         accelerations[i].x = 0;	accelerations[i].y = 0; accelerations[i].z = 0;
  7.         for(j=0;j<bodies;j++){
  8.             if(i!=j){
  9.                 //accelerations[i] = addVectors(accelerations[i],scaleVector(GravConstant*masses[j]/pow(mod(subtractVectors(positions[i],positions[j])),3),subtractVectors(positions[j],positions[i])));
  10.                 vector sij = {positions[i].x-positions[j].x,positions[i].y-positions[j].y,positions[i].z-positions[j].z};
  11.                 vector sji = {positions[j].x-positions[i].x,positions[j].y-positions[i].y,positions[j].z-positions[i].z};
  12.                 double mod = sqrt(sij.x*sij.x + sij.y*sij.y + sij.z*sij.z);
  13.                 double mod3 = mod * mod * mod;
  14.                 double s = GravConstant*masses[j]/mod3;
  15.                 vector S = {s*sji.x,s*sji.y,s*sji.z};
  16.                 accelerations[i].x+=S.x;accelerations[i].y+=S.y;accelerations[i].z+=S.z;
  17.             }
  18.         }
  19.     }
  20. }
  22. 我尝试做类似的事情：
  24. ```c
  25. void computeAccelerations_static(int num_of_threads){
  26. int i,j;
  27. #pragma omp parallel for num_threads(num_of_threads) schedule(static)
  28. for(i=0;i<bodies;i++){
  29.     accelerations[i].x = 0; accelerations[i].y = 0; accelerations[i].z = 0;
  30.     for(j=0;j<bodies;j++){
  31.         if(i!=j){
  32.             //accelerations[i] = addVectors(accelerations[i],scaleVector(GravConstant*masses[j]/pow(mod(subtractVectors(positions[i],positions[j])),3),subtractVectors(positions[j],positions[i])));
  33.             vector sij = {positions[i].x-positions[j].x,positions[i].y-positions[j].y,positions[i].z-positions[j].z};
  34.             vector sji = {positions[j].x-positions[i].x,positions[j].y-positions[i].y,positions[j].z-positions[i].z};
  35.             double mod = sqrt(sij.x*sij.x + sij.y*sij.y + sij.z*sij.z);
  36.             double mod3 = mod * mod * mod;
  37.             double s = GravConstant*masses[j]/mod3;
  38.             vector S = {s*sji.x,s*sji.y,s*sji.z};
  39.             accelerations[i].x+=S.x;accelerations[i].y+=S.y;accelerations[i].z+=S.z;
  40.         }
  41.     }
  42. }

自然而然地添加#pragma omp parallel for num_threads(num_of_threads) schedule(static)，但这并不正确。

我认为accelerations[i]存在一些虚假共享问题，但我不知道如何解决它。我感激任何帮助。谢谢。

  2. <details>
  3. <summary>英文:</summary>
  4. So I have this function that I have to parallelize with OpenMP static scheduling for n threads
  5. void computeAccelerations(){
  6. int i,j;
  7. for(i=0;i&lt;bodies;i++){
  8. accelerations[i].x = 0;	accelerations[i].y = 0; accelerations[i].z = 0;
  9. for(j=0;j&lt;bodies;j++){
  10. if(i!=j){
  11. //accelerations[i] = addVectors(accelerations[i],scaleVector(GravConstant*masses[j]/pow(mod(subtractVectors(positions[i],positions[j])),3),subtractVectors(positions[j],positions[i])));
  12. vector sij = {positions[i].x-positions[j].x,positions[i].y-positions[j].y,positions[i].z-positions[j].z};
  13. vector sji = {positions[j].x-positions[i].x,positions[j].y-positions[i].y,positions[j].z-positions[i].z};
  14. double mod = sqrt(sij.x*sij.x + sij.y*sij.y + sij.z*sij.z);
  15. double mod3 = mod * mod * mod;
  16. double s = GravConstant*masses[j]/mod3;
  17. vector S = {s*sji.x,s*sji.y,s*sji.z};
  18. accelerations[i].x+=S.x;accelerations[i].y+=S.y;accelerations[i].z+=S.z;
  19. }
  20. }
  21. }
  22. }
  23. I tried to do something like:

void computeAccelerations_static(int num_of_threads){
int i,j;
#pragma omp parallel for num_threads(num_of_threads) schedule(static)
for(i=0;i<bodies;i++){
accelerations[i].x = 0; accelerations[i].y = 0; accelerations[i].z = 0;
for(j=0;j<bodies;j++){
if(i!=j){
//accelerations[i] = addVectors(accelerations[i],scaleVector(GravConstantmasses[j]/pow(mod(subtractVectors(positions[i],positions[j])),3),subtractVectors(positions[j],positions[i])));
vector sij = {positions[i].x-positions[j].x,positions[i].y-positions[j].y,positions[i].z-positions[j].z};
vector sji = {positions[j].x-positions[i].x,positions[j].y-positions[i].y,positions[j].z-positions[i].z};
double mod = sqrt(sij.xsij.x + sij.ysij.y + sij.zsij.z);
double mod3 = mod * mod * mod;
double s = GravConstantmasses[j]/mod3;
vector S = {ssji.x,ssji.y,ssji.z};
accelerations[i].x+=S.x;accelerations[i].y+=S.y;accelerations[i].z+=S.z;
}
}
}

  2. It comes naturally to just add the ```#pragma omp parallel for num_threads(num_of_threads) schedule(static)``` but it isn&#39;t correct. 
  3. I think there is some kind of false sharing with the ``accelerations[i]`` but I don&#39;t know how to approach it. I appreciate any kind of help. Thank you.
  4. </details>
  5. # 答案1
  6. **得分**: 2
  7. 在你的循环嵌套中，只有外部循环的迭代被并行化。因为`i`是循环控制变量，每个线程都有自己的私有副本，但从代码风格的角度来看，最好在循环控制块中声明`i`。
  8. `j`则不同。它在并行区域之外声明，*不是*并行化循环的控制变量。因此，它在线程之间是共享的。因为执行`i`循环迭代的每个线程都操作共享变量`j`，所以会出现数据竞争的严重问题。这可以通过将`j`的声明移到并行区域内来解决（除其他替代方案之外），最好放在其关联循环的控制块中。
  9. 总之，代码应该是这样的：
  10. ```cpp
  11. // int i, j;
  12. #pragma omp parallel for num_threads(num_of_threads) schedule(static)
  13. for (int i = 0; i < bodies; i++) {
  14. accelerations[i].x = 0;
  15. accelerations[i].y = 0;
  16. accelerations[i].z = 0;
  17. #pragma omp parallel for num_threads(num_of_threads) schedule(static)
  18. for (int j = 0; j < bodies; j++) {
  19. if (i != j) {
  20. vector sij = { positions[i].x - positions[j].x,
  21. positions[i].y - positions[j].y,
  22. positions[i].z - positions[j].z };
  23. double mod = sqrt(sij.x * sij.x + sij.y * sij.y + sij.z * sij.z);
  24. double mod3 = mod * mod * mod;
  25. double s = GravConstant * masses[j] / mod3;
  26. accelerations[i].x -= s * sij.x;
  27. accelerations[i].y -= s * sij.y;
  28. accelerations[i].z -= s * sij.z;
  29. }
  30. }
  31. }

还要注意，计算sji似乎是浪费的，因为在数学上它只是-sij，而且既不sji也不sij被修改。我建议将上述代码简化为以下形式：

  1. #pragma omp parallel for num_threads(num_of_threads) schedule(static)
  2. for (int i = 0; i < bodies; i++) {
  3.     accelerations[i].x = 0;
  4.     accelerations[i].y = 0;
  5.     accelerations[i].z = 0;
  7.     for (int j = 0; j < bodies; j++) {
  8.         if (i != j) {
  9.             vector sij = { positions[i].x - positions[j].x,
  10.                            positions[i].y - positions[j].y,
  11.                            positions[i].z - positions[j].z };
  12.             double mod = sqrt(sij.x * sij.x + sij.y * sij.y + sij.z * sij.z);
  13.             double mod3 = mod * mod * mod;
  14.             double s = GravConstant * masses[j] / mod3;
  16.             accelerations[i].x -= s * sij.x;
  17.             accelerations[i].y -= s * sij.y;
  18.             accelerations[i].z -= s * sij.z;
  19.         }
  20.     }
  21. }
英文:

In your loop nest, only the iterations of the outer loop are parallelized. Because i is the loop-control variable, each thread gets its own, private copy, but as a matter of style, it would be better to declare i in the loop control block.

j is another matter. It is declared outside the parallel region and it is not the control variable of a parallelized loop. As a result, it is shared among the threads. Because each of the threads executing i-loop iterations manipulates shared variable j, you have a huge problem with data races. This would be resolved (among other alternatives) by moving the declaration of j into the parallel region, preferrably into the control block of its associated loop.

Overall, then:

  1. // int i, j;
  2. #pragma omp parallel for num_threads(num_of_threads) schedule(static)
  3. for (int i = 0; i &lt; bodies; i++) {
  4. accelerations[i].x = 0;
  5. accelerations[i].y = 0;
  6. accelerations[i].z = 0;
  7. for (int j = 0; j &lt; bodies; j++) {
  8. if (i != j) {
  9. //accelerations[i] = addVectors(accelerations[i],scaleVector(GravConstant*masses[j]/pow(mod(subtractVectors(positions[i],positions[j])),3),subtractVectors(positions[j],positions[i])));
  10. vector sij = { positions[i].x - positions[j].x,
  11. positions[i].y - positions[j].y,
  12. positions[i].z - positions[j].z };
  13. vector sji = { positions[j].x - positions[i].x,
  14. positions[j].y - positions[i].y,
  15. positions[j].z - positions[i].z };
  16. double mod = sqrt(sij.x * sij.x + sij.y * sij.y + sij.z * sij.z);
  17. double mod3 = mod * mod * mod;
  18. double s = GravConstant * masses[j] / mod3;
  19. vector S = { s * sji.x, s * sji.y, s * sji.z };
  20. accelerations[i].x += S.x;
  21. accelerations[i].y += S.y;
  22. accelerations[i].z += S.z;
  23. }
  24. }
  25. }

Note also that computing sji appears to be wasteful, as in mathematical terms it is just -sij, and neither sji nor sij is modified. I would probably reduce the above to something more like this:

  1. #pragma omp parallel for num_threads(num_of_threads) schedule(static)
  2. for (int i = 0; i &lt; bodies; i++) {
  3. accelerations[i].x = 0;
  4. accelerations[i].y = 0;
  5. accelerations[i].z = 0;
  6. for (int j = 0; j &lt; bodies; j++) {
  7. if (i != j) {
  8. vector sij = { positions[i].x - positions[j].x,
  9. positions[i].y - positions[j].y,
  10. positions[i].z - positions[j].z };
  11. double mod = sqrt(sij.x * sij.x + sij.y * sij.y + sij.z * sij.z);
  12. double mod3 = mod * mod * mod;
  13. double s = GravConstant * masses[j] / mod3;
  14. accelerations[i].x -= s * sij.x;
  15. accelerations[i].y -= s * sij.y;
  16. accelerations[i].z -= s * sij.z;
  17. }
  18. }
  19. }

huangapple
  • 本文由 发表于 2023年2月16日 04:05:40
  • 转载请务必保留本文链接：https://go.coder-hub.com/75464941.html
  • c
  • false-sharing
  • multithreading
  • openmp
  • parallel-processing
匿名

发表评论

匿名网友

:?: :razz: :sad: :evil: :!: :smile: :oops: :grin: :eek: :shock: :???: :cool: :lol: :mad: :twisted: :roll: :wink: :idea: :arrow: :neutral: :cry: :mrgreen:

确定