対流項の発散スキームの比較

2019年6月26日

はじめに

対流項の発散スキームを比較する。

使用バージョン

OpenFOAM v6

ファイル

対流項の発散スキームの比較

scalarTransportFoam で 1 次元問題を解き、対流項 div(phi,T) のスキームを比較する。以下のスキームを比較する。

内挿スキーム

  • linear
  • cubic

風上差分スキーム

  • upwind
  • linearUpwind
  • linearUpwind with faceLimited
  • QUICK
  • quadraticUpwindFit
  • cubicUpwindFit

TVD スキーム

  • Minmod
  • SuperBee
  • vanLeer
  • vanAlbada
  • UMIST
  • MUSCL
  • limitedLinear
  • limitedCubic

NVD スキーム

  • Gamma

制限関数適用 (limitWith)

  • QUICK with Minmod
  • QUICK with SuperBee
  • QUICK with vanAlbada

勾配制限

  • faceLimited
  • cellLimited Minmod (Barth and Jespersen)
  • cellLimited Venkatakrishnan
  • cellLimited cubic (Michalak and Ollivier-Gooch)

1 次元領域のスカラー値を初期値 0 とし、左端のスカラー値を 1 に設定、左から右に向かって流れを作る。スカラーの前面が領域の中央に来た時のスカラーの分布形状を比較する。速度は 1m/s、格子幅は 0.01 m、時間刻み幅は 1e-4 s、時間のスキームは Euler である。

中心差分・風上差分スキームの結果を示す。

upwind の結果ががかなりなまっている (数値拡散が大きい) のがわかる。また linear や cubic といった単純な内挿スキームは大きな振動が起こっていることが見て取れる。

差がよく見えるように表示範囲を絞る。

linearUpwind で少しオーバーシュートやアンダーシュートが起こっている。"linearUpwind with faceLimited" は勾配計算に faceLimited を使用したものだが、制限なしの場合に比べてオーバーシュートやアンダーシュートが抑えられている。QUICK と quadraticUpwindFit (セル界面の値を 2 次多項式で風上側に補間する) の結果はほぼ一致しているが、これは QUICK の定式化の仕方を考えると当然であろう。cubicUpwindFit が拡散が小さくて精度がよさそうであるが、QUICK もそれに近い結果を与えている。

つぎに、代表的な TVD スキーム (2 次精度 TVD スキーム) の結果を示す。

ちょうど Minmod と SuperBee の間に収まっているのがわかる。制限付き linearUpwind は MUSCL とまったく同じ結果を示している。

Gamma スキームの結果を示す。

これは TVD スキームと少し異なった傾向を示している。

続いて、limitedLinear スキームの結果を示す。

これは次の図に示しているように、Gamma スキームの結果に近い。パラメタを 0 にしたものは両者の結果が一致している (重なってしまって見えないが)。安定側であるパラメタ 1 は Gamma のほうが拡散が大きい。

limitedCubic スキームの結果を示す。

ちょうどオーバーシュートを抑えた QUICK のような挙動を示している。

limitedLinear と limitedCubic の比較を次に示す。limitedCubic のほうが若干拡散が小さいのがわかる。

limitWith を使って QUICK に制限関数を適用した結果を示す。

拡散が大きくなったりアンダーシュートが起こったり、あまりいい結果は得られていない。UMIST は QUICK に vanAlbada 制限関数を適用したものに近い。QUICK の制限版である UMIST よりもむしろ limitedCubic のほうが QUICK の結果に近い。

linearUpwind に提供する勾配制限のタイプを比較した結果を次に示す。拡散の様子はあまり変わらない。cubic については、パラメータを大きくすると拡散が大きくなっている。

おまけ: WENO

Tobias Martin 氏により WENO スキームのライブラリ WENOEXT が公開されている。libWENOEXT を wmake でコンパイルする (最新の OpenFOAM で使うためには多少修正が必要: 修正版)。使う時には controlDict でライブラリの指定がいる。

libs ( "libWENOEXT.so" );

実行時にメッシュから constant に WENOBase2 あるいは WENOBase3 ディレクトリが作成される。あれば読み込むようになるが、情報がメッシュと関連しているため、メッシュが変わったら消して実行する。

結果を以下に示す。少し振動している。

どれを選ぶべきか

スキームを選ぶ観点は、以下のものがあげられる。

  • 安定性
  • 精度 (数値拡散の小ささ)
  • オーバーシュート・アンダーシュートを許容できるか
  • 計算速度

安定性と精度は相反する。また、オーバーシュート・アンダーシュートを抑えるとどうしても数値拡散が導入されてしまう。安定性を最優先するなら upwind を選べばよい。振動を抑制しつつ精度の高いものが欲しいなら TVD/NVD スキームを使えばよい。TVD/NVD スキームは、安定側から高精度側までいくつも提案されており、性能もさほど変わらないので、どれを選ぶかは好みの問題のように見える。ただし、計算時間の問題がある。以下に、計算時間を計測した結果を示す。

各スキームを 10 回ずつ計算して平均を取り、upwind の時間で割っている。つまり、upwind の何倍の時間がかかるかを示している。平均を取っても計測のたびに結果が結構入れ替わるので、ひとつの目安である。だいたい計算手順の複雑さに比例して時間がかかっている。2次精度 TVD は SuperBee が高精度な上に計算速度も速いことになる。UMIST は TVD の中では時間がかかっている。ざっくりいえば、除算があると時間がかかるようである。linearUpwind は TVD より時間がかかっている。勾配制限が入るとさらにかかる。勾配制限付き linearUpwind は "MUSCL" と結果的には同じもののはずだが、計算時間は結構違う。この結果からすると、わざわざ linearUpwind を使う理由はないかもしれない。withLimited 系は時余計な処理を負荷するわけだから、当然時間がかかる。

ほんとか？ simpleFoam/pitzDaily ケースでの性能を調べてみた。チュートリアルケースのまま最大ステップ数を 500 に設定し、収束判定で止まるものは止まるままとした。速度の対流項のみスキームを変更する。まず、ステップ数。

"QUICK with SuperBee" は発散したので除いた。主要な TVD スキームのほとんどは収束せずに最大ステップ数に達している。結果を見ると、速度分布がゆらゆらとゆらめいており、落ち着く様子はない。linearUpwind は upwind についで少ないステップ数で収束している。勾配制限の差はほとんどなさそうである。基本的には勾配制限付き linearUpwind と同一のスキームであると考えられる "MUSCL" が異なる挙動を示している。メッシュの不均一性により差が出たか？ つづいては計算時間。

ステップ数が少ないだけあってか、linearUpwind 系が上位を占めている。1 ステップあたりの計算時間が多くても、トータルの時間でスピードを稼いでいるわけである。勾配制限の差はあまり大きくない。強いて言えば、cellLimited では Venkatakrishnan が一番速いが、そう言い切ってよいのだろうか？ 後続に limitedLinear および Gamma が来ている。挙動が似ているためか、性能も似ている。quadraticUpwindFit、cubicUpwindFit は、途中で収束しているにも関わらず時間がかかっている。

1 次元のスカラー輸送での検討結果だけからは予想できない結果になった。問題によって変わりそうでもある。以上の検討結果からなんとなく判断すると、スキーム選択の指針は次のようになろうか。

  • 安定性重視なら upwind。
  • 数値粘性を抑えたいなら、速度については linearUpwind、limitedLinear、Gamma あたりにしておくのが無難か。
  • 速度以外も基本は同じ？ 濃度や体積分率のようなものは TVD 系がよいのかも。
  • ちなみに、LES のケースだと LUST (linear と linearUpwind の混合スキーム) と limitedLinear などが使われていたりする。