MOOSEフレームワーク
理論:残差で物理を書く仕組み
概要
MOOSEは、ほかのオープンソースの有限要素ソフトと何が違うのですか?
「解析ソフト」ではなく「解析ソフトを作るための枠組み」である点だ。MOOSE(Multiphysics Object Oriented Simulation Environment)は米国アイダホ国立研究所が開発していて、メッシュ処理はlibMesh、線形・非線形ソルバーはPETScに任せ、利用者は支配方程式の「残差」だけを書く。熱伝導なら拡散項と熱源項を別々のKernelとして書き、入力ファイルで並べると、MOOSEがそれらを足し合わせて解く。物理を足すときもKernelを追加するだけで、熱・構造・化学反応などすべての未知数を1つの非線形方程式として連立させ、Newton法でまとめて解く。原子炉燃料の挙動を解くBISONをはじめ、多くの研究用ソルバーがこの上に作られている。
残差とNewton法
ψᵢ は試験関数、k(u) は解に依存する係数、f は源項、J はヤコビアン。MOOSEはJを近似的に作るか、行列を作らずに残差の差分で J·v を計算するJFNK(Jacobian-Free Newton-Krylov)法で解く。
製造解法(MMS)による検証
解をあらかじめ決め、それを満たす源項 f を逆算して与える。要素数 N を2倍にしたときの誤差の比から収束次数 p を求め、理論値(1次要素なら関数値で2、勾配で1)と比べる。
原子炉燃料のために生まれた枠組み
MOOSEの開発は2008年ごろ、アイダホ国立研究所で原子炉燃料の挙動を解くために始まった。燃料ペレットの中では熱、変形、核分裂ガスの放出、材料の変化が互いに影響し合い、別々のコードを順番に回す方法では強い結合を扱いきれなかった。そこで「物理は残差として書き、連立とソルバーは枠組みが受け持つ」という設計が選ばれ、2014年にオープンソースとして公開された。
計算例
例1:製造解 sin(πx) で収束次数を確かめる
区間 [0, 1] で −u″ = f を1次要素で解く。両端は u = 0、厳密解は sin(πx)。
| 要素数 N | L2誤差 | 勾配の誤差 | L2の次数 | 勾配の次数 |
|---|---|---|---|---|
| 4 | 3.93×10⁻² | 0.499 | — | — |
| 8 | 9.92×10⁻³ | 0.251 | 1.99 | 0.99 |
| 16 | 2.49×10⁻³ | 0.126 | 2.00 | 1.00 |
| 32 | 6.22×10⁻⁴ | 0.0630 | 2.00 | 1.00 |
要素を半分にするたびにL2誤差は4分の1、勾配の誤差は2分の1になり、理論どおりの次数が出ている。
節点の値を厳密解と比べたら、4要素でも誤差が10⁻⁷以下でした。もう十分正確ということですか?
それは1次元の1次要素に特有の性質で、だまされやすい。この問題では節点の値だけは厳密解とほぼ一致する(超収束)が、節点の間の値や勾配は、表のとおり4要素ではL2誤差4×10⁻²、勾配の誤差0.5もある。節点値だけを見ていると、要素の実装が間違っていても気づけないことがある。MOOSEでは、誤差をL2ノルムやH1半ノルムのような積分量で評価し、要素を2倍、4倍と細かくしたときの次数で判定するのが基本だ。次数が理論値から外れたら、Kernelの符号や積分点の扱い、境界条件の与え方を疑う。
例2:非線形拡散でNewton法とPicard反復を比べる
−(k(u)u′)′ = 10、k = 1 + u²、両端 u = 0、32要素。残差のノルムが10⁻¹⁰を下回るまでの反復回数を比べる。
| 解き方 | 反復回数 | 最大値 u |
|---|---|---|
| Newton法(ヤコビアンを正しく組む) | 5 | 0.957 |
| Picard反復(k を前回の値で固定) | 12 | 0.957 |
| 参考:k = 1(線形) | 1 | 1.250 |
Newton法は収束の終盤で残差が 7×10⁻³ → 4×10⁻⁶ のように2乗で減るのに対し、Picard反復は1回ごとにほぼ10分の1ずつしか減らない。k が温度とともに増えるため、最大値は線形の場合の1.25から0.957に下がる。
Newton法がいつも有利なら、なぜ分離して順番に解く方法も使われるのですか?
Newton法の弱点は、ヤコビアンを正しく作る手間と、初期値が解から遠いと発散しやすいことだ。物理ごとに既存のコードを使い回したい場合や、結合が弱く数回のやり取りで収束する場合は、分離型のほうが簡単で速い。MOOSEはこの両方に対応していて、強く結合した変数は1つのNewton系にまとめ、弱い結合や時間スケールの違う物理はMultiAppという仕組みで別々に解いてデータをやり取りできる。ヤコビアンが不正確だと二次収束が失われ、反復回数がPicard並みに増えるので、MOOSEには差分で作ったヤコビアンと比べてKernelの微分の誤りを見つける機能もある。
始め方の流れ
- condaなどで環境を作り、MOOSEとテストを一度通す。
- 付属の熱伝導などの例題を入力ファイルで動かし、Kernel・境界条件・Executionerの書き方を覚える。
- 自分のKernelを書いたら、製造解で収束次数を確かめる。
- 非線形問題ではNewton法の反復回数を見て、ヤコビアンの正しさを確かめる。
- 確かめた計算をテストとして登録し、変更のたびに自動で回す。
「収束次数が1.5で止まった」
新しい境界条件をKernelとして実装し、製造解で確かめると、L2の次数が2ではなく1.5前後で止まる、ということがよく起きる。原因の多くは、境界の項を要素の内部と同じ積分点で評価していない、あるいは源項を節点値で近似している、といった小さな不整合だ。結果の図だけ見ればもっともらしく見えるので、こうした誤りは次数を測らないと見つからない。MOOSEの開発で検証試験が重視されるのは、このためだ。
よくある間違い
間違いと対策
| 間違い | 起きること | 対策 |
|---|---|---|
| 節点値だけで検証 | 実装の誤りを見逃す | L2・H1ノルムで次数を見る |
| ヤコビアンの誤り | Newton法の回数が増える | 差分ヤコビアンと比較 |
| Kernelの符号の誤り | 発散や逆向きの解 | 製造解で確認 |
| 初期値が解から遠い | Newton法が発散 | 荷重増分・時間発展で近づける |
| 弱い結合まで1つの系に | 計算が重い | MultiAppで分離 |
関連する内容も知りたいです。
関連するページはMMS(製造解法)の概要と手順、MMS収束次数の算出方法、FEniCSプロジェクト、Elmerマルチフィジックス、Code_Aster入門が参考になる。
関連トピック
なった
詳しく
報告