Conversion from point-value flux to numerical flux

通常、有限差分ではフラックスに対して再構築を行うため、HLLD法といった解像度の高い近似Riemann解法を用いることができない。 密度や圧力の正値性などの観点からも、フラックスよりも基本量に対していろいろ操作したほうが都合がいいことが多いので、このフラックス再構築は有限体積と比してデメリットであると言える。 しかも純粋な有限差分ではCT法とカップルさせるのも難しい。 でも有限差分の利点である dimension-by-dimension による多次元化の恩恵は受けたい。

そういう諸々の要請をまとめて解決するのが、ここで採用している「基本量補間」+「フラックス変換」の手法。

フラックスと数値フラックス

3次以上の高次精度有限差分スキーム(仮に5次精度としよう)を用いるときには、セル中心で求めたフラックスをセル境界へ再構築(reconstruction)するのが常套手段であった。 この時、セル境界で求まるのは数値フラックス \(\hat{\mathbf{F}}\) であって、次のような特性を持つ。 \begin{align*} \left.\frac{\partial \mathbf{F}}{\partial x}\right|_i = \frac{\hat{\mathbf{F}}_{i+1/2} - \hat{\mathbf{F}}_{i-1/2}}{\Delta x} + \mathcal{O}\left(\Delta x^5\right) \end{align*} 要するに、フラックス値自体は5次精度ではないが、二点差分を計算するとセル中心でちょうど5次精度の微分値が得られる、という代物。 大事なのは、ここで使った再構築は、有限体積で普通に使う再構築と全く同じであるということ。 再構築ではなく補間(interpolation)をしてしまうと、5次精度のフラックスは求まるが、5次精度の数値フラックスにはならない。

ところがこのフラックスと数値フラックスの間には、変換公式が存在することが知られている。 5次精度のフラックス変換に必要となるのは6次精度の変換公式で、次のようなもの。 \begin{align*} \hat{\mathbf{F}} = \mathbf{F} - \frac{1}{24}\frac{\partial^2 \mathbf{F}}{\partial x^2}\Delta x^2 + \frac{7}{5760}\frac{\partial^4 \mathbf{F}}{\partial x^4}\Delta x^4 + \mathcal{O}\left(\Delta x^6\right) \end{align*} 2階微分と4階微分を評価する必要があるが、2階微分には\(\Delta x^2\)が、4階微分には\(\Delta x^4\)が掛かっているので、それぞれ4次精度と2次精度で評価すれば適切な精度が保たれることになる。

再構築から補間へ

フラックスから数値フラックスへの変換公式がわかったので、もはや再構築に縛られる必然性はなくなった。 我々は、好きなところから好きなところへ、好きな物理量を、好きな精度で「補間」すればよい。 そして微分値が必要になったらその都度、補間された物理量を数値フラックスに変換して二点差分を評価すればよい。

こうなってしまうとフラックス補間を用いる道理はない。 基本量をセル中心からセル境界へ単に補間して、好きな近似Riemann解法を使ってフラックスを計算する。 フラックスを数値フラックスに変換してから、保存量を更新する。 保存量や基本量はあくまで point value なので、dimension-by-dimension の多次元化が許される。 しかも有限体積では許されなかった、掛け算や割り算も自由に行える (一般に積の平均と平均の積は異なるので、有限体積では平均値の積は適切な処置をしないと精度が落ちる)。

さらに、セル境界で求めたフラックスはあくまで point value なので、こいつも好きに補間してよい。 セル境界からセルコーナーに電場を補間してから数値フラックスに変換すれば、容易にCT法とカップルできる。 point value としていろいろな物理量を保持しているため、あらゆる操作に対して非常に自由度が高く、 さまざまな拡張が容易になっている。 なんと夢のようなスキームか。

まとめると

  1. セル中心からセル境界へ基本量を補間
  2. セル境界で近似Riemann解法を噛ませてフラックスを計算
  3. セル境界でフラックスを数値フラックスに変換
  4. 二点差分でセル中心の微分値を評価、保存量更新
となる。必要ならCT用に同じ要領で電場を補間し、face-centerの磁場を更新。

欠点

一応欠点らしい欠点を述べておくと、まず必要なステンシルが増える。 これは並列化した際の通信量の増大に直結するので、苦しいようなら隠蔽とかで頑張る。 計算量の増加はそれほど問題にならないと思う。 積が自由に計算できたりするので、有限体積と比べると計算量はむしろ少ない。 有限差分ベースなので、多次元だとことさらに少ない。

あとはまぁ、なんだろね。 めんどくさいとか?