11. CUDA+MPI

11.1 CUDAとMPIの併用

CUDAとMPIを併用すると、複数のグラフィックスボードを持ったコンピュータ、 またはグラフィックスボードを持った複数台のコンピュータで並列計算することができます。
グラフィックスボードではCUDAを用いて高速に計算を行い、 プロセス間の通信はMPIを用いて行います。

CUDA+MPIのプログラミング方法
CUDA+MPIではMPIを用いて問題をプロセスに分割しCUDAで計算を行います。
host-device memoryでは以下のような手順になります。

  1. CPUで前処理を行います。
  2. CPUからGPUにメモリーをコピーします。
  3. GPUで計算を行います。
  4. GPUの計算の途中でプロセス間の通信が必要になったときはGPU間では通信ができないため、 GPUからCPUにメモリーコピーし、MPIを用いてCPU間で通信を行い、 結果をGPUにメモリーコピーします。
  5. 計算が終了したらGPUからCPUにメモリーコピーします。
  6. CPUで後処理を行います。
CUDA6から追加されたunified memory(managed memory)を使用すると以下のようになります。
  1. CPUで前処理を行います。
  2. GPUで計算を行います。
  3. GPUの計算の途中でプロセス間の通信が必要になったときはMPIを用いて通信を行います。
  4. CPUで後処理を行います。
host-device memoryの2.と5.を省略することができ、4.のメモリーコピーも不要になります。 これからわかるように、unified memoryを用いるとCUDA+MPIのプログラムが簡素化されます。
なお、CUDA対応のMPIを用いるとhost-device memoryの4.のメモリーコピーは不要になります。[14]
以下ではhost-device memoryとunified memoryをプログラムの引数で使い分けるプログラムを説明します。

CUDA+MPIの作業手順
CUDA+MPIは2段階の並列計算になりますので、以下の手順で開発することが必要です。 2.と3.は逆でも構いません。

  1. 逐次版を作成し十分動作確認を行います。
  2. CUDA版を作成し十分動作確認を行います。
  3. MPI版を作成し十分動作確認を行います。
  4. CUDA+MPI版を作成します。

11.2 CUDA+MPIプログラミング例

リスト11-1にベクトルの和をCUDA+MPIで並列計算するプログラムを示します。


リスト11-1 CUDA+MPIプログラム (vadd_cuda_mpi.cu)
     1	/*
     2	add two vectors (CUDA + MPI)
     3	
     4	Compile,Link:
     5	> nvcc -O3 -DMPI -o vadd_cuda_mpi vadd_cuda_mpi.cu cuda_memory.cu msmpi.lib
     6	
     7	Usage:
     8	> mpiexec -n <proc> vadd_cuda_mpi [-gpu|-cpu] [-hdm|-um] <n> <loop>
     9	> mpiexec -hosts <hosts> <host> <proc> [<host> <proc> ...] vadd_cuda_mpi [-gpu|-cpu] [-hdm|-um] <n> <loop>
    10	*/
    11	
    12	// GPU/CPU
    13	__host__ __device__
    14	static void vadd_calc(float a, float b, float *c)
    15	{
    16		*c = a + b;
    17	}
    18	
    19	// GPU
    20	__global__
    21	static void vadd_gpu(int n, const float *a, const float *b, float *c)
    22	{
    23		int tid = threadIdx.x + (blockIdx.x * blockDim.x);
    24		if (tid < n) {
    25			vadd_calc(a[tid], b[tid], &c[tid]);
    26		}
    27	}
    28	
    29	// CPU
    30	static void vadd_cpu(int n, const float *a, const float *b, float *c)
    31	{
    32		for (int i = 0; i < n; i++) {
    33			vadd_calc(a[i], b[i], &c[i]);
    34		}
    35	}
    36	
    37	// GPU/CPU
    38	static void vadd(int gpu, int n, const float *a, const float *b, float *c)
    39	{
    40		if (gpu) {
    41			int block = 256;
    42			int grid = (n + block - 1) / block;
    43			vadd_gpu<<<grid, block>>>(n, a, b, c);
    44		}
    45		else {
    46			vadd_cpu(n, a, b, c);
    47		}
    48	}
    49	
    50	#include <stdlib.h>
    51	#include <stdio.h>
    52	#include <string.h>
    53	#include <time.h>
    54	#ifdef MPI
    55	#include <mpi.h>
    56	#endif
    57	
    58	static int device_number(int, int, int []);
    59	extern void cuda_malloc(int, int, void **, size_t);
    60	extern void cuda_free(int, void *);
    61	extern void cuda_memcpy(int, void *, const void *, size_t, cudaMemcpyKind);
    62	
    63	int main(int argc, char **argv)
    64	{
    65		int    gpu = 1;
    66		int    um = 0;
    67		int    comm_size = 1;
    68		int    comm_rank = 0;
    69		int    n = 1000;
    70		int    nloop = 1000;
    71		int    nhost = 1;
    72		int    *ndevice;
    73		float  *a, *b, *c;
    74		clock_t t0 = 0, t1 = 0;
    75	
    76		// initialize (MPI)
    77	#ifdef MPI
    78		MPI_Init(&argc, &argv);
    79		MPI_Comm_size(MPI_COMM_WORLD, &comm_size);
    80		MPI_Comm_rank(MPI_COMM_WORLD, &comm_rank);
    81	#endif
    82	
    83		// arguments
    84		while (--argc) {
    85			argv++;
    86			if      (!strcmp(*argv, "-hosts")) {
    87				if (--argc) {
    88					nhost = atoi(*++argv);
    89					if (nhost < 1) nhost = 1;
    90					ndevice = (int *)malloc(nhost * sizeof(int));
    91					for (int ihost = 0; ihost < nhost; ihost++) {
    92						if (argc > 1) {
    93							ndevice[ihost] = atoi(*++argv);
    94							argc--;
    95						}
    96						else {
    97							ndevice[ihost] = 1;
    98						}
    99					}
   100				}
   101			}
   102			else if (!strcmp(*argv, "-gpu")) {
   103				gpu = 1;
   104			}
   105			else if (!strcmp(*argv, "-cpu")) {
   106				gpu = 0;
   107			}
   108			else if (!strcmp(*argv, "-hdm")) {
   109				um = 0;
   110			}
   111			else if (!strcmp(*argv, "-um")) {
   112				um = 1;
   113			}
   114			else if (argc == 2) {
   115				n = atoi(*argv);
   116			}
   117			else if (argc == 1) {
   118				nloop = atoi(*argv);
   119			}
   120		}
   121	
   122		// GPU device
   123		if (gpu) {
   124			// rank -> device number
   125			int device = device_number(comm_rank, nhost, ndevice);
   126	
   127			// GPU info
   128			cudaDeviceProp prop;
   129			cudaGetDeviceProperties(&prop, device);
   130			printf("Rank-%d GPU-%d : %s, C.C.%d.%d, U.M.%s\n", comm_rank, device, prop.name, prop.major, prop.minor, (um ? "ON" : "OFF"));
   131			fflush(stdout);
   132	
   133			// set device
   134			cudaSetDevice(device);
   135		}
   136	
   137		// global size : n
   138		// local size : l_n
   139		int l_n = (n + (comm_size - 1)) / comm_size;
   140	
   141		// alloc device memory
   142		size_t size = l_n * sizeof(float);
   143		cuda_malloc(gpu, um, (void **)&a, size);
   144		cuda_malloc(gpu, um, (void **)&b, size);
   145		cuda_malloc(gpu, um, (void **)&c, size);
   146	
   147		// alloc host memory
   148		float *h_a = (float *)malloc(size);
   149		float *h_b = (float *)malloc(size);
   150	
   151		// setup problem
   152		for (int i = 0; i < l_n; i++) {
   153			int gid = (comm_rank * l_n) + i;
   154			if (gid < n) {
   155				h_a[i] = gid;
   156				h_b[i] = gid + 1;
   157			}
   158		}
   159	
   160		// copy host to device
   161		cuda_memcpy(gpu, a, h_a, size, cudaMemcpyHostToDevice);
   162		cuda_memcpy(gpu, b, h_b, size, cudaMemcpyHostToDevice);
   163	
   164		// timer
   165	#ifdef MPI
   166		MPI_Barrier(MPI_COMM_WORLD);
   167	#endif
   168		if (comm_rank == 0) {
   169			t0 = clock();
   170		}
   171	
   172		// calculation
   173		for (int loop = 0; loop < nloop; loop++) {
   174			vadd(gpu, l_n, a, b, c);
   175		}
   176		if (gpu) cudaDeviceSynchronize();
   177	
   178		// timer
   179	#ifdef MPI
   180		MPI_Barrier(MPI_COMM_WORLD);
   181	#endif
   182		if (comm_rank == 0) {
   183			t1 = clock();
   184		}
   185	
   186		// copy device to host
   187		float *h_c = (float *)malloc(size);
   188		cuda_memcpy(gpu, h_c, c, size, cudaMemcpyDeviceToHost);
   189	
   190		// local sum
   191		double sum = 0;
   192		for (int i = 0; i < l_n; i++) {
   193			sum += h_c[i];
   194		}
   195	
   196		// reduction (MPI)
   197	#ifdef MPI
   198		double l_sum = sum;
   199		MPI_Reduce(&l_sum, &sum, 1, MPI_DOUBLE, MPI_SUM, 0, MPI_COMM_WORLD);
   200	#endif
   201	
   202		// output
   203		if (comm_rank == 0) {
   204			double exact = (double)n * n;
   205			double sec = (double)(t1 - t0) / CLOCKS_PER_SEC;
   206			printf("nproc=%d n=%d nloop=%d %e(%e) %s[sec]=%.3f\n",
   207				comm_size, n, nloop, sum, exact, (gpu ? "GPU" : "CPU"), sec);
   208			fflush(stdout);
   209		}
   210	
   211	#ifdef MPI
   212		MPI_Finalize();
   213	#endif
   214	
   215		// free
   216		free(h_a);
   217		free(h_b);
   218		free(h_c);
   219		cuda_free(gpu, a);
   220		cuda_free(gpu, b);
   221		cuda_free(gpu, c);
   222	
   223		return 0;
   224	}
   225	
   226	// rank -> device number
   227	static int device_number(int comm_rank, int nhost, int ndevice[])
   228	{
   229		int device = 0;
   230	
   231		if (nhost <= 1) {
   232			// single node
   233			device = comm_rank;
   234		}
   235		else {
   236			// cluster
   237			device = -1;
   238			int rank = -1;
   239			for (int ihost = 0; ihost < nhost; ihost++) {
   240				for (int idevice = 0; idevice < ndevice[ihost]; idevice++) {
   241					if (++rank == comm_rank) {
   242						device = idevice;
   243						break;
   244					}
   245				}
   246				if (device >= 0) {
   247					break;
   248				}
   249			}
   250			if (device < 0) device = 0;
   251		}
   252	
   253		int num_device;
   254		cudaGetDeviceCount(&num_device);
   255		if (device >= num_device) device = num_device - 1;
   256	
   257		return device;
   258	}

ソースコードの説明
86-101行目, 123-135行目, 227-258行目: MPIではプログラムの知っている情報はプロセス数と自分のランクのみです。 これでは自分の割り当てるべきデバイス番号がわかりませんので、 引数の"-hosts"オプションによってランクからデバイス番号が求められるようにしています。
139行目:全体の問題サイズnから各プロセスの問題サイズl_nを求めます。
142-145行目: 各プロセスは部分配列を持ちます。 これによってプロセス数を増やせばそれに比例して扱える問題サイズが大きくなります。
153-154行目: 各プロセスのランクと問題サイズから全体のインデックスを求めます。
全体の問題サイズはプロセス数の倍数とは限りませんので全体のインデックスが全体の問題サイズを超えないように条件判定が必要です。

コンパイル・リンク方法
コンパイル・リンク方法は以下の通りです。
> nvcc -O3 -DMPI -o vadd_cuda_mpi -arch=sm_30 vadd_cuda_mpi.cu cuda_memory.cu msmpi.lib
VC++では多数の"warning C4819"が出ることがあります。 そのときは以下のようにコンパイルオプションを追加してください。
> nvcc -O3 -DMPI -Xcompiler "/wd4819" -o vadd_cuda_mpi -arch=sm_30 vadd_cuda_mpi.cu cuda_memory.cu msmpi.lib

プログラムの実行方法
プログラムの実行方法は以下の通りです。
1台に複数のグラフィックスボードがあるとき:
> mpiexec -n プロセス数 vadd_cuda_mpi [-gpu|-cpu] [-hdm|-um] 配列の大きさ 繰り返し回数
複数ノードで計算するとき:
> mpiexec -hosts ホスト数 ホスト名 プロセス数 [[ホスト名 プロセス数] ...] vadd_cuda_mpi -hosts ホスト数 プロセス数 ... [-gpu|-cpu] [-hdm|-um] 配列の大きさ 繰り返し回数
mpiexecの引数の"ホスト名 プロセス数"はホストの数だけ必要です。
"プロセス数"には通常は各ノードのグラフィックスボードの数を指定します。 デバイス番号の小さい順に使用されます。
vadd_cuda_mpiの引数の"ホスト数 プロセス数 ..."にはmpiexecの引数の"ホスト数 ホスト名 プロセス数 ..."からホスト名を除いて数字だけを取り出したものを入力してください。
例えば以下のようになります。
> mpiexec -n 2 vadd_cuda_mpi 100000000 1000 (1台に2GPUがあるとき)
> mpiexec -n 2 vadd_cuda_mpi -um 100000000 1000 (unified memoryを使用するとき)
> mpiexec -n 2 vadd_cuda_mpi -cpu 100000000 1000 (CPUで計算するとき、デバッグ用)
> mpiexec -hosts 2 localhost 1 PC2 1 vadd_cuda_mpi -hosts 2 1 1 100000000 1000 (2台で計算するとき)
> mpiexec -hosts 2 localhost 2 PC2 1 vadd_cuda_mpi -hosts 2 2 1 100000000 1000 (2台で計算しlocalhostに2GPUがあるとき)
繰り返し回数は計算時間の測定誤差を小さくするためです。
mpiexecの使い方と実行環境の設定については7.5を参考にしてください。

11.3 CUDA+MPIの計算時間

表11-1に計算時間を示します。
1台のコンピュータに2個のGPUが実装されている場合です。
2GPUで計算すると1GPUのほぼ2倍速くなります。

表11-1 CUDA+MPIの計算時間(Windows, GTX 1070 X 2, ()内は1GPUとの速度比)
No.配列の大きさN繰り返し回数L1GPU2GPU
1100,000,0001,0005.78秒 (1.0)2.89秒 (2.00)
210,000,00010,0005.80秒 (1.0)2.92秒 (1.99)