MMS: 2次元定常熱伝導 — トラブルシューティングガイド
より充実した内容を mms-heat-2d.html でご覧いただけます。
「最初の1問」としての2次元定常熱伝導MMS
定常熱伝導のMMSは、コード検証の入門に最適な問題です。理由はソース項が手計算できるから——例えば単位正方形領域・熱伝導率 \( k \) 一定で、製造解に \( T_m = \sin(\pi x)\sin(\pi y) \) を選ぶと
$$ s = -k\,\nabla^2 T_m = 2\pi^2 k \,\sin(\pi x)\sin(\pi y) $$
と1行で出ます。境界は全周 \( T_m = 0 \)(Dirichlet)にでき、シンボリック計算すら不要。この「答え合わせが完全にできる1問」で検証の型(メッシュ系列・誤差ノルム・観測次数)を身につけてから、弾性・流体へ進むのが最短の学習路です。とはいえこの単純な問題でも詰まりどころは複数あり、本ページで症状別に解消します。
5ステップの標準手順
- ソース項 \( s(x,y) \) を体積熱源(単位W/m³:\( k \) の単位と整合させる)としてソルバーに定義
- 全周にDirichlet境界 \( T = 0 \)(一般には \( T_m \) の境界値)
- メッシュ3水準(例:10×10、20×20、40×40)で解く
- 各水準で誤差ノルム \( E_h = \|T_h - T_m\|_{L^2} \) を計算
- 観測次数 \( p = \ln(E_{2h}/E_h)/\ln 2 \) を要素の理論次数(1次要素で2、2次要素で3——\( L^2 \) ノルムの場合)と比較
症状1: 誤差がいつまでも大きい・明らかに別の解
| チェック | 説明 |
|---|---|
| ソース項の符号 | 方程式の書き方(\( -\nabla\cdot(k\nabla T) = s \) か \( \nabla\cdot(k\nabla T) + s = 0 \) か)で符号規約が変わる。逆符号だと解が上下反転する——温度コンターを見れば一目瞭然 |
| 係数の入れ忘れ | \( 2\pi^2 k \) の \( k \) を落とすと、解の形は同じで振幅が \( 1/k \) 倍ずれる。「形は合うが振幅が違う」はこのパターン |
| 境界の与え忘れ | FEMで無指定の辺は自動的に断熱(Neumannゼロ)。1辺でも指定漏れがあると解は大きく変わる。境界の割り当てを可視化で確認 |
| 単位系 | mm系モデルにW/m³でソースを入れると10⁹倍ずれる。手計算のオーダーと比較 |
症状2: 収束するが次数が理論値に届かない
- ソース項の空間分布の与え方 — ツールの体積熱源を「領域一定値」でしか与えられず、要素ごとの定数近似にしていると1次精度に律速される。座標依存の式・ユーザー関数・テーブル補間(十分細かい)で与える
- 誤差の測り方 — 節点値と厳密解を同位置で比較しているか、\( L^2 \) の体積重みを入れているか(詳細は次数算出のトラブルシュート)
- Neumann辺を混ぜた場合 — 熱流束 \( q = -k\,\partial T_m/\partial n \) の符号(外向き正か流入正か:ツール規約に合わせる)と、辺に沿って変化する値の積分精度
- 2次要素で誤差が「ゼロ」になる — 多項式の製造解(2次以下)は2次要素が厳密再現するため次数が測れない。三角関数系を使う(本ページの \( T_m \) なら問題なし)
症状3: 温度依存 k(T) にしたら合わなくなった
入門の次の一歩として \( k(T) \) を入れると、ソース項の導出が変わることを見落としがちです。正しくは
$$ s = -\nabla\cdot\left(k(T_m)\,\nabla T_m\right) = -k(T_m)\nabla^2 T_m - k'(T_m)\,|\nabla T_m|^2 $$
で、第2項(\( k' \) 項)の入れ忘れが定番ミスです。ここからはシンボリック導出(SymPy等)に切り替えるのが安全です。また非線形になるため反復(Picard/Newton)の収束を離散化誤差より2桁深くまで落とすこと——収束の浅い非線形解析は、線形で出ていた次数を簡単に壊します。
症状4: 細かいメッシュで誤差が下げ止まる
この単純な問題で最初に「誤差フロア」を体験することになります。原因は①線形ソルバーの収束判定(反復法の場合)、②丸め誤差(誤差が \( 10^{-10} \) 台に達した場合)、③ソース項・境界値の評価に混ぜた近似(テーブル補間の粗さ等)。診断は誤差の絶対値を見る+収束判定を1桁締めて変化を見る、の2手です。「フロアに当たるまで細分して、フロアの原因を特定して報告する」ところまでやると、検証レポートとして完成度が上がります。
症状5: 商用ツールで座標依存ソースが入れられない
MMSの実施可否はソルバー本体より「座標依存の体積熱源を定義できるか」で決まります。主要ツールの経路は:Ansys Mechanical=関数定義(Function/APDLの%table%・座標関数)、Abaqus=解析場(Analytical Field)またはDFLUXサブルーチン、COMSOL=式をそのまま入力(最も簡単)、OpenFOAM=codedSource/fvOptions、Fluent=UDF(DEFINE_SOURCE)。GUIから式入力できないツールでも、細かい空間テーブル+補間で代用できます(補間誤差が次数を汚さない細かさにする——症状2参照)。この「自分のツールでのMMS実施経路」を一度確立しておくと、以後のあらゆる検証で再利用できます。
チェックリストと発展
- ソース項の符号規約・係数・単位をソルバーの方程式表記と照合したか
- 全境界の割り当てを可視化で確認したか(無指定=断熱に注意)
- 誤差ノルムの測り方(同位置・体積重み)を確認したか
- 3水準以上で観測次数を計算し、理論次数と比較したか
- 誤差フロアの原因(収束判定/丸め)を特定したか
- 発展: Neumann/Robin境界を混ぜた版、k(T)非線形版、非定常版(時間次数の分離測定)へ段階的に拡張したか
正直、こんな簡単な問題をわざわざやる意味があるのかなと思っていました…。
この問題の真の成果物は「解けた」ではなく、検証の全手順を通した経験と、再利用できるスクリプト一式なんだ。メッシュ系列の作り方、誤差ノルムの計算、次数のプロット、フロアの見極め——ここで作った道具は、弾性でもNSでも、そのまま使い回せる。しかも問題が単純だから、手順のバグと物理のバグが混ざらない。楽器の音階練習と同じで、簡単な課題で「型」を作るから、難しい曲が弾けるようになる。1日で終わる投資としては、コード検証の世界で最も利回りが高い1問だよ。
関連記事:MMSの概要(統合版)、収束次数算出のトラブルシュート、2次元弾性MMS(次のステップ)。
関連トピック
なった
詳しく
報告