OpenMPで入れ子を使う

C言語のOpenMPで入れ子を使う際のメモ。既存のOpenMPでスレッド並列化されている関数をsections構文内で呼び出すときや再帰処理で使える。
このページの説明で使っているコンパイラは、

```sh
$ gcc --version
gcc (GCC) 9.1.0
Copyright (C) 2019 Free Software Foundation, Inc.
This is free software; see the source for copying conditions.  There is NO
warranty; not even for MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
```

OpenMPのバージョンは、

```sh
$ echo |cpp -fopenmp -dM |grep -i open
#define _OPENMP 201511
```
なので OpenMP4.5。入れ子自体はOpenMP2.5から使えるようになっているらしい。

## 使い方

```c
omp_set_nested(1);
```
引数に1を指定してnestedフラグをオンにする。0にすればフラグをオフにする。現在のフラグがどうなっているかは、

```c
int omp_get_nested()
```
で取得可能。1ならオン、0ならオフ。オフの場合のときの入れ子領域は1スレッドし、parallel regionは無視される。


### 動作確認
下記の2段階入れ子のコードで動作を確認してみる。

```c:omp_nest2.c
#include <stdio.h>
#include <stdlib.h>
#include <omp.h>

void print_omp_nested()
{
  int omp_active_level_l0 = omp_get_active_level(); // 現在の入れ子レベルの取得
  int omp_max_threads_l0  = omp_get_max_threads();　// 最大スレッド数(parallel regionを作ったときのデフォルト並列数)の取得　
  int omp_num_threads_l0  = omp_get_num_threads();  // 今入っているparallel regionでの並列数の取得
  int omp_thread_l0       = omp_get_thread_num();   // 今入っているparallel regionでのスレッドIDの取得

  printf("# level : max,num_threads : omp_thread(l0,l1,l2)\n");
  printf("# %2d    :   %3d,%3d       : %3d\n", omp_active_level_l0,
	 omp_max_threads_l0, omp_num_threads_l0, omp_thread_l0);
  fflush(stdout);

#pragma omp parallel
  {
    int omp_active_level_l1 = omp_get_active_level();
    int omp_max_threads_l1  = omp_get_max_threads();
    int omp_num_threads_l1  = omp_get_num_threads();
    int omp_thread_l1       = omp_get_thread_num();

    printf("# %2d    :   %3d,%3d       : %3d,%3d\n", omp_active_level_l1,
	   omp_max_threads_l1,omp_num_threads_l1,
	   omp_thread_l0, omp_thread_l1);
    fflush(stdout);

#pragma omp parallel
    {   
        int omp_active_level_l2 = omp_get_active_level();
        int omp_max_threads_l2  = omp_get_max_threads();
        int omp_num_threads_l2  = omp_get_num_threads();
        int omp_thread_l2       = omp_get_thread_num();

	printf("# %2d    :   %3d,%3d       : %3d,%3d,%3d\n", omp_active_level_l2,
	       omp_max_threads_l2, omp_num_threads_l2,
	       omp_thread_l0, omp_thread_l1, omp_thread_l2);
	fflush(stdout);

    }
  }
}

int main(int argc,char** argv)
{
  printf("# omp_get_num_procs         = %d\n", omp_get_num_procs());
  printf("# omp_get_max_threads       = %d\n", omp_get_max_threads());
  printf("# omp_get_num_threads       = %d\n", omp_get_num_threads());
  printf("# omp_get_max_active_levels = %d\n", omp_get_max_active_levels());
  printf("# omp_get_thread_limit      = %d\n", omp_get_thread_limit());

  // ここで入れ子をオンにする。デフォルトはオフ。
  omp_set_nested(1);

  print_omp_nested();

  // 入れ子処理が終わったら一応デフォルトに戻しておく。
  omp_set_nested(0);

  return EXIT_SUCCESS;
}
```

上記をコンパイルし、実行する。実行環境では4コア8スレッドなのでデフォルトで8スレッドが使われる。
最初に `omp_set_nested(1);` を `omp_set_nested(0);` にして実行してみると、

```sh
$ ./a.out 
# omp_get_num_procs         = 8
# omp_get_max_threads       = 8
# omp_get_num_threads       = 1
# omp_get_max_active_levels = 2147483647
# omp_get_thread_limit      = 2147483647
# level : max,num_threads : omp_thread(l0,l1,l2)
#  0    :     8,  1       :   0
#  1    :     8,  8       :   0,  0
#  1    :     8,  8       :   0,  2
#  1    :     8,  8       :   0,  4
#  1    :     8,  8       :   0,  3
#  1    :     8,  8       :   0,  5
#  1    :     8,  8       :   0,  7
#  1    :     8,  8       :   0,  1
#  1    :     8,  1       :   0,  0,  0
#  1    :     8,  1       :   0,  2,  0
#  1    :     8,  1       :   0,  5,  0
#  1    :     8,  1       :   0,  1,  0
#  1    :     8,  1       :   0,  4,  0
#  1    :     8,  1       :   0,  3,  0
#  1    :     8,  8       :   0,  6
#  1    :     8,  1       :   0,  6,  0
#  1    :     8,  1       :   0,  7,  0
```
レベル0はparallel regionに入る前。レベル1で8スレッド並列になるが、2つ目のparallel regionが無視されており、新たにレベル2でのスレッドが生成されることはない。nested flagを1にして実行すると、

```sh
$ gcc -fopenmp omp_nest2.c
$ ./a.out                        
...
# level : max,num_threads : omp_thread(l0,l1,l2)
#  0    :     8,  1       :   0
#  1    :     8,  8       :   0,  0
#  1    :     8,  8       :   0,  7
#  1    :     8,  8       :   0,  1
#  1    :     8,  8       :   0,  3
#  1    :     8,  8       :   0,  6
#  1    :     8,  8       :   0,  5
#  1    :     8,  8       :   0,  2
#  1    :     8,  8       :   0,  4
#  2    :     8,  8       :   0,  4,  1
#  2    :     8,  8       :   0,  4,  2
#  2    :     8,  8       :   0,  4,  3
...省略...
#  2    :     8,  8       :   0,  5,  5
#  2    :     8,  8       :   0,  1,  5
```
レベル1では8スレッド、レベル2ではレベル1での各スレッドがそれぞれ8スレッド生成していることがわかる。今はレベル2までなので最大8*8=64スレッド生成されているがもう1レベル増えると512スレッドとなりスレッド数が膨大になるため、入れ子を使わないとき以上にスレッド数のコントロールを意識する必要がある。

### 環境変数`OMP_NUM_THREADS`によるスレッド数の指定

実行時に全てのレベルで一括にスレッド数を変更するには`OMP_NUM_THREADS=x`を指定する。また、各レベルでスレッド数を変更したい場合は`OMP_NUM_THREADS=x,y,z`のようにカンマ区切りで数値を指定すると各レベルで使われるスレッド数を制限できる。実行コードは上記と同じ。

```sh
$ OMP_NUM_THREADS=3 ./a.out        
# omp_get_num_procs         = 8
# omp_get_max_threads       = 3
...
# level : max,num_threads : omp_thread(l0,l1,l2)
#  0    :     3,  1       :   0
#  1    :     3,  3       :   0,  0
#  1    :     3,  3       :   0,  1
#  1    :     3,  3       :   0,  2
#  2    :     3,  3       :   0,  0,  0
#  2    :     3,  3       :   0,  0,  2
...省略...
#  2    :     3,  3       :   0,  2,  2
#  2    :     3,  3       :   0,  1,  2
全てのレベルで3スレッドづつ生成。

$ OMP_NUM_THREADS=2,4 ./a.out
# omp_get_max_threads       = 2
# level : max,num_threads : omp_thread(l0,l1,l2)
#  0    :     2,  1       :   0
#  1    :     4,  2       :   0,  0
#  1    :     4,  2       :   0,  1
#  2    :     4,  4       :   0,  0,  0
#  2    :     4,  4       :   0,  0,  1
#  2    :     4,  4       :   0,  0,  2
#  2    :     4,  4       :   0,  0,  3
#  2    :     4,  4       :   0,  1,  1
#  2    :     4,  4       :   0,  1,  2
#  2    :     4,  4       :   0,  1,  3
#  2    :     4,  4       :   0,  1,  0
レベル0の`omp_get_max_threads`が2なのでレベル1で2スレッド生成。
レベル1の`omp_get_max_threads`が4なのでレベル2で4スレッドづつ生成。

$ OMP_NUM_THREADS=4,2 ./a.out
# omp_get_num_procs         = 8
# omp_get_max_threads       = 4
# level : max,num_threads : omp_thread(l0,l1,l2)
#  0    :     4,  1       :   0
#  1    :     2,  4       :   0,  0
#  1    :     2,  4       :   0,  1
#  1    :     2,  4       :   0,  2
#  2    :     2,  2       :   0,  0,  0
#  2    :     2,  2       :   0,  0,  1
#  1    :     2,  4       :   0,  3
#  2    :     2,  2       :   0,  1,  0
#  2    :     2,  2       :   0,  3,  0
#  2    :     2,  2       :   0,  2,  0
#  2    :     2,  2       :   0,  2,  1
#  2    :     2,  2       :   0,  3,  1
#  2    :     2,  2       :   0,  1,  1
レベル0の`omp_get_max_threads`が4なのでレベル1で4スレッド生成。
レベル1の`omp_get_max_threads`が2なのでレベル2で2スレッドづつ生成。

$ OMP_NUM_THREADS=4,1 ./a.out
# omp_get_num_procs         = 8
# omp_get_max_threads       = 4
#  0    :     4,  1       :   0
#  1    :     1,  4       :   0,  0
#  1    :     1,  1       :   0,  0,  0
#  1    :     1,  4       :   0,  2
#  1    :     1,  4       :   0,  3
#  1    :     1,  4       :   0,  1
#  1    :     1,  1       :   0,  2,  0
#  1    :     1,  1       :   0,  3,  0
#  1    :     1,  1       :   0,  1,  0
レベル0の`omp_get_max_threads`が4なのでレベル1で4スレッド生成。
レベル1の`omp_get_max_threads`が1なのでレベル2のスレッド階層は生成されないで
レベル1のスレッド階層のまま次のparallel regionを無視していることがわかる。
```
シェルで`$ export OMP_NUM_THREADS=x`で指定してもいい。`omp_get_max_threads()`で取得できるスレッド数は今の入れ子レベルでの(次の入れ子レベルが並列実行される)スレッド数`omp_get_num_threads()`であることがわかる。

### コード内でスレッド数の指定

環境変数でなくてもコード内で使うスレッド数を変更できるようにする。`print_omp_nested()`関数をparallel構文でスレッド数を指定するように変更。スレッド数の優先度は環境変数指定よりこちらの方が高い。

```c
void print_omp_nested()
{
   // (変数の定義などは省略)
  ...
  printf("# level : max,num_threads : omp_thread(l0,l1,l2)\n");
  printf("# %2d    :   %3d,%3d       : %3d\n", omp_active_level_l0,
         omp_max_threads_l0, omp_num_threads_l0, omp_thread_l0);
  fflush(stdout);

#pragma omp parallel num_threads(4)
  {
    ...
    printf("# %2d    :   %3d,%3d       : %3d,%3d\n", omp_active_level_l1,
           omp_max_threads_l1,omp_num_threads_l1,
           omp_thread_l0, omp_thread_l1);
    fflush(stdout);

#pragma omp parallel num_threads(2)
    {
        ...
        printf("# %2d    :   %3d,%3d       : %3d,%3d,%3d\n", omp_active_level_l2,
               omp_max_threads_l2, omp_num_threads_l2,
               omp_thread_l0, omp_thread_l1, omp_thread_l2);
        fflush(stdout);

    }
  }
}
```

上記コードを実行すると

```sh
# omp_get_num_procs         = 8
# omp_get_max_threads       = 8
# level : max,num_threads : omp_thread(l0,l1,l2)
#  0    :     8,  1       :   0
#  1    :     8,  4       :   0,  0
#  1    :     8,  4       :   0,  2
#  1    :     8,  4       :   0,  3
#  1    :     8,  4       :   0,  1
#  2    :     8,  2       :   0,  0,  0
#  2    :     8,  2       :   0,  0,  1
...省略...
#  2    :     8,  2       :   0,  1,  1
#  2    :     8,  2       :   0,  3,  1
```
想定通りの振る舞いをするが`omp_get_max_threads`がデフォルトのスレッド数の値になっている。これは入れ子関係なくparallel構文でスレッド数を指定したらこうなるので不具合というわけではない。

`omp_set_num_threads()`でも指定ができる。これを実行すると、

```c
 omp_set_num_threads(4);
#pragma omp parallel
  {
   ...    
   omp_set_num_threads(2);
   #pragma omp parallel
    {
```

```sh
./a.out                        
# omp_get_num_procs         = 8
# omp_get_max_threads       = 8
## level : max,num_threads : omp_thread(l0,l1,l2)
#  0    :     8,  1       :   0
#  1    :     4,  4       :   0,  0
#  1    :     4,  4       :   0,  1
#  1    :     4,  4       :   0,  2
#  2    :     2,  2       :   0,  0,  0
#  2    :     2,  2       :   0,  0,  1
#  1    :     4,  4       :   0,  3
#  2    :     2,  2       :   0,  1,  0
...省略...
#  2    :     2,  2       :   0,  3,  1
```
`omp_get_max_threads`にもスレッド数が反映される。注意点としては入れ子を抜けたあとも`omp_get_max_threads`の値が最初にレベル0でセットした値に上書きされているため、場合によっては自分でデフォルトの値に戻す必要がある。

少し前までは入れ子領域内のparallel構文でスレッド数を指定する方法ではうまく機能してくれなかった(指定したスレッド数が無視されデフォルトのスレッド数で並列される)が今回はその再現ができなかったのでgccかOpenMPのバージョンが上がって直ったのかもしれない。(そもそも勘違いだったかもしれない)

### 入れ子のレベル / OpenMP Active Level
再帰処理などを使い、入れ子がどこまで深くなるかわからないような場合には最も深いレベルを指定してあげることができる。そのレベルに達したら新たなスレッドは作られずに1スレッドで処理される。現在の入れ子のレベル(深さ)は、

```c
int omp_get_active_level()
```
で取得可能。実行環境でどのレベルまでネストできるかは

```c
int omp_get_max_active_levels()
```
で確認できる。デフォルトでは制限なし(2147483647レベル)でネストできるようになっているので最大レベルを制限するにはコード内で

```c
omp_set_max_active_levels(4);
```
を指定、またはプログラム実行時に環境変数でレベルの最大値を与える。

```sh
$ OMP_MAX_ACTIVE_LEVELS=4 ./a.out
```

下記のコードで説明。

```c:omp_nest_recursive.c  
#include <stdio.h>
#include <stdlib.h>
#include <omp.h>

void print_thread()
{
  int omp_active_level = omp_get_active_level();
  int omp_max_threads  = omp_get_max_threads();
  int omp_num_threads  = omp_get_num_threads();
  int omp_thread       = omp_get_thread_num();

  printf("# %2d    :   %3d,%3d       : %3d\n", omp_active_level,
	 omp_max_threads, omp_num_threads, omp_thread);
}

// 2スレッドづつ生成する再帰的な入れ子処理
void print_omp_nested(int active_level)
{
  if(active_level == omp_get_active_level())  return;

  print_thread();

  active_level=omp_get_active_level();

#pragma omp parallel num_threads(2)
  {
    print_omp_nested(active_level);
  }
}

int main(int argc,char** argv)
{
  //  omp_set_max_active_levels(4);

  printf("# omp_get_num_procs         = %d\n", omp_get_num_procs());
  printf("# omp_get_max_threads       = %d\n", omp_get_max_threads());
  printf("# omp_get_num_threads       = %d\n", omp_get_num_threads());
  printf("# omp_get_max_active_levels = %d\n", omp_get_max_active_levels());
  printf("# omp_get_thread_limit      = %d\n", omp_get_thread_limit());

  omp_set_nested(1);

  printf("# level : max,num_threads : omp_thread\n");

  print_omp_nested(-1);

  omp_set_nested(0);
  return EXIT_SUCCESS;
}
```
このコードをそのまま実行すると2^2147483647までスレッドを生成しようとして破綻する。環境変数を設定するか`omp_set_max_active_levels();`のコメントアウトを外すことによって`omp_get_max_active_levels()`に設定した値が入り、それ以上深いレベルの入れ子スレッドを生成しなくなる。

```sh
OMP_MAX_ACTIVE_LEVELS=2 ./a.out  # この場合は2レベルの入れ子 
# omp_get_max_active_levels = 2  # <- ここが適切な値になる 
# level : max,num_threads : omp_thread
#  0    :     8,  1       :   0
#  1    :     8,  2       :   0
#  1    :     8,  2       :   1
#  2    :     8,  2       :   0
#  2    :     8,  2       :   1
#  2    :     8,  2       :   0
#  2    :     8,  2       :   1
```
今は`omp_set_max_active_levels(4);`はコメントアウトしている。コメントアウトを外すとこちらの値が使われて、環境変数に関わらず4レベルの入れ子を作る。**実行時の環境変数よりコード内に書いた数値が優先される。**

### 同時実行最大スレッド数 / OpenMP Thread Limit

入れ子を使ったOpenMPでは`omp_get_max_threads()`で取得できるスレッド数にかかわらず際限なくスレッドが生成されるおそれがある
(たいていはOpenMPスタックサイズの関係上、ネストを繰り返すと落ちるが)。現在の同時実行の最大スレッド数は

```c
int omp_get_thread_limit()
```
で取得可能。デフォルトでは制限なし(2147483647レベルまで)で同時にスレッドが発行できるようになっている。最大スレッド数を制限するにはプログラム実行時に環境変数で与える。

```sh
$ OMP_THREAD_LIMIT=20 ./a.out
```
最大スレッド数を与えることにより、そのスレッド数を超えた新たな入れ子領域は生成されなくなる。プログラム内で

```c
omp_set_thread_limit(20);
```

を指定すると最大スレッド数を指定できそうだが、このような関数はないようだ。
再帰版の確認コードだとわかりにくいので以下の最初の`omp_nest2.c`を2スレッドづつ5レベル入れ子するコードで動作確認をする。

```sh
$ OMP_THREAD_LIMIT=4 ./a.out
# omp_get_num_procs         = 8
# omp_get_max_threads       = 8
# omp_get_num_threads       = 1
# omp_get_max_active_levels = 2147483647
# omp_get_thread_limit      = 4
# level : max,num_threads : omp_thread(l0,l1,l2,l3,l4,l5)
#  0    :     8,  1       :   0
#  1    :     8,  2       :   0,  0
#  1    :     8,  2       :   0,  1
#  2    :     8,  2       :   0,  0,  0
#  2    :     8,  1       :   0,  0,  0,  0
#  2    :     8,  2       :   0,  1,  0
#  2    :     8,  1       :   0,  1,  0,  0
#  2    :     8,  1       :   0,  1,  0,  0,  0
#  2    :     8,  1       :   0,  1,  0,  0,  0,  0
#  2    :     8,  2       :   0,  0,  1
#  2    :     8,  2       :   0,  1,  1
#  2    :     8,  1       :   0,  0,  0,  0,  0
#  2    :     8,  1       :   0,  0,  0,  0,  0,  0
#  2    :     8,  1       :   0,  0,  1,  0
#  2    :     8,  1       :   0,  0,  1,  0,  0
#  2    :     8,  1       :   0,  0,  1,  0,  0,  0
#  2    :     8,  1       :   0,  1,  1,  0
#  3    :     8,  2       :   0,  1,  1,  0,  0
#  3    :     8,  1       :   0,  1,  1,  0,  0,  0
#  3    :     8,  2       :   0,  1,  1,  0,  1
#  3    :     8,  1       :   0,  1,  1,  0,  1,  0

$ OMP_THREAD_LIMIT=8 ./a.out
# omp_get_num_procs         = 8
# omp_get_max_threads       = 8
# omp_get_num_threads       = 1
# omp_get_max_active_levels = 2147483647
# omp_get_thread_limit      = 8
# level : max,num_threads : omp_thread(l0,l1,l2,l3,l4,l5)
#  0    :     8,  1       :   0
#  1    :     8,  2       :   0,  0
#  1    :     8,  2       :   0,  1
#  2    :     8,  2       :   0,  0,  0
#  2    :     8,  2       :   0,  0,  1
#  2    :     8,  2       :   0,  1,  0
#  3    :     8,  2       :   0,  0,  0,  0
#  2    :     8,  2       :   0,  1,  1
#  3    :     8,  2       :   0,  0,  0,  1
#  3    :     8,  2       :   0,  0,  1,  0
#  2    :     8,  1       :   0,  1,  0,  0
#  2    :     8,  1       :   0,  1,  0,  0,  0
#  3    :     8,  2       :   0,  0,  1,  1
#  3    :     8,  1       :   0,  0,  1,  0,  0
#  3    :     8,  1       :   0,  0,  1,  0,  0,  0
#  2    :     8,  1       :   0,  1,  0,  0,  0,  0
#  3    :     8,  1       :   0,  0,  1,  1,  0
#  3    :     8,  1       :   0,  0,  1,  1,  0,  0
#  4    :     8,  2       :   0,  0,  0,  0,  0
#  5    :     8,  2       :   0,  0,  0,  0,  0,  0
#  2    :     8,  1       :   0,  1,  1,  0
#  2    :     8,  1       :   0,  1,  1,  0,  0
#  2    :     8,  1       :   0,  1,  1,  0,  0,  0
#  5    :     8,  2       :   0,  0,  0,  0,  0,  1
#  4    :     8,  2       :   0,  0,  0,  0,  1
#  4    :     8,  2       :   0,  0,  0,  1,  0
#  4    :     8,  1       :   0,  0,  0,  1,  0,  0
#  4    :     8,  2       :   0,  0,  0,  1,  1
#  4    :     8,  1       :   0,  0,  0,  1,  1,  0
#  5    :     8,  2       :   0,  0,  0,  0,  1,  0
#  5    :     8,  2       :   0,  0,  0,  0,  1,  1
```
使用する最大スレッド数に応じてなんとなく入れ子の深さも深くなることがわかるがどれくらい深くなるか、スレッド生成されるかは実行するたびに変化する。全てのスレッドが同期をとって並行で動いてるわけではなくバラバラで動くので速いスレッドが次の入れ子に入り、遅いスレッドは入れ子レベルが浅いところで処理を終えたりするため。うまく同期を取るようにすればもう少しわかり易い動作確認ができると思うが割愛する。

**`OMP_THREAD_LIMIT`, `OMP_MAX_ACTIVE_LEVELS` を同時に指定することも可能**

## まとめ
- OpenMPの入れ子機能の基本的な使い方の説明をした。
- 入れ子領域の度にスレッド生成するのでコストに注意。
- 入れ子の深さがわからない場合は`OMP_THREAD_LIMIT, OMP_MAX_ACTIVE_LEVELS (omp_set_max_active_levels())`で制御。
- nested lock機能などもあるのでまとめたら紹介しようと思う。[ここ](https://qiita.com/csms/items/c991db46591b8d8f27ad)にまとめた。