数学とか語学とか楽しいよね

フランス語、ドイツ語、ロシア語、アラビア語、オランダ語、英語、スペイン語、ラテン語とか数学とか数値計算(有限要素法、有限体積法、差分法、格子ボルツマン法、数理最適化、C++コード付き)とか勉強したことをまとめます。右のカテゴリーから興味のある記事を探してください。最近はクラシックの名演も紹介しています。noteにも書いています。https://note.mu/kumag_min

【Navier-Stokes方程式】フラクショナルステップ法によるNavier-Stokes方程式の離散化

ラクショナルステップ法(部分段階法)によるNavier-Stokes方程式の離散化について説明していきます。この手法は流速と圧力を直接未知数として計算する手法です。SMAC法では中間流速(予測子)を計算する際に圧力を使いましたが、フラクショナルステップ法では圧力を使わず「仮の流速」を計算します。最終的にはMAC法のように圧力に関するPoisson方程式が導かれます。

以下の内容は『流体解析の基礎』のpp.75-76を参考にしています。

流体解析の基礎

流体解析の基礎


簡単のため二次元直交座標系で考えます。まず保存型のNavier-Stokes方程式

 \displaystyle
\begin{eqnarray} 
\frac{\partial u}{\partial t} + \frac{\partial (u^2)}{\partial x} + \frac{\partial (uv)}{\partial y} &=& -\frac{\partial p}{\partial x} + \frac{1}{Re} \left( \frac{\partial^2 u}{\partial x^2} + \frac{\partial^2 u}{\partial y^2}\right) \\
\frac{\partial v}{\partial t} + \frac{\partial (uv)}{\partial x} + \frac{\partial (v)^2}{\partial y} &=& -\frac{\partial p}{\partial y} + \frac{1}{Re} \left( \frac{\partial^2 v}{\partial x^2} + \frac{\partial^2 v}{\partial y^2}\right)
\end{eqnarray}

です。連続式は

 \displaystyle \frac{\partial u}{\partial x} + \frac{\partial v}{\partial y} = 0

です。Navier-Stokes方程式と連続式を以下のように時間に対して離散化します。

 \displaystyle
\begin{eqnarray} 
\frac{u^{n+1}-u^n}{\Delta t} + \frac{\partial (u^n)^2}{\partial x} + \frac{\partial (u^nv^n)}{\partial y} &=& -\frac{\partial p^{n+1}}{\partial x} + \frac{1}{Re} \left( \frac{\partial^2 u^n}{\partial x^2} + \frac{\partial^2 u^n}{\partial y^2}\right) \\
\frac{v^{n+1}-v^n}{\Delta t} + \frac{\partial (u^nv^n)}{\partial x} + \frac{\partial (v^n)^2}{\partial y} &=& -\frac{\partial p^{n+1}}{\partial y} + \frac{1}{Re} \left( \frac{\partial^2 v^n}{\partial x^2} + \frac{\partial^2 v^n}{\partial y^2}\right) \\
\frac{\partial u^{n+1}}{\partial x} + \frac{\partial v^{n+1}}{\partial y} &=& 0
\end{eqnarray}

ここまではSMACと同じです。次に、以下のように圧力を省いてnステップの値を用いて仮の流速  u^* v^* を計算します。

 \displaystyle
\begin{eqnarray} 
\frac{u^{*}-u^n}{\Delta t} + \frac{\partial (u^n)^2}{\partial x} + \frac{\partial (u^nv^n)}{\partial y} &=& \frac{1}{Re} \left( \frac{\partial^2 u^n}{\partial x^2} + \frac{\partial^2 u^n}{\partial y^2}\right) \\
\frac{v^{*}-v^n}{\Delta t} + \frac{\partial (u^nv^n)}{\partial x} + \frac{\partial (v^n)^2}{\partial y} &=& \frac{1}{Re} \left( \frac{\partial^2 v^n}{\partial x^2} + \frac{\partial^2 v^n}{\partial y^2}\right)
\end{eqnarray}

「仮の流速」と呼ばれているのは、この流速はNavier-Stokes方程式を満たさないからです。対応する運動方程式をそれぞれ引き算すると

 \displaystyle
\begin{eqnarray} 
\frac{u^{n+1}-u^*}{\Delta t} &=& -\frac{\partial p^{n+1}}{\partial x} \\
\frac{v^{n+1}-v^*}{\Delta t} &=& -\frac{\partial p^{n+1}}{\partial y}
\end{eqnarray}

を得ます。これが次のタイムステップの流速を決める式です。ここから  u^{n+1} v^{n+1} を連続式に代入します(n+1ステップにおいて連続式を満たすとする)。

 \displaystyle
\begin{eqnarray} 
\frac{\partial^2 p^{n+1}}{\partial x^2} + \frac{\partial^2 p^{n+1}}{\partial y^2} = \frac{1}{\Delta t} \left( \frac{\partial u^{*}}{\partial x} + \frac{\partial v^{*}}{\partial y} \right)
\end{eqnarray}

これがフラクショナルステップ法における圧力のPoisson方程式です。

計算手順は以下のように1~3の繰り返しとなります。

1. nステップでの流速  \boldsymbol{(u^n, v^n)} を用いて仮の流速  \boldsymbol{(u^*, v^*)} を得る
2. 仮の流速  \boldsymbol{(u^*, v^*)} を用いて圧力のPoisson方程式を解き、次のタイムステップの圧力  \boldsymbol{p^{n+1}} を求める
3. 次のタイムステップの圧力  \boldsymbol{p^{n+1}} を用いて、流速を補正(流速の式を使う)し、次のタイムステップの流速  \boldsymbol{(u^{n+1}, v^{n+1})} を得る

空間方向の離散化には、保存型をそのまま用いれば有限体積法を、非保存型を用いれば有限要素法や差分法を用いることができます。