∑ 流体解析における新しい並列解法の提案

流体解析における新しい並列解法の提案

東京大学大学工学部     学生員  Pham Van Phuc 東京大学大学院工学系研究科 正会員  石原孟

東京大学大学院工学系研究科 フェロー 藤野陽三 1. はじめに

キーワード: 並列計算、領域分割、並列解法、残差切除法

連絡先  : 〒113-8656 東京都文京区本郷7-3-1 tel. 03-5841-6099; fax. 03-5841-7454 風による構造物のギャロッピング振動を予測する際

には空気力係数が必要となる。しかし、空気力係数の計 算には3次元非定常流体解析が必要なため、多くの計算 時間がかかり、並列計算が不可欠である。しかし、並列 計算では多数のコンピュータを用いれば、必ず速くなる とは限らない。

本研究では並列計算を高速するために従来の並列解 法の問題点を明らかにし、線形連立方程式における高速 かつ安定な新しい並列解法を提案する。

2. 領域分割

複数台のコンピュータを用いて並列計算を行う場合 には、領域分割手法がよく用いられる¹⁾。すなわち、全 体の計算領域をいくつかの小さな領域に分割し、それぞ れの領域を近隣領域の境界値を境界条件として独立に 解く。分散メモリ型の並列計算機を用いる場合にマシン 間に通信が必要となる。

図 1 には本研究で開発した並列計算コードを用いて 計算された2次元キャビティ内の流れ場を示す。8つの 領域間の流れ関数コンターはなめらかにつながってい る。図２には1つ領域と8つの領域に分割されたときの 計算結果の比較を示し、その差は0.006％となり、よく 一致した解が得られることが分かる。

‑0.6

‑0.4

‑0.2 0 0.2 0.4

0 0.25 0.5 0.75 1

1Domain 8Domains

Vy

X

図-2  1領域と8領域で計算されたAB断面での速度 成分Vyの比較

3. 線形連立方程式の並列解法

流体解析においては、Poisson方程式を解く部分が計 算時間の大部分を占める。従来の並列解法の問題点を明 らかにするために式（1）に示すようなポワソン方程式 を例として用いた。境界条件はDirichlet境界条件とし て、計算格子は256×256とした。

(

) (

y

x i

iπ π

φ

φ ₂ sin sin

2 2 2



 



=

∂ +∂

∂

∂

∑ ⁾

20 , 8 , 2 ,

=1 ki

3.1 従来の並列解法

並 列 数 値 解 法 と し て は 非 構 造 格 子 用 の GS  (Guass-Seidel)¹⁾,構造格子用の SIP(Strongly Implicit

Procedure)¹⁾を用いて、シリアル計算（単一領域）と領

域分割による並列計算を行った。0.001残差でシリアル 計算ではSIP法はGS法より10倍も速い。

図3には領域分割よる並列計算の収束状況を示し、シ リアル計算による収束までの反復回数で無次化した。領 域の数が増えてもGS法の反復回数はほとんど変わらな い。一方、SIP法の反復回数は領域の数が増えるにつれ、

大幅に増大する。以上のことから、シリアル計算での速 い数値解法はそのまま並列計算に適用することができ ないことが分かる。この理由としてGS法では、データ の依存関係がローカルのため、解析速度は領域分割の影 Re=1000

160x160CVs 誤差：0.0001

図-1 二次元キャビティ内の流れ関数のコンター

土木学会第58回年次学術講演会（平成15年9月）

‑253‑

I‑127

響を殆ど受けない。一方、SIP法ではデータの依存関係 がグローバルのため、領域分割によりその関係が破壊さ れてしまう。従って、領域の数が増えると、反復回数が 大きくなると考えられる。

0

3.2 新しい並列解法

図４には本研究で提案した新しい並列解法のフロ-チ ャートを示す。まずGS法あるいはSIP法などを用いて、

ローカル通信（LC）により各領域における近似解を求 める。そして、グロ

ーバルな通信（GC）

に よ り 各 マ シ ン の 残差を1台マシンに 集め、最小２乗法に よ り 全 体 の 誤 差 を 最 小 と な る よ う に 係数を決定する。最 後に、これらの係数 を グ ロ ー バ ル な 通 信 に よ り 各 マ シ ン に配送し、新しい推 測値を求める。この よ う な プ ロ セ ス を 収 束 条 件 を 満 足 す るまでに繰り返す。

新 し い 並 列 解 法 で は 従 来 の 並 列 解

法のようにローカルな残差を消去するだけではなく、残 差切除法（RCM：Residual Cutting Method）を組み込 むことによりグローバルな残差も消去できる²⁾。これに より、高速かつ安定な並列解法を実現した。

図５には新しい並列解法による収束状況を示す。SIP 法とRCM法とを組合せた場合には領域の数が増えても

反復の回数がほぼ一定となった。

0.0001 0.001 0.01 0.1 1

0 0.5 1 1.5

1Domain 2Domains 4Domains 8Domains 16Domains

2 0.0001

0.001 0.01 0.1 1

0 0.5 1 1.5 2

3.3 加速率の比較

図６には各種の並列解法によるCPU数と加速率の関 係を示す。新しい並列解法を用いた場合に、計算の加速 率とCPU数との関係がほぼ直線となる。すなわち、CPU の数が多くなると、計算が高速になることが分かる。

4. まとめ

本研究では従来の並列解法の問題点を明らかにし、グ ローバル残差を消去できる残差切除法を組み込むこと により、高速かつ安定な並列解法を実現した。

参考文献

１ ） J.H.Ferziger and M.Peric: Computational Methods for Fluid Dynamics. Springer, 2002。

２）石原孟、山口敦、藤野陽三：複雑地形における局所 風況の数値予測と大型風洞実験による検証、土木学会論 文集、2003.4（掲載予定）。

図-3 従来の並列解法による収束の様子

0.001 0.01 0.1 1

0 0.5 1 1.5 2

相対残差のノルム

SIPRCM

a) b)

.001 0.01 0.1 1

0 0.5 1 1.5 2

1Domain 2Domains 4Domains 8Domains 16Domains

相対残差のノルム

GSRCM GS SIP

無次元化反復回数 無次元化反復回数 b)

a) 図-5 新しい並列解法による収束の様子

無次元化反復回数 無次元化反復回数

2 4 6 8 10

12 GS

SIP GSRCM SIPRCM

1 2 4 8 16

1

12.1 12.5

GS,SIP法等により近似解 を求め：LC 残差を収集：GC

全体の誤差が最小になる ように係数を決定

係数を配送：GC 新しい推測値を求め

収束条件

終了 開始

NO

YES

9.5

加速率

5.6

CPU数

図-6各種の並列解法によるCPU数と加速率の関係

図-4 新しい並列解法の流れ

土木学会第58回年次学術講演会（平成15年9月）

‑254‑

I‑127