3.2 OpenMPによる並列化

3.2.1 OpenMPによる並列化

リスト3-2-1にSOR法のCPU版の反復計算部分を示します。
8行目でOpenMPを用いてループiを並列化しています。 またreduction句を用いて残差の和を求めています[6]。
14行目でred-black法の判定を行っています。

リスト3-2-1 OpenMPを用いた反復計算部の並列化（update.cの一部, CPU版, novectorモード）

     1	real_t update_no_vector(int oe)
     2	{
     3		const real_t omega = (real_t)Solver.omega;
     4	
     5		real_t res2 = 0;
     6		int i;
     7	#ifdef _OPENMP
     8	#pragma omp parallel for reduction (+: res2)
     9	#endif
    10		for (    i = iMin; i <= iMax; i++) {
    11		for (int j = jMin; j <= jMax; j++) {
    12		for (int k = kMin; k <= kMax; k++) {
    13			const int64_t n = D(i, j, k);
    14			if (((i + j + k) % 2 == oe) && (idVolt[n] == 0)) {
    15				const real_t e = Epsr[idEpsr[n]];
    16	
    17				const real_t axp = (e + Epsr[idEpsr[n + Ni]]) * RXp[i];
    18				const real_t axm = (e + Epsr[idEpsr[n - Ni]]) * RXm[i];
    19				const real_t ayp = (e + Epsr[idEpsr[n + Nj]]) * RYp[j];
    20				const real_t aym = (e + Epsr[idEpsr[n - Nj]]) * RYm[j];
    21				const real_t azp = (e + Epsr[idEpsr[n + Nk]]) * RZp[k];
    22				const real_t azm = (e + Epsr[idEpsr[n - Nk]]) * RZm[k];
    23				const real_t asum = axp + axm + ayp + aym + azp + azm;
    24	
    25				const real_t res = omega * (
    26					((axp * V[n + Ni]) + (axm * V[n - Ni]) +
    27					 (ayp * V[n + Nj]) + (aym * V[n - Nj]) +
    28					 (azp * V[n + Nk]) + (azm * V[n - Nk])) / asum - V[n]);
    29	
    30				V[n] += res;
    31				if (((i == 0) || (i > iMin)) &&
    32				    ((j == 0) || (j > jMin)) &&
    33				    ((k == 0) || (k > kMin))) {  // MPI
    34					res2 += res * res;
    35				}
    36			}
    37		}
    38		}
    39		}
    40	
    41		return res2;
    42	}

3.2.2 OpenMPの計算時間

表3-2-1にOpenMPの計算時間を示します。
vectorモードはスレッド数が少ないときは速いですが、 スレッド数が多くなると遅くなります。
これからCPUでは、使用メモリーが少なく、多スレッド時に速いnovectorモードを推奨します。
また、ハイパースレッディングの効果はほとんどないので物理コア数を推奨します。

表3-2-2に各ベンチマークの計算時間を示します。