MPI並列の使い方

# 基本的な仕組み
例えば３つのコアで並列化することを考えます。
すると次の図のように実行されます。

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2445529/21a47c4c-9f1c-56b9-09a0-b2bae6e5eb11.png)

ただ並列化をしても同じプログラムが３つ同時に実行されるだけで、全く同じ処理が３回行われるだけです。

MPIのパッケージでは使用するコアに番号を振って自分の番号を読み取る関数があります。そこでこの番号に応じて異なる処理を行うプログラムを書いて、それによって処理を分担させるのがMPI並列化のやり方です。

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2445529/e7eb2f1a-e292-cd1a-0105-96efbf695783.png)

# コーデイングとコンパイル

mpi並列化を使うにはパッケージmpi.hのインストールが必要です。

```test_mpi.cpp
#include <mpi.h>
:
:
```

また、コンパイルもg++やi++ではなく使用するマシンに応じた並列化用のコンパイラが必要です。例えば

```
mpic++ test_mpi.cpp
```

などのようにコンパイルします。(使うコンパイラは適宜確認が必要です。)
実行の際もmpirunなど、マシンごとに並列化してプログラムを実行するコマンドを使用する必要があります。

# 並列化の宣言

MPI並列処理をする際には次の部分が必要です。

```mpi.cpp
#include <mpi.h> //MPIのパッケージの読み込み

int main(int argc, char **argv){
    :
    :
    int procs, rank;
    MPI_Init( &argc, &argv ); //並列化開始の宣言(初期化や番号の割り当てなどを行う)
    MPI_Comm_size( MPI_COMM_WORLD, &procs ); //並列化に使用するコア数
  　MPI_Comm_rank( MPI_COMM_WORLD, &rank ); //自分が何番目のコアかを表す番号(0　~　procs-1番目)
    ///////////////////////
    // 並列処理を含む処理 //
    ///////////////////////
    MPI_Finalize(); //並列化終了の宣言
    :
    :
    return 0;
}
```

# 割り振り

全部で$N$回の操作を繰り返すとき、各ノードが何番目から何番目までの計算を担当するかを決めなければいけません。(なるべく均等になるように分けます。そうでないと並列による時間短縮の効果が不十分になってしまいます。)

それを計算するための関数を作りました。一つ目は$r$番目のノードが担当する回数$n_r$を計算し、２つ目はその累積和を求めて、担当する最初のループのインデックスを計算しています。

```mpi.cpp
int get_num_of_points_for_proc( int size, int procs, int rank ){
  int n_points = size/procs;
  if ( rank < size%procs ) n_points++;
  return n_points;
}

int get_start_point_for_proc( int size, int procs, int rank ){
  int start_point = 0;
  for( int i=0; i<rank; i++ ) start_point += get_num_of_points_for_proc( size, procs, i );
  return start_point;
}
```

# 集計と出力の注意点

### 集計

それぞれのノードは勝手に計算をしているだけなので、他のノードがどこまで計算を進めているのかを気にしません。
各ノードが部分的に計算した結果を集計する際に他のノードの処理がまだ終わっていないのに、データを他のノードに送受信したり勝手に集計を始めたりします。
そこで、次の「他のノードが同じところに到達するまでその先に進まない」という命令が必要です。

```mpi.cpp
MPI_Barrier( MPI_COMM_WORLD );
//集計orデータの送受信
```

### 出力

何もしないと出力結果も並列化されて各ノードが順番も気にせずに好き勝手に出力をしてしまうので、1つのノード(ふつうは0番ノード)にデータを送って0番ノードだけが出力を行うようにします。新規ファイルやディレクトリの作成も各ノードが独立してやるので、ノードを指定しないと変なタイミングで初期化されたり上書きされたりするので注意が必要です。

```mpi.cpp
int rank;
MPI_Init( &argc, &argv ); //並列化開始の宣言
MPI_Comm_rank( MPI_COMM_WORLD, &rank );
:
//並列処理
:
//集計やデータの送受信
:
if(rank==0){ //0番ノードを指定
    //データの出力
}
:
MPI_Finalize(); //並列化終了の宣言
```

# 集計

各コアでのデータを通信でやりとりしてそれぞれのコア(or指定したコア)で所望のデータを共有する関数です。様々な関数がありますが、使用例の一部をご紹介します。

# 例:結果の総和

次のコードで
```math
\mathrm{sum} = \sum_{x=0}^{99} x = \sum_{r=0}^{\mathrm{procs}-1} \sum_{x = 0}^{n_\mathrm{assign}-1} (x + i_\mathrm{start}) = 4950
```

を計算しました。

```mpi_sum.cpp
#include <mpi.h>
#include <stdio.h>
#include <stdlib.h>

#include <sstream>
#include <fstream>
#include <iostream>

using namespace std;

int get_num_of_points_for_proc( int size, int procs, int rank ){
  int n_points = size/procs;
  if ( rank < size%procs ) n_points++;
  return n_points;
}

int get_start_point_for_proc( int size, int procs, int rank ){
  int start_point = 0;
  for( int i=0; i<rank; i++ ) start_point += get_num_of_points_for_proc( size, procs, i );
  return start_point;
}

int calc_sum(int i_start, int n_resign){
  int sum = 0;
  for( int i=i_start; i-i_start<n_sign; ++i) sum += i;
  return sum;
}

int main(int argc, char** argv) {
  int procs, rank;
  int sum, size = 100;
  
  MPI_Init( &argc , &argv );
  MPI_Comm_size( MPI_COMM_WORLD , &procs );
  MPI_Comm_rank( MPI_COMM_WORLD , &rank );

  int i_start = get_start_point_for_proc(size, procs, rank);
  int n_resign = get_num_of_points_for_proc(size, procs, rank);
//部分和の計算
  int subsum = calc_sum(i_start, n_resign);

  MPI_Barrier( MPI_COMM_WORLD );
// 各部分和の総和
  MPI_Allreduce(&subsum, &sum, 1, MPI_INT, MPI_SUM, MPI_COMM_WORLD);
//0番ノードが出力
  if(rank==0) cout << procs << " " << sum << endl;
  MPI_Finalize();
  return 0;
}
```

出力結果は

```
10 4950
```

となり、10並列で計算した結果は確かに合っています。

## 配列の各要素を集計

各コアで

```math
\text{data} = (\text{rank},\text{rank},\cdots,\text{rank})
```

という配列を用意して、各要素で全コアの総和をとった配列

```math
\text{sum} = \frac{\text{procs}(\text{procs}-1)}{2} \\
\text{data} = (\text{sum},\text{sum},\cdots,\text{sum})
```

```mpi_test.cpp
#include <mpi.h>
#include <stdio.h>
#include <stdlib.h>

#include <sstream>
#include <fstream>
#include <iostream>

using namespace std;

int main(int argc, char** argv) {
  int procs, rank;
  int size = 10;
  int *data = new int[size]; 
  MPI_Init( &argc , &argv );
  MPI_Comm_size( MPI_COMM_WORLD , &procs );
  MPI_Comm_rank( MPI_COMM_WORLD , &rank );  
  for(int i=0; i<size; ++i) data[i] = rank;
  MPI_Barrier( MPI_COMM_WORLD );
  MPI_Allreduce( MPI_IN_PLACE, data, size, MPI_INT, MPI_SUM, MPI_COMM_WORLD );
  if(rank==0){
    cout << procs << endl;
    for(int i=0; i<size; ++i) cout << data[i] << " ";
    cout << endl;
  }
  MPI_Finalize();
  delete[] data;
  return 0;
}
```
出力は次のようになります。

```
10
45 45 45 45 45 45 45 45 45 45 
```

$\text{procs}=10$なので、結果は$45$で合っています。