📂流体力学

ナビエ・ストークス方程式の導出

定理

$$ \mathbf{u} = \mathbf{u} \left( t ; \mathbf{x} \right) = \left( u_{1} \left( t ; \mathbf{x} \right) , u_{2} \left( t ; \mathbf{x} \right) , u_{3} \left( t ; \mathbf{x} \right) \right) $$ 特に、3次元空間で視点 $t$ と空間座標 $\mathbf{x} = \left( x_{1} , x_{2} , x_{3} \right)$ の速度場を上のような速度ベクトルで表すとする。それに類似して、$p : \mathbb{R}^{3} \to \mathbb{R}$ は各座標でかかる圧力 $p = p \left( \mathbf{x} \right)$ を表す。$\mathbf{u}$ が非圧縮性のニュートン流体の速度であれば、次の支配方程式に従う。 $$ {\frac{ \partial \mathbf{u} }{ \partial t }} + \left( \mathbf{u} \cdot \nabla \right) \mathbf{u} = - \nabla w + \nu \nabla^{2} \mathbf{u} + \mathbf{g} $$ ここで $\nabla \cdot$ は発散、$\nu = \mu / \rho$ は動粘性係数 $\mu$ と密度 $\rho$ の比、$\nabla w = \nabla p / \rho$ は熱力学的仕事^{thermodynamic work}、$\mathbf{g}$ は重力加速度である。

説明

ナビエ–ストークス方程式^{Navier-Stokes equation}は流体力学におけるオイラー方程式にニュートン流体を仮定して粘性項 $\nabla^{2} \mathbf{u}$ が追加された支配方程式である。

3次元空間でこの偏微分方程式の滑らかな解が常に存在するか、もし存在しないならば反例を見つけられるかという問題はナビエ–ストークスの存在性と滑らかさ^{Navier–Stokes existence and smoothness}と呼ばれ、ミレニアム懸賞問題の一つとして残されている。

時折「解くだけで100万ドルもらえる」といった言い方をされることがあるが、難易度を考えると「解くだけで」という表現は適切ではない。流体力学の方程式に関する人類の知見が深まるということは、ほとんど天の機を読むようなものだ。もちろん非常に難しく、その導出を理解するだけでもかなりの学習が必要である。

導出

流体力学におけるオイラー方程式: $$ {\frac{ \partial \mathbf{u} }{ \partial t }} + \left( \mathbf{u} \cdot \nabla \right) \mathbf{u} = - {\frac{ 1 }{ \rho }} \nabla p + \mathbf{g} $$

オイラー方程式から始める。物質微分 $D$ を用い、両辺に $\rho$ を掛けると次のようになる。 $$ \rho {\frac{ D \mathbf{u} }{ D t }} = - \nabla p + \rho \mathbf{g} $$ このままでは粘性による反発力を考慮していないので、右辺で作用する $- \nabla p$ を修正する。

コーシー応力テンソル: 主に物理学で、ある点において各方向およびせん断として作用する応力を成分に持ち、次のように定義される正方行列 $\mathbf{\sigma} \in \mathbb{R}^{3 \times 3}$ をコーシー応力テンソル^{Cauchy stress tensor}という。 $$ \sigma = \begin{bmatrix} \sigma_{11} & \sigma_{12} & \sigma_{13} \\ \sigma_{21} & \sigma_{22} & \sigma_{23} \\ \sigma_{31} & \sigma_{32} & \sigma_{33} \end{bmatrix} $$

オイラー方程式ではコーシー応力テンソルは $\sigma = - p I$ のように圧力項のみを含んでいた。体積が $V$ の流体塊に作用する力 $\mathbf{f}$ は、 $V$ を囲む表面 $\partial V$ に対する閉曲面積分 $\mathbf{f} = \oint_{\partial V} \sigma \cdot d S$ で表すことができる。

発散定理（ガウスの発散定理）: 3次元ベクトル関数 $\mathbf{F}$ に対して次が成り立つ。 $$ \int_{\mathcal{V}} \nabla \cdot \mathbf{F} dV = \oint_{\mathcal{S}} \mathbf{F} \cdot d \mathbf{S} $$ ここで $\nabla \cdot \mathbf{F}$ はダイバージェンス（発散）、$\int_{\mathcal{V}}$ は体積積分、$\oint_{\mathcal{S}}$ は閉曲面積分である。

発散定理により体積積分に変えると $\nabla \cdot \left( - p I \right) = - \nabla p$ となるので次を得る。 $$ \begin{align*} \mathbf{f} =& \oint_{\partial V} \sigma \cdot d S \\ =& \int_{V} \nabla \cdot \sigma d V \\ =& \int_{V} \nabla \cdot \left( - p I \right) d V \\ =& \int_{V} - \nabla p d V \\ \implies {\frac{ d \mathbf{f} }{ d V }} =& - \nabla p \end{align*} $$ ここで混同しやすいので少し $d \mathbf{f} / d V$ と $\rho \mathbf{g}$ の次元解析で単位を確認してみよう。 $$ \begin{align*} {\frac{ d \mathbf{f} }{ d V }} :& {\frac{ \text{force} }{ \text{volume} }} = \left[ {\frac{ kg \cdot m / s^{2} }{ m^{3} }} \right] = \left[ {\frac{ kg }{ m^{2} \cdot s^{2} }} \right] \\ \rho \mathbf{g} :& {\frac{ \text{mass} }{ \text{volume} }} \cdot \text{acceleration} = \left[ {\frac{ kg }{ m^{3} }} \right] \cdot \left[ {\frac{ m }{ s^{2} }} \right] = \left[ {\frac{ kg }{ m^{2} \cdot s^{2} }} \right] \end{align*} $$

要点は、オイラー方程式の導出はコーシー応力テンソル $\sigma$ に対する閉曲面積分によっても可能だということだ。より一般的には、重力による項である $\rho \mathbf{g}$ を除く影響は $d \mathbf{f} / d V = \nabla \cdot \sigma$ のようにコーシー応力テンソルのダイバージェンスを加えることで反映できる。ここで $\sigma$ に次のようにせん断応力 $\tau$ を追加しよう。 $$ \sigma = - p I + \tau $$

ニュートンの粘性法則: 流体力学で、非圧縮性かつ等方的な流体に作用する応力は速度の対称化された勾配に比例するというものをニュートンの粘性法則^{Newton’s law of viscosity}という。式としては応力テンソル $\tau \in \mathbb{R}^{3 \times 3}$ と速度場 $\mathbf{u}$ のヤコビアン $\nabla \mathbf{u}$ に対して次のように表される。 $$ \tau = \mu \left( \nabla \mathbf{u} + \left( \nabla \mathbf{u} \right)^{T} \right) $$

方程式の大前提でニュートン流体を仮定したので、$\nabla \cdot \sigma$ はニュートンの粘性法則に従って次のように書ける。 $$ \begin{align*} \nabla \cdot \sigma =& \nabla \cdot \left( - p I + \tau \right) \\ =& \nabla \cdot \left( - p I + \mu \left( \nabla \mathbf{u} + \left( \nabla \mathbf{u} \right)^{T} \right) \right) \\ =& - \nabla p + \mu \nabla \cdot \left( \nabla \mathbf{u} + \left( \nabla \mathbf{u} \right)^{T} \right) \\ =& - \nabla p + \mu \nabla \cdot \left( \nabla \mathbf{u} \right) + \mu \nabla \cdot \left( \nabla \mathbf{u} \right)^{T} \end{align*} $$

$\nabla \mathbf{u}$ のダイバージェンス: $$ \nabla \cdot \nabla \mathbf{u} = \nabla^{2} \mathbf{u} $$ $\left( \nabla \mathbf{u} \right)^{T}$ のダイバージェンス: $$ \nabla \cdot \left( \nabla \mathbf{u} \right)^{T} = \nabla \left( \nabla \cdot \mathbf{u} \right) $$

一方、流体は非圧縮性と仮定したので $\nabla \cdot \mathbf{u} = 0$ であり、最後の項を消すことができる。 $$ \begin{align*} \nabla \cdot \sigma =& - \nabla p + \mu \nabla^{2} \mathbf{u} + \mu \nabla \left( \nabla \cdot \mathbf{u} \right) \\ =& - \nabla p + \mu \nabla^{2} \mathbf{u} \end{align*} $$

これで元の方程式に戻ると、次のように二通りの表現を得る。 $$ \begin{align*} \rho {\frac{ D \mathbf{u} }{ D t }} =& \nabla \cdot \sigma + \rho \mathbf{g} \\ =& - \nabla p + \mu \nabla^{2} \mathbf{u} + \rho \mathbf{g} \\ \overset{\div \rho}{\implies} {\frac{ D \mathbf{u} }{ D t }} =& - \nabla w + \nu \nabla^{2} \mathbf{u} + \mathbf{g} \end{align*} $$

■