オイラー方程式(圧縮性)
オイラー方程式(圧縮性)の理論基礎
概要
先生、圧縮性のオイラー方程式ってどういうものですか? 非圧縮のオイラー方程式とは別物なんですか?
オイラー方程式は粘性を無視した流体の運動方程式だけど、圧縮性の場合は密度が変化するから、エネルギー方程式も含めた連立保存則の形になる。これが圧縮性CFDの出発点なんだ。
保存形オイラー方程式
1次元の場合、保存変数ベクトル $\mathbf{U}$ とフラックス $\mathbf{F}$ を使ってこう書ける。
ここで
3つの式が連立してるんですね。それぞれ質量保存、運動量保存、エネルギー保存ですか?
そのとおり。3次元では運動量方程式が3成分になるから5元連立系になる。
これを閉じるために状態方程式が必要で、理想気体なら
を使う。空気なら比熱比 $\gamma = 1.4$ だ。
双曲型方程式としての性質
オイラー方程式が双曲型って聞いたことがあるんですけど、それはどういう意味ですか?
重要な概念だね。1次元オイラー方程式のヤコビアン $\mathbf{A} = \partial\mathbf{F}/\partial\mathbf{U}$ の固有値が全て実数であることが双曲型の条件だ。固有値は
ここで $a = \sqrt{\gamma p/\rho}$ は音速。これらは特性速度と呼ばれ、情報が伝播する速度を表す。
$\lambda_1$ と $\lambda_3$ は音波の前進波・後退波、$\lambda_2$ はエントロピー波(接触不連続面)に対応する。超音速流れ($|u| > a$)では全固有値が同符号になるから、情報は一方向にしか伝わらない。これがCFDの境界条件設定に直結するんだ。
超音速だと下流の情報が上流に伝わらないってことですね。境界条件の設定がマッハ数によって変わる理由がやっと分かりました。
亜音速入口では圧力などの情報が上流に伝播するから、1つの物理量は外部から指定し残りは内部から外挿する。超音速入口では全ての特性が流入方向だから、全変数を指定する。この原則を守らないと非物理的な反射波が発生するんだ。
18世紀の天才が書いた式が今もCFDの心臓部に
レオンハルト・オイラーが流体の運動方程式を導いたのは1757年のこと。粘性も熱伝導もない「理想流体」を仮定したこの方程式は、当時は抽象的すぎると批判された。ところが270年後、超音速・極超音速のCFDではこのオイラー方程式が粘性の支配するナビエ・ストークス方程式よりも先に使われる。なぜなら粘性の効果が慣性力に比べて極めて小さくなる高速流では、オイラー方程式が現実をよく近似するからだ。時代を超えて生き続ける数式の強さを感じる瞬間だ。
オイラー方程式(圧縮性)の数値計算手法
数値解法
オイラー方程式をコンピュータで解くには、衝撃波の扱いがポイントになるんですよね?
そのとおり。オイラー方程式の解には衝撃波や接触不連続面といった不連続解が含まれる。これらを数値的に正しく捕獲するための手法が、圧縮性CFDの核心技術だ。
Godunovの方法とリーマンソルバー
Godunovの方法では、セル界面で局所的なリーマン問題を解いて数値フラックスを求める。
ここで $\mathbf{U}^*$ はリーマン問題の解、$\mathbf{U}_L, \mathbf{U}_R$ はセル界面の左右の状態だ。厳密リーマンソルバーは計算コストが高いので、近似リーマンソルバーが広く使われている。
| ソルバー | 特徴 | 長所 | 短所 |
|---|---|---|---|
| Roe (1981) | 線形化リーマン解 | 接触不連続の解像度高い | エントロピー修正が必要 |
| HLL (1983) | 2波近似 | 堅牢、実装が簡単 | 接触不連続が拡散 |
| HLLC (1994) | 3波近似(接触波含む) | 接触不連続も解像 | Roeより若干拡散 |
| AUSM+ (1996) | 質量流束分離 | 全マッハ数対応 | パラメータ調整必要 |
Roeのリーマンソルバーが有名ですけど、具体的にはどう計算するんですか?
Roeは左右の状態からRoe平均を計算して、ヤコビアンを線形化する。Roe平均密度は $\hat{\rho} = \sqrt{\rho_L \rho_R}$ で定義される。数値フラックスは
第2項が数値的な散逸を加える風上項だ。
高次精度化: MUSCL法とリミッター
1次精度のGodunovスキームでは衝撃波が数十セルにまたがってしまう。高次精度化にはMUSCL(Monotone Upstream-centered Schemes for Conservation Laws)法を使う。
ここで $\phi(r)$ はスロープリミッター関数。リミッターを使わないと衝撃波近傍でGibbs振動が発生する。代表的なリミッターには minmod、van Leer、superbee がある。
リミッターの選び方で結果が変わるんですか?
かなり変わる。minmod は最も拡散的だけど安定、superbee は鮮鋭だけど疎密波を階段状にしてしまう傾向がある。実用的には van Leer か van Albada が良いバランスだ。WENO(Weighted Essentially Non-Oscillatory)スキームは5次精度で衝撃波捕獲も優れているが、計算コストが高い。
時間積分はどうするんですか?
陽解法ならTVD Runge-Kutta法(Shu-Osher 3段3次)が定番だ。CFL条件は
CFL $\leq 1$ が安定条件だ。陰解法ではLU-SGS(Lower-Upper Symmetric Gauss-Seidel)法が効率的で、CFL数を大きく取れるから定常計算の収束が速いんだ。
Godunov法の誕生——1959年のソ連が生んだ衝撃波数値解法
オイラー方程式を数値的に解く上で革命的だったのが、1959年にソ連の数学者Sergei Godunovが提案した「Godunov法」です。それまでの差分法は衝撃波まわりで数値振動(ギブス現象)が抑えられず苦労していたのですが、Godunovは「各セル界面でRiemann問題を解く」というアイデアを持ち込んだ。実は元の論文にはGodunov自身による「この1次精度スキームはあまり精度が高くない」という正直なコメントが書いてあったといわれています。それでも現代のHLLC・Roeスキームに至る衝撃波計算の系譜はすべてこのアイデアの延長線上にある。
オイラー方程式(圧縮性)の実務適用
実践ガイド
先生、オイラー方程式ソルバーって実務のどこで使うんですか? 粘性を無視してるのに意味あるんですか?
いい疑問だね。オイラーソルバーは次のような場面で今でも重要だ。
- 初期設計段階での空力評価(翼型の圧力分布、衝撃波位置)
- ミサイルや発射体の抗力推算
- Navier-Stokesソルバーの検証用ベンチマーク
- 非粘性効果が支配的な問題(強い衝撃波、爆風波)
ベンチマーク問題
オイラーソルバーの検証ってどうやるんですか?
標準的なベンチマーク問題がいくつかあるよ。
| 問題 | マッハ数 | 検証ポイント |
|---|---|---|
| Sod衝撃波管 | 内部流 | 衝撃波・接触不連続・膨張波の捕獲 |
| 双楔反射 | M=2〜4 | 衝撃波の反射と交差 |
| NACA0012超音速 | M=0.8〜1.2 | 遷音速での衝撃波位置 |
| Shu-Osher問題 | M=3 | 衝撃波-渦干渉の解像度 |
| Sedov爆風波 | 高圧比 | 球対称爆風の自己相似解との比較 |
Sod問題は必ず最初にやるべきだ。解析解が存在するから、スキームの拡散特性やリミッターの影響を定量的に評価できる。Sod問題の初期条件は
で、$t = 0.2$ での衝撃波位置、接触不連続のシャープさ、膨張波の滑らかさを比較する。
実際にメッシュはどのくらい細かくすればいいんですか?
Sod問題なら100セルでスキームの基本特性が見える。400セルで2次精度の効果が明確になる。実問題では衝撃波近傍にAMR(Adaptive Mesh Refinement)を使うのが効率的だ。衝撃波の数値的厚みは概ね3〜5セルになるから、衝撃波の構造に興味がなければそれ以上の細分化は不要だよ。
境界条件の設定
境界条件の設定って、マッハ数で変えなきゃいけないんですよね?
そう、特性理論に基づいて正しく設定する必要がある。
| 境界 | 亜音速 (M<1) | 超音速 (M>1) |
|---|---|---|
| 入口 | 全圧・全温指定、静圧は外挿 | 全変数を指定 |
| 出口 | 静圧指定、他は外挿 | 全変数を外挿 |
| 壁面 | 法線速度=0(slip条件) | 法線速度=0 |
| 遠方場 | 特性ベースNRBC | 全変数指定(流入側) |
遠方場境界の設定ミスで反射波が出るのは、特性の取り扱いが原因だったんですね。
そのとおり。Fluentのfar-field境界条件は特性ベースで実装されているから、Mach numberとflow directionを正しく指定すれば問題ない。OpenFOAMでは freestreamPressure / freestreamVelocity boundary conditionが対応するよ。
オイラー方程式で翼を計算するとき、境界層はどこへ行った?
オイラー方程式は粘性項を含まないため、翼まわりを計算しても「境界層」という概念が生まれません。それなのに実務ではオイラー計算で揚力係数の見積もりをすることがあります。なぜ使えるかというと、揚力は主に圧力分布で決まり、粘性の寄与は比較的小さいから。ただし抗力(特に摩擦抗力)はオイラーでは正確に出ません。航空機の初期コンセプト設計で「とにかく揚力とおおまかな衝撃波位置を速く掴みたい」というケースにオイラー計算が今でも現役で使われる理由がここにあります。
オイラー方程式(圧縮性)のソフトウェア比較
商用ツールでのオイラーソルバー
商用CFDソフトでオイラー方程式を解く場合、各ソフトの違いってありますか?
主要ソルバーでの設定方法と特徴を比較してみよう。
Ansys Fluent
Fluentでは粘性モデルを「Inviscid」に設定するとオイラー方程式が解かれる。
- Solver: Density-Based を推奨(圧力ベースでも可能だがM>1では密度ベースが安定)
- Flux Type: Roe-FDS がデフォルト。AUSM は低マッハ数に強い
- Spatial Discretization: 2nd Order Upwind 以上
- 2024R2以降、Mosaic meshingとAMRの組み合わせが強力
Simcenter STAR-CCM+
STAR-CCM+では物理モデルの選択で「Inviscid」を指定。結合型ソルバーが標準で、陰的時間積分のため大きなCFL数が使える。
- Convection scheme: 2nd Order, MUSCL 3rd Order Reconstruction が選択可能
- Flux Function: Roe, AUSM+, HLLC
- Mesh morphing と overset mesh が超音速飛翔体解析で有用
OpenFOAM
OpenFOAMでオイラー方程式を解くにはどのソルバーを使うんですか?
rhoCentralFoam が最適だ。KNP(Kurganov-Noelle-Petrova)スキームベースで、人工粘性不要のcentral-upwindスキームを採用している。粘性項を計算しないようにtransportPropertiesでmuを0に設定すればオイラーソルバーとして動作する。
チュートリアルの shockTube がSod問題に対応しているから、まずはこれを動かしてみるといい。
オイラーソルバーとNavier-Stokesソルバーの使い分けの判断基準はありますか?
粘性効果が結果に与える影響の大きさで判断する。
- オイラーで十分: 衝撃波の位置・強度、圧力波の伝播、爆風解析
- N-Sが必要: 抗力(摩擦抗力成分)、熱伝達、剥離を伴う流れ、境界層効果
設計初期段階でオイラー解析を回して衝撃波パターンを把握し、その後N-S解析で詳細評価するのがワークフローの定石だ。Eulerで出した圧力分布は、衝撃波が支配的な流れ場ではN-S解とかなり一致するからね。
なぜ宇宙機開発ではオイラー専用ソルバーが今も生き残るのか
粘性ありのNavier-Stokes計算が主流の今でも、宇宙機・ミサイル設計の初期フェーズではオイラー専用コードが使われ続けています。理由はシンプルで「速い」から。NASAのCART3DやEuler3Dのようなオイラーソルバーは、Navier-Stokesコードの10〜50倍速く計算できます。設計初期に100以上の形状バリエーションを流したいとき、1形状数分で終わるオイラー計算は非常に重宝します。完全な精度は求めず「衝撃波の位置と大まかな圧力分布を掴む」という用途に絞れば、今でも現役で十分な道具なのです。
オイラー方程式(圧縮性)の先端研究
先端トピック
オイラー方程式の数値解法って、もう完成された分野じゃないんですか?
そう思われがちだけど、実はまだ活発に研究されている分野だよ。特に高次精度スキーム、全マッハ数対応、AMR(適応格子細分化)の3つが注目されている。
高次精度スキーム
WENO(Weighted ENO)スキームは5次精度で衝撃波近傍でも振動のない解を得られる。最近はWENO-Z、TENO(Targeted ENO)といった改良版が提案されていて、より少ない数値散逸で衝撃波と渦を同時に高精度に捕捉できる。
$\omega_k$ が非線形重みで、滑らかな領域では最適な高次精度を回復し、不連続近傍では自動的にENO的な選択を行う。
DG法は有限要素法の枠組みで高次精度を実現する手法で、要素間の不連続を許容するからGodunovスキームとの親和性が高い。$p$-refinement(多項式次数の向上)が容易で、h-p適応が可能。ただし自由度数が多くメモリと計算コストが課題だ。最近はFR(Flux Reconstruction)法というDGの効率的な再定式化が広まりつつある。
全マッハ数対応スキーム
従来のGodunovスキームは低マッハ数($M \ll 1$)で精度が低下する問題があった。圧力の数値散逸が $O(1/M)$ でスケールするため、非圧縮極限で正しい解が得られない。
これを解決するのが全マッハ数スキーム(all-speed scheme)で、AUSM+-up、SLAU(Simple Low-dissipation AUSM)、Low-Mach preconditioned Roeなどがある。
全マッハ数対応って、遷音速の翼まわりの流れみたいに局所的にM=0に近い領域とM>1の領域が混在する場合に重要ですか?
まさにそのとおり。旅客機の巡航条件($M_\infty = 0.85$)では、翼上面で局所マッハ数が1.2程度まで上がる一方、よどみ点近傍ではM→0になる。全マッハ数スキームならこの両方を正しく解ける。
AMRとimmersogeometric解析
AMR(Adaptive Mesh Refinement)はoctreeベースの格子で衝撃波近傍を動的に細分化する技術だ。AMReXやp4estといったライブラリが公開されていて、数十億セルの計算も可能になっている。
また、immersed boundary法やcut-cell法を組み合わせることで、複雑形状に対するメッシュ生成の手間を大幅に削減できる。Cartesian grid + AMR + cut-cellの組み合わせは「meshless」に近いワークフローを実現するんだ。
メッシュ生成が楽になるのは嬉しいです。今後の方向性として、GPUの活用はどうですか?
GPU対応は急速に進んでいる。陽解法のオイラーソルバーは並列性が高いからGPUとの相性が非常に良い。JAXAのFaSTARやNASAのFUN3DもGPU対応版が開発されているよ。
オイラー方程式で「エントロピー生成」が起きる矛盾
オイラー方程式は粘性なし・熱伝導なしの「理想流体」を記述する。本来なら可逆過程だから、エントロピーは変化しないはずだ。ところが衝撃波を含む解を計算すると、衝撃波の通過後にエントロピーが不連続に増加する——これは物理的には正しいが、数式的には矛盾に見える。この謎の答えは「衝撃波がオイラー方程式の弱解(weak solution)であり、古典解が存在しない」という数学的事実にある。CFDの先端研究では今もエントロピー条件と弱解の扱いが議論のテーマになっており、数値スキームの設計に直結する。
オイラー方程式(圧縮性)のトラブル対応
トラブルシューティング
先生、オイラーソルバーで計算がうまくいかないときの対処法を教えてください。
よくあるトラブルをパターン別に整理しよう。
1. 膨張衝撃波(Expansion Shock)の発生
症状: 物理的にありえない膨張衝撃波が解に含まれる
原因: Roeスキームはエントロピー条件を自動的には満たさない。固有値がゼロに近い場合にentropy violationが起こる。
対策: エントロピー修正(Harten-Hyman fix)を適用する。
$\delta$ は小さな正の値(例: $0.1 \times (|u| + a)$)。Fluentではデフォルトでentropy fixが適用されているから通常は問題ない。
2. Carbuncle現象
Carbuncle現象って何ですか?
症状: bow shock(離脱衝撃波)の前面にニキビ状の凸凹が発生する。特に格子が衝撃波に対して正対している場合に起きやすい。
原因: Roe等のリーマンソルバーの固有の不安定性。接触不連続に対する散逸が不足。
対策:
- HLLCソルバーに切り替え(接触波の散逸が適度にある)
- Roeソルバーに人工的なshear viscosityを追加
- メッシュを衝撃波に対してわずかに斜めにする(実用的な回避策)
3. 負の密度・圧力
症状: 計算中に密度や圧力が負値になって発散
原因: 高次精度再構築で、衝撃波近傍でオーバーシュート/アンダーシュートが発生。リミッターが不十分。
対策:
- より散逸的なリミッター(minmod)に切り替え
- positivity-preserving limiter を使用
- CFL数を下げる(0.3〜0.5程度に)
- failsafe機構: 負値が検出されたら局所的に1次精度に落とす
圧力が負になるって怖いですね。どうやって事前に防げますか?
強い衝撃波(圧力比10以上)の問題では、最初から保守的な設定(1次精度、小さいCFL)で計算を開始し、安定したら徐々に高次精度に切り替えるのが安全だ。
4. 境界からの非物理的反射
症状: 計算領域の境界から衝撃波や圧力波が反射して、定常解が得られない
対策:
- 遠方場境界にcharacteristic-based NRBC(Non-Reflecting Boundary Condition)を使用
- 計算領域を十分大きく取る(翼弦長の20〜50倍)
- スポンジ層(ダンピングゾーン)を境界近傍に追加
オイラーソルバーは粘性がない分、数値的な問題が目立ちやすいんですね。
その通り。Navier-Stokes方程式では粘性項が自然な散逸を提供するから安定化されるが、オイラー方程式では数値スキームの散逸だけが安定化要因になる。だからスキームの選択とパラメータ調整がより重要なんだ。
「計算が爆発する」——オイラー解析あるある発散の原因トップ3
オイラー方程式の計算が突然発散するとき、原因の大半は「初期条件の設定ミス」「CFL数の過大設定」「衝撃波とメッシュの位置関係」の3つに集約されます。特に多いのが初期条件問題で、「フリーストリーム一様流で初期化したらノズル喉部付近で即発散」というケース。喉部近傍でマッハ数が1に近づくと特性線が鋭くなり、わずかな不連続が増幅します。対策として現場では「マッハ数を0.1くらいから徐々に上げるランプアップ手法」や「局所的に圧力・密度を滑らかにするピーシング初期化」が使われます。CFL数は0.3〜0.5から始めるのが無難です。
なった
詳しく
報告