ナビエ-ストークス方程式

glass

これまで、流体の粘性を考慮しないで運動方程式を導いてきました。

ここでは、流体の粘性も考慮して流体における支配方程式を導出します。

粘性

固定された板の上に液体があり、更にその上に板を置きます。

そして、この上の板を水平に引っ張るのにある程度の力が必要です。

これが液体の粘性です。

粘性力

粘性力は板の表面積や引く速度に比例すると予測されます。よって

 \displaystyle F=\mu S\left(\frac{du}{dy}\right)

が固定された板にかかる力となります。

 \displaystyle\tau =\frac{F}{S}

のような単位面積当たりの力を応力(おうりょく)と呼び、\mu粘性係数と呼びます。

上記のように流体の速度勾配に比例して粘性力がかかる液体をニュートン流体と呼びます。

そして、ニュートン流体でない流体を非ニュートン流体と呼びます。

流体にかかる応力

では流体にかかる応力を計算していきます。

以下のように流体表面にかかる力を設定します(少し複雑ですが、規則性を見つけてください)。

設定

せん断力

そして、回転における運動方程式を立てると、

 \tau_{xy}-\tau_{yx}+\tau_{xy}-\tau_{yx}=0

 \tau_{yz}-\tau_{zy}+\tau_{yz}-\tau_{zy}=0

 \tau_{zx}-\tau_{xz}+\tau_{zx}-\tau_{xz}=0

となります。ここでは限りなく小さい微小要素を考えているので慣性項や体積項は無視され、結局モーメントの釣り合いの式になっています。

これを解くことで

 \tau_{xy}=\tau_{yx}

 \tau_{yz}=\tau_{zy}

 \tau_{zx}=\tau_{xz}

が得られます。

このことからわかるのですが、粘性力の説明で使用した状況では、y軸に平行な面にも同様にせん断応力がかかっていることになります。

せん断力

このことを一般化すれば

 \displaystyle\tau_{xy}=\tau_{yx}=\mu\left(\frac{\partial v}{\partial x}+\frac{\partial u}{\partial y}\right)・・・①

 \displaystyle\tau_{yz}=\tau_{zy}=\mu\left(\frac{\partial w}{\partial y}+\frac{\partial v}{\partial z}\right)

 \displaystyle\tau_{zx}=\tau_{xz}=\mu\left(\frac{\partial u}{\partial z}+\frac{\partial w}{\partial x}\right)

が得られます。

軸方向の力

軸方向の力は圧力pの効果と体積変化の効果があり

 \displaystyle\sigma_{xx}=-p+2\mu\frac{\partial u}{\partial x}+\lambda\Theta・・・②

 \displaystyle\sigma_{yy}=-p+2\mu\frac{\partial v}{\partial y}+\lambda\Theta

 \displaystyle\sigma_{zz}=-p+2\mu\frac{\partial w}{\partial z}+\lambda\Theta

となります。ここで、

 \displaystyle\Theta=\frac{\partial u}{\partial x}+\frac{\partial v}{\partial y}+\frac{\partial w}{\partial z}

であり、この項がゼロでないと体積膨張しているということになります。

非圧縮性流体の場合はゼロになります。

ここで圧力pは周りの流体から押される応力\sigma_{xx},\sigma_{yy},\sigma_{zz}の平均だと考えると

 \displaystyle p=-\frac{\sigma_{xx}+\sigma_{yy}+\sigma_{zz}}{3}

となり、\sigma_{xx},\sigma_{yy},\sigma_{zz}に式②を代入すると

 \displaystyle\lambda=-\frac{2}{3}\mu

が得られます。

ナビエ-ストークス方程式

では、流体の微小要素にかかる力がわかったので運動方程式を立てていきます。

設定

図を見ながらx成分の運動方程式をたてると

 \displaystyle\rho\Delta x\Delta y\Delta z\frac{Du}{Dt}=\rho\Delta x\Delta y\Delta z F_{x}\\+\left(\sigma_{xx}(x+\Delta x,y,z)-\sigma_{xx}(x,y,z)\right)\Delta y\Delta z\\+\left(\tau_{yz}(x,y+\Delta y,z)-\tau_{yz}(x,y,z)\right)\Delta z\Delta x\\+\left(\tau_{zx}(x,y,z+\Delta z)-\tau_{zx}(x,y,z)\right)\Delta x\Delta y

となります。テイラー展開で高次の微小項を無視して計算すると

 \displaystyle\rho\frac{Du}{Dt}=\rho F_{x}+\left(\frac{\partial\sigma_{xx}}{\partial x}+\frac{\partial\tau_{yx}}{\partial y}+\frac{\partial\tau_{zx}}{\partial z}\right)

が得られます。あとは式①と②を代入すれば

 \displaystyle\frac{Du}{Dt}=F_{x}-\frac{1}{\rho}\frac{\partial p}{\partial x}+\frac{1}{3}\frac{\mu}{\rho}\frac{\partial \Theta}{\partial x}+\frac{\mu}{\rho}\left(\frac{\partial^{2}u}{\partial x^{2}}+\frac{\partial^{2}u}{\partial y^{2}}+\frac{\partial^{2}u}{\partial z^{2}}\right)

が導かれます。同様に、y、z成分を計算すると

 \displaystyle\frac{Dv}{Dt}=F_{y}-\frac{1}{\rho}\frac{\partial p}{\partial y}+\frac{1}{3}\frac{\mu}{\rho}\frac{\partial \Theta}{\partial y}+\frac{\mu}{\rho}\left(\frac{\partial^{2}v}{\partial x^{2}}+\frac{\partial^{2}v}{\partial y^{2}}+\frac{\partial^{2}v}{\partial z^{2}}\right)

 \displaystyle\frac{Dw}{Dt}=F_{z}-\frac{1}{\rho}\frac{\partial p}{\partial z}+\frac{1}{3}\frac{\mu}{\rho}\frac{\partial \Theta}{\partial z}+\frac{\mu}{\rho}\left(\frac{\partial^{2}w}{\partial x^{2}}+\frac{\partial^{2}w}{\partial y^{2}}+\frac{\partial^{2}w}{\partial z^{2}}\right)

となります。これらの方程式をナビエ-ストークス方程式と呼び、流体力学の基礎方程式となります。

ベクトルを使って簡潔に書くと

 \displaystyle\frac{D\vec{v}}{Dt}=\vec{F}-\frac{1}{\rho}\nabla p+\frac{\mu}{3}\nabla\Theta+\nu\nabla^{2}\vec{v}

となります。ここで、

 \displaystyle\nu=\frac{\mu}{\rho}

動粘性係数といいます。

ちなみに非圧縮性流体の場合は、\Thetaがゼロなので

 \displaystyle\frac{D\vec{v}}{Dt}=\vec{F}-\frac{1}{\rho}\nabla p+\nu\nabla^{2}\vec{v}

となります。

著者:安井 真人(やすい まさと)