OpenMP 並列プログラミング2の続きです。

実際にOpenMPを使いこなしていきましょう。

グレゴリー・ライプニッツ級数でπを求めます。

π / 4 = arctan1 
 = (1 / 1) - (1 / 3) + (1 / 5) - (1 / 7) + (1 /9) - (1 / 11) ...無限に続く
 
グレゴリー・ライプニッツ級数は収束が遅いので実際にπを求める際には使われていません。

式がシンプルなforループで記述できるので並列プログラミング向きなのです。

並列化なし pi-leibniz.c
#include <stdio.h>
#include <time.h>

#define MAX 1000000000

int main()
{
int i;
double is = 1, a = 0;

clock_t start,end;
  start = clock();
for (i = 0; i < MAX; i++) {
if ((i % 2) == 0) {
is = 1;
} else {
is = -1;
}
a += is / (2 * i + 1); 
}
printf("%.15f\n", 4 *a);
end = clock();
  printf(" time [sec.] = %f \n", (double)(end - start) / CLOCKS_PER_SEC);
return 0;
}

$ gcc pi-leibniz.c -o pi-leibniz
$ ./pi-leibniz
3.141592652588050
 time [sec.] = 6.853473 
6.8秒かかりました。

OpenMPで並列化 pi-leibniz2.c
#include <stdio.h>
#include <time.h>
#include <omp.h>

#define MAX 1000000000

int main()
{
int i;
double is = 1, a = 0;

double start,end;
  start = omp_get_wtime();

#pragma omp parallel for private(is) reduction(+:a)
for (i = 0; i < MAX; i++) {
if ((i % 2) == 0) {
is = 1;
} else {
is = -1;
}
a += is / (2 * i + 1); 
}
printf("%.15f\n", 4 *a);
end = omp_get_wtime();
  printf(" time [sec.] = %lf \n", end - start);

return 0;
}

$ gcc pi-leibniz2.c -o pi-leibniz2 -fopenmp
$ ./pi-leibniz2
3.141592652589210
 time [sec.] = 1.899626 
 
並列化をすると1.8秒。 すごい効果が出ています。コア４つなので、早く終わります。

ところがpiの結果が異なります。

並列なし　 3.141592652588050
並列あり　 3.141592652589210

微妙に異なります。MAX 10にして調べました。

並列なし　ループ毎に加算しています。
0 : 1.000000000000000
1 : 0.666666666666667
2 : 0.866666666666667
3 : 0.723809523809524
4 : 0.834920634920635
5 : 0.744011544011544
6 : 0.820934620934621
7 : 0.754267954267954
8 : 0.813091483679719
9 : 0.760459904732351
3.041839618929403

並列化あり　
6 : 0.076923076923077
7 : -0.073777203188968
8 : 0.135746606334842
9 : -0.126408782136336
3 : -0.007110536522301
4 : -0.015297671025225
5 : -0.106206761934316
0 : 0.893793238065684
1 : 0.560459904732351
2 : 0.760459904732351
3.041839618929402
スレッド毎に加算して最後にまとめ上げています。
当たり前といえば当たり前か。そういう仕組みだからしょうがないのか。

その辺りに誤差の原因がありそうですね。

誤差のレベルなので無視してもいいのかもしれませんが、わかりません。

並列なしと並列ありで突き合せて、確かめた方がいいですね。注意が必要です。