CFDの時間積分法
理論：陽解法と陰解法
概要
円柱の後ろのカルマン渦をCFDで解いたら、渦がすぐに消えて定常の流れになってしまいました。何が悪いのでしょう？
時間積分法を確かめてみよう。多くのソフトの初期設定は1次の陰解法（後退オイラー法）で、とても安定な代わりに、振動を少しずつ弱める『数値減衰』が強い。大きな時間刻みと組み合わせると、本来続くはずの渦の放出まで弱めてしまうことがある。陽解法は時間刻みがCFL条件で厳しく制限され、陰解法は大きな刻みでも安定だけど、振動を正しく保てるかは方法と刻みで決まる。簡単な振動のモデルで比べると、その違いがはっきり見えるよ。
1ステップの増幅率
振動の式に各方法を当てはめたときの、1ステップあたりの増幅率 $G$。$|G| = 1$ なら振幅が保たれ、$|G| < 1$ なら減衰、$|G| > 1$ なら発散する。
時間刻みの制限
陽解法では、1ステップで流体がセルを越えない（CFL ≤ 1程度）ことと、拡散の制限が必要になる。
天気予報から始まった時間の刻み
1922年、リチャードソンは手計算で天気を予報しようと試み、6時間の予報に何週間もかかったうえに、結果は大きく外れた。原因の一つは、時間刻みと格子の関係が計算の安定性を決めることが知られていなかったことだった。1928年、クーラン、フリードリッヒス、レヴィが、陽解法では時間刻みを格子と波の速さで決まる値以下にしないと解が発散することを示した（CFL条件）。この条件を避けるために陰解法が発達し、大きな時間刻みの代わりに数値減衰という別の問題と向き合うことになった。
計算例：振動のモデルで比べる
減衰のない振動（振幅は本来一定）を、1周期あたりのステップ数を変えて解き、10周期後に残る振幅と、1周期あたりの位相の誤差を比べた。
| 1周期のステップ数 | 陽解法（前進オイラー） | 陰解法（後退オイラー） | BDF2 | クランク–ニコルソン |
|---|---|---|---|---|
| 10 | 発散 | 6×10⁻⁸ | 0.14 | 1.00 |
| 20 | 発散 | 8×10⁻⁵ | 0.68 | 1.00 |
| 40 | 発散（131倍） | 0.0076 | 0.944 | 1.00 |
| 80 | 発散（11.7倍） | 0.086 | 0.993 | 1.00 |
位相の誤差（1周期あたり）は、40ステップで後退オイラー法とBDF2が2.9°、クランク–ニコルソン法が0.7°、80ステップではそれぞれ0.7°と0.2°。
1次の陰解法だと、40ステップでも10周期で振幅が1%以下になるんですね。
そう。1次の陰解法は1ステップごとに振幅を少しずつ削るから、長く計算するほど効いてくる。実際のカルマン渦は流れの不安定さで自分から振動を生み出すので、完全には消えないこともあるけれど、強さや周波数が変わってしまう。BDF2は2次精度で減衰がずっと小さく、実務の非定常計算の標準になっている。クランク–ニコルソン法は振幅を保つけれど、急な変化では振動が残りやすい。例えば直径0.1 mの円柱に10 m/sの流れなら渦の周波数は約20 Hzで、1周期40ステップなら時間刻みは1.25 ms。このときセル1 mmならCFLは12.5になり、陽解法では解けないから陰解法を使うことになるんだ。
選び方の目安
- 解きたい最も速い変化（渦の放出、脈動など）の周期を見積もる。
- 非定常の結果を見るなら、1次の陰解法ではなくBDF2などの2次精度の方法を使う。
- 1周期に20〜40ステップ以上を目安に時間刻みを決め、CFLも確かめる。
- 時間刻みを半分にして、周波数や振幅が変わらないことを確かめる。
- 1ステップ内の反復（外側の反復）が十分に収束しているかも確かめる。
「渦の周波数は合うが大きさが半分」
橋の断面の空力振動を調べるCFDで、渦の放出の周波数は風洞実験と合ったのに、揚力の変動の大きさが半分ほどしかなかった——という例がある。時間積分は1次の陰解法で、時間刻みは1周期あたり15ステップ程度だった。BDF2に変えて1周期40ステップにすると、変動の大きさは実験の1割以内に入った。周波数は位相の誤差、大きさは振幅の減衰で決まるので、周波数が合っているからといって大きさも正しいとは限らない。
よくある間違い
間違いと対策
| 間違い | 起きること | 対策 |
|---|---|---|
| 非定常で1次の陰解法 | 振動が消える・弱まる | BDF2などの2次 |
| 時間刻みが粗い | 位相・振幅の誤差 | 1周期20〜40ステップ以上 |
| 陽解法でCFLを無視 | 発散 | CFLを確認 |
| 周波数だけで検証 | 振幅の誤差を見落とす | 振幅も比べる |
| ステップ内の反復が不十分 | 時間精度が落ちる | 残差を確認 |
関連する内容も知りたいです。
関連トピック
なった
詳しく
報告