OpenMP超入門その1 #C++

OpenMPとは?

OpenMPのことをWikipediaさんで調べると...

OpenMPは、並列コンピューティング環境を利用するために用いられる標準化された基盤。OpenMPは主に共有メモリ型並列計算機で用いられる。OpenMP:Wikipedia

とかかれてますね...

なんのこっちゃい, というと,
並列計算するともっとはやくなるお!! 拙者を使えばメモリも共有して並列計算できるから便利だお
というわけです!!

対応言語としてはC/C++, fortranが対応しています!!

OpenMPを使う羽目になったきっかけ

実はPythonで並列化できるだろうと思っていたら, メモリ共有に関して自分で適宜したClassではむりそうだったからです...

ホントは全部Pythonで書いてJupyterにまとめようと思ったのですが僕の見通しが悪かったみたいです.

環境

Windows
- Cygwinで確認できてます!!
Mac
- 悪いこといわないからbrew install gccとかしてgccに切り替えよう(なぜかMacではgccがclangのエイリアスになってる...)
Linux系統
- 確認中(たぶんいける)

すべての環境に関していえること

gccは4.2以降

というか, 新しいバージョンであれば大丈夫, という認識でOKです...

Step1 : OpenMPが動くか確認

今回, 実行についてC++, g++ 5,4,0でやっております...

Hello World with OpenMP

比較的どこのサイトでもやっているやつ

#include <iostream>
#include <omp.h> // OpenMPのヘッダをinclude

int main(){
//ここが並列処理される
  #pragma omp parallel
  {
    std::cout << "Hello World!\n";
  }
}

これを実行すると

Hello World!
Hello World!
Hello World!
Hello World!

と表示される. このままだと, 異なるスレッドに全く同じ処理を任せていることになる.

配列を操作してみる

#include <iostream>
#include <omp.h>

int main(){
  int A[10];
  #pragma omp parallel for
  for(int i = 0; i < 10; i++){
    A[i] = i;
  }

  for(int i = 0; i < 10; i++){
    std::cout << A[i] << "\n";
  }
}

上のようなプログラムがあったとすると,

と表示してくれる. もちろん, 並列化されている.

単純に並列化してみる

今から, 1,...1000を足してみる. 500500が答えだ.
普通に書くと

#include <iostream>

int main(){
  int sum = 0;
  for(int i = 0; i < 1001; i++){
    sum += i;
  }

  std::cout << "sum = " << sum <<"\n";

}

とかくはずだ. もちろん, これは正しい. しかし, これを並列化する際に

#include <iostream>
#include <omp.h>

int main(){
  int sum = 0;
  #pragma omp parallel for
  for(int i = 0; i < 1001; i++){
    sum += i;
  }

  std::cout << "sum = " << sum <<"\n";

}

と書いてしまうと, 安定しない. 500500にならないときもある.
正しくは,

#include <iostream>
#include <omp.h>

int main(){
  int sum = 0;
  #pragma omp parallel for reduction(+:sum)
  for(int i = 0; i < 1001; i++){
    sum += i;
  }

  std::cout << "sum = " << sum <<"\n";

}

ひとつの変数に足しこむときはreductionでその変数を設定する必要がある.

もっと実用的っぽいコードを書いてみる...

ということで数値積分をやってみます...
eを求めてみましょう!!

e - 1 = \int_{0}^1 e^x dx

なので, これを元に求めてみます.

逐次計算だと...

#include <iostream>
#include <cmath>

#define N 10000

int main(){
  double step = 1.0 / N;
  double temp;
  double ans;
  for(int i = 0; i < N; i++){
    temp = (i + 0.5) * step;
    ans += exp(temp)/N;
  }
  std::cout << "ans = " << ans + 1 << "\n";
}

とかけます. 結果は

ans = 2.71828

です.

並列計算でやると...

いちばんやっちゃいけないやつ

#include <iostream>
#include <cmath>
#include <omp.h>
#define N 10000

int main(){
  double step = 1.0 / N;
  double temp;
  double ans;
  #pragma omp parallel for
  for(int i = 0; i < N; i++){
    temp = (i + 0.5) * step;
    ans += exp(temp)/N;
  }
  std::cout << "ans = " << ans + 1 << "\n";
}

動かすと

ans = 2.0528

もちろん, 間違ってます.

惜しい!!って感じの間違え

#include <iostream>
#include <cmath>
#include <omp.h>
#define N 10000

int main(){
  double step = 1.0 / N;
  double temp;
  double ans;
  #pragma omp parallel for reduction(+:ans)
  for(int i = 0; i < N; i++){
    temp = (i + 0.5) * step;
    ans += exp(temp)/N;
  }
  std::cout << "ans = " << ans + 1 << "\n";
}

まぁ, 動くには動くけど, お作法的にはお行儀がよくない. ちゃんと, Loop内で書き換えられる変数(用は使い捨ての変数とか)はちゃんとprivateしてあげたほうがいい. (@yohhoy 様の指摘ありましたので変更します. )

実は, このままだと変数tempはshared変数ではないため, 実行に際して動作の保証ができないようです. そのため, 以下のようにしましょう...

ベストアンサー

#include <iostream>
#include <cmath>
#include <omp.h>
#define N 10000

int main(){
  double step = 1.0 / N;
  double temp;
  double ans;
  #pragma omp parallel for reduction(+:ans) private(temp)
  for(int i = 0; i < N; i++){
    temp = (i + 0.5) * step;
    ans += exp(temp)/N;
  }
  std::cout << "ans = " << ans + 1 << "\n";
}

こう書いてあげたほうがいい. (パフォーマンスに響くんだとか)

最後に

ほんまにさわりだけやりました. たぶんこれより先の内容は別の記事にします!!

OpenMP超入門 その1