3rd Symposium

学際大規模情報基盤共同利用・共同研究拠点 公募型共同研究 平成２３年度採択課題

11-NA25

研究代表者氏名 (所属) 佐藤雅彦（核融合科学研究所）

研究課題名 陰的時間積分法による核融合プラズマの非線形MHDシミュレーションの高速化に関する研究
研究体制

目的
核融合炉の実現のためには、高温高密度プラズマを安定に閉じ込めることが必要であり、巨視的不安定性である磁気流体力学（MHD）不安定性の非線形現象の解明が重要である。
これまでに多くのMHD非線形シミュレーションが行われてきたが、今後は、・現実的なパラメータでの高精度MHD非線形シミュレーション・自由境界を取り扱ったMHD非線形シミュレーション
・二流体モデル等の拡張MHD非線形シミュレーション等を行っていくことが求められる。
このようなシミュレーションを陽解法で行うことは、クーラン条件の制約から非常に計算コストが大きくなる。
・クーラン条件を緩和するため、完全陰解法を採用したMHDコードの開発する。
・大規模連立一次方程式を高速に処理することが鍵。
→ Xabclibの使用、MHDモデルに対する適切な前処理の研究により、高速なMHDコードの開発を目指す。

MIPS code (MHD Infrastructure for Plasma Simulation)
核融合科学研究所で開発された基盤MHDコード
・一流体MHDモデル・4次精度の有限差分による離散化
・4次精度ルンゲクッタ法による時間積分・MPIによる並列化

佐藤雅彦(NIFS) : MIPSコードの陰解法化
片桐孝洋（東大）：Xabclibの開発
藤堂泰(NIFS)：陽的MIPSコードの開発、モデル拡張
中島徳嘉(NIFS)：MHD平衡コードの整備

プラズマシミュレータ

研究計画
1.
2.
3.
4.

Schur分解を利用したPhysics Based Preconditioningの導入
Schur分解を利用したStatic Condensation Method（静的縮約法）の導入
XabclibのMHDコードへの導入
拡張MHDモデルへの適用

Schur分解æ D
çç 1çççè L
æ Dçç 1çççè LU ö÷÷æç u ö÷ æçç r1 ö÷÷÷÷ç÷÷÷÷ = ççç÷÷çD2 ÷øè v ø çè r2 ÷÷ø(1) u* = D1-1r1(2) v* = r2 - Lu*(3) v = S-1v*ö÷æö÷
ö÷æ
-1U ö÷÷ æç I0D0IDUç
ç÷÷÷
÷
11
÷÷ = çç÷÷çç
÷÷çç-1÷I ÷øçè 0 S ÷øççè 0D2 ÷÷ø ççè LD1I ÷ø(4) u = u* - D1-1Uv

-1
LDU
2
1

ここで、 S = D¶r= -·(rv) + D2(r - req )¶t
æ v 2 ö÷4¶v= -rw ´ v - r ççç ÷÷ - P + j´ v + [nr(·v)]-  ´[nrw]r
çè 2 ÷ø3¶t¶B= - ´ E¶t4¶p= -·(pv) -(g -1)p·v + (g -1)[nrw 2 + nr(·v)2 + h j·(j - jeq )]3¶t
E = -v´ B + h(j - jeq )
1j=
´B
m0

「物理変数」で分割し、Schur分解を適用「 領域 」で分割し、Schur分解を適用

Static condensation method

Physics-based preconditioning
理想MHDモデルにおけるヤコビアン行列

æ D
çç r
çç
çç 0¶F·du = çççç 0¶uçççç Lrv
è00

DP00

DBLPv

LBv

Uvr ö÷÷æ÷÷ç dr ö÷÷çUvP ÷÷÷çç dP ÷÷÷÷÷çç÷÷÷UvB ÷÷çç dB ÷÷÷÷÷çç÷÷dv
èø÷Dv ÷ø双曲型の性質をもたらす1階の微分演算子 v·, B· は、ヤコビアン行列の中のL, Uに含まれる。
近似的なSchur補行列の作成により放物型問題へ変換

w = ´v

Physics‐based preconditioning
Static condensation method-1 *-11-1v = S u = (D2 - LD U) uæ Lçç IIçç
çè LGI

LI G ö÷÷æçç u I ö÷÷ æçç rI ö÷÷÷÷ç
÷÷ = ç÷÷çç
LGG ÷÷øçè u G ÷÷ø çè rG ÷÷ø
I : Interior
G : Interface: プロセスごとに独立に計算可能
L-1
II

S-1 : グローバル行列だが、オリジナル行列よりも次元が小さい。

*(D1-1 ~ 1 / Dt)(1 / Dt)LU = (1 / Dt)(Lrv Uvr + LPv UvP + LBv UvB )拡散型演算子,対角優位行列LとUの積をとることは、1階の微分演算子を2回作用させることであり、2階の微分演算子を作用させることと同一。
→ 双曲型問題を放物型問題へ変換·U(v·) · (v·)(v·)(v·)(B·) · (B·)(B·)(B·)L

行列LUは対角優位な行列であるため、クリロフ法の収束を早めることが期待できる。

陰解法の必要性
アルフヴェン波〜10‐6秒 << 抵抗性MHD不安定性〜10‐3秒
陽的な時間積分法では、時間刻み幅がアルフヴェン波の伝播時間に制限される。
MHDモデルから二流体モデルへ拡張した場合、さらに速い時間スケールの現象が現れ、クーラン条件による制約はより厳しくなる。
完全陰解法より、大きな時間刻み幅が可能に¶u= f(u)¶tun +1 - un
= f(un +1 )DtF(u) = f(uNewton法

Xabclib

n +1un +1 - un)=0DtF：ヤコビアン行列
u
次元〜108
1行あたり非零要素数=13〜54高速な線形ソルバーが必要

4次の有限差分における分割方法

自動チューニング機能付き数値計算ライブラリパッケージ

・汎用的自動チューニングAPIの提供
・自動チューニング機能の再利用性を考慮した、汎用的APIのOpenATLibを提供
・世界に先駆けた先進的AT機能の提供
1. 自動チューニングポリシ設定機能
実行速度、メモリ量、演算精度の方針設定機能
2. マルチコア向きの負荷バランスの改善が可能な疎行列演算方式
非零要素均等化方式,Branchless Segmented Scan (BSS) 方式
3. 直交化アルゴリズム動的選択機能
関数名

完全陰解法

n +1

LU

OpenATI_INIT

機能の説明
OpenATLib と Xabclibのパラメタを設定する

OpenATI_DAFRT

Krylov部分空間のリスタート周期を増加するか判定

OpenATI_DSRMVCRS形式において、倍精度演算 対称行列用の疎行列－ベクトル積で、最もよい実装方式の選択
CRS形式において、倍精度演算 非対称行列用の疎行列－ベクトル積で、最もよい実装方式の選択

OpenATI_DURMVOpenATI_DSRMV_Setup

OpenATI_DSRMVのための初期データ設定処理

OpenATI_DURMV_Setup

OpenATI_DURMVのための初期データ設定処理

OpenATI_DAFGS

Gram‐Schmidt直交化処理において、４種類の実装から最もよい実装方式の選択
数値計算ポリシを適用する連立一次方程式ソルバ用
メタ・インタフェース
数値計算ポリシを適用する固有値ソルバ用
メタ・インタフェース

OpenATI_LINEARSOLVEOpenATI_EIGENSOLVE

学際大規模情報基盤共同利用・共同研究拠点 第３回シンポジウム
2011年 7月 14日，15日

THE GRAND HALL （品川）