MMS: Spalart-Allmaras乱流モデル

カテゴリ: V&V・製造解の方法 | 2026-10-01 改稿
CAE visualization for mms turbulence sa theory - technical simulation diagram
MMS: Spalart-Allmaras乱流モデル

理論:答えを先に決めて源項を作る

概要

🙋

乱流モデルの実装が正しいかは、実験と比べればわかるのではないですか?

🎓

実験と合わないとき、それがモデルの限界なのか、プログラムの誤りなのか区別できない。たまたま合っても、誤りが別の誤差と打ち消し合っているだけかもしれない。プログラムが方程式を正しく解いているかは、実験とは別に確かめる必要がある。これが検証(verification)で、その強力な方法が製造解の方法だ。解を自分で決めて方程式に代入し、余った分を源項として足せば、決めた解が厳密解になる。あとは格子を細かくしたときに誤差が理論どおりの速さで減るかを見るだけでよい。SAモデルのように非線形な項がいくつもある方程式では、項の入れ忘れや係数の誤りがよく起き、この方法で見つかる。

SAモデルを簡略化した1次元の方程式

$$ \frac{d}{dy}\!\left[\frac{\nu + \tilde\nu}{\sigma}\frac{d\tilde\nu}{dy}\right] + \frac{c_{b2}}{\sigma}\left(\frac{d\tilde\nu}{dy}\right)^2 + c_{b1} S\,\tilde\nu - c_{w1}\left(\frac{\tilde\nu}{d}\right)^2 + f(y) = 0 $$

SAモデルの ν̃ の輸送方程式から、対流を除き、消滅項の補正関数を fw = 1、渦度の大きさを S = 1 とした。定数は σ = 2/3、cb1 = 0.1355、cb2 = 0.622、cw1 = cb1/κ² + (1 + cb2)/σ = 3.239(κ = 0.41)。区間 0.1 ≤ y ≤ 1、壁からの距離 d = y、分子動粘度 ν = 0.1(無次元)。

製造解と観測次数

$$ \tilde\nu_m(y) = 1 + 0.5\sin(2\pi y),\qquad f = -\mathcal{L}(\tilde\nu_m),\qquad p = \frac{\ln(e_h / e_{h/2})}{\ln 2} $$

L は左辺の演算子。源項 f は記号計算(sympy)で求めた。e_h は格子幅 h でのL2誤差。

Coffee Break よもやま話

SAモデルと検証の文化

Spalart-Allmarasモデルは1992年に航空機の外部流れ向けに発表された1方程式の乱流モデルで、扱いやすさから広く使われている。一方で、実装によって細部の扱いが少しずつ違い、同じ問題でも結果が合わないことが問題になった。NASAのTurbulence Modeling Resourceは、標準的な式の形と検証用の問題をまとめ、複数のプログラムで格子を細かくして同じ答えに収束するかを確かめられるようにしている。製造解による検証も、こうした実装の確認に使われている。

計算例

例:格子を細かくしたときの誤差(2次精度中心差分、ニュートン法で解く)

実装20分割40分割80分割160分割観測次数
正しい実装3.40×10⁻³8.38×10⁻⁴2.08×10⁻⁴5.19×10⁻⁵2.02 / 2.01 / 2.00
cb2項の入れ忘れ1.37×10⁻¹1.34×10⁻¹1.33×10⁻¹1.33×10⁻¹0.03 / 0.01 / 0.01
面の拡散係数を片側の値で計算7.98×10⁻³3.32×10⁻³1.55×10⁻³7.52×10⁻⁴1.26 / 1.10 / 1.04

観測次数は、隣り合う2つの格子の誤差の比から求めた。正しい実装は理論どおり2に収束する。

🙋

cb2項を忘れると誤差がまったく減らないのに、片側の拡散係数だと誤差は減っていくんですね。

🎓

2つは性質の違う誤りだ。cb2項の入れ忘れは、プログラムが解いている方程式そのものが違うので、格子をいくら細かくしても別の方程式の解に近づくだけで、誤差は0.13あたりで止まる。これは見つけやすい。一方、面の拡散係数を片側の値で計算する誤りは、方程式は正しく、離散化の精度だけが2次から1次に落ちている。格子を細かくすれば正しい解に近づくので、1つの格子で結果を見ても、実験との比較でも気づきにくい。でも観測次数を見れば1.04で、理論の2と合わないので一目でわかる。次数の検証は、誤差の大きさだけでなく、減り方を見ることに意味がある。

検証の進め方

  1. 滑らかで、すべての項がはたらく製造解を選ぶ(どの微分も0にならないもの)。
  2. 源項は記号計算で作り、手計算での転記ミスを避ける。
  3. 3〜4段階の格子で誤差を計算し、観測次数を求める。
  4. 次数が理論と合わなければ、項を1つずつ外して原因を絞る。
  5. 壁関数や制限関数など、精度を下げる仕組みは検証中は切っておく。
Coffee Break よもやま話

「平板の摩擦は合っていたのに」

自作のソルバーにSAモデルを入れ、平板の境界層で壁面摩擦を実験と比べたところ、数%以内で一致した。ところが製造解で検証すると、観測次数が2にならず、誤差が途中で止まった。調べると、cb2項の係数で σ で割るのを忘れていた。平板の境界層ではこの項の影響が比較的小さく、ほかの誤差に紛れて見えなかった。はく離のある流れでは結果が大きく変わり、修正後は別のプログラムの結果とよく合うようになった。

よくある間違い

間違いと対策

間違い起きること対策
実験との一致で検証済みとする実装の誤りが残るMMSで次数を確認
一部の項が0になる製造解その項の誤りを見逃す全項がはたらく解
源項を手で導出源項自体の誤り記号計算
格子が粗すぎる漸近域に入らない4段階以上で確認
制限関数を入れたまま次数が下がる検証中は切る
🙋

関連する内容も知りたいです。

関連シミュレーター

この分野のインタラクティブシミュレーターで理論を体感しよう

シミュレーター一覧

関連する分野

構造解析流体解析熱解析
この記事の評価
ご回答ありがとうございます!
参考に
なった
もっと
詳しく
誤りを
報告
参考になった
0
もっと詳しく
0
誤りを報告
0