ナビエ–ストークス方程式の復習 — 二つの式、四つの力

非圧縮性流体の運動を決める連続の式と運動量方程式を一項ずつ読み解き、粘性係数とレイノルズ数を改めて定義する。

はじめに

第2章でテンソル表記を整えたので、いよいよその表記を用いて実際の流体がどう動くかを語る二本の方程式を書ける。本章を終えると、読者はナビエ–ストークス方程式の各項を指で押さえながら「これは加速度、これは圧力、これは粘性が生む摩擦」と日本語で説明できるようになるはずだ。式を丸暗記する必要はない — 各項がどんな物理的効果を意味するのかを覚えておく方が、ずっと長く残る。

本論 1 — 質量保存:非圧縮流体は漏れない

非圧縮の仮定の下では、流体はどこかで生成も消滅もしない。ある一点に入ってくる量と出ていく量がぴったり等しくなければならない。これをベクトル形式で書けば

u=0\nabla \cdot \vec{u} = 0

である。ここで u\vec{u}(u-ベクトル)は速度場、\nabla \cdot(発散)は「ある点での湧き出し量」を測る演算子だ。同じ式を第2章で導入したテンソル(添字)表記で書くと

uixi=0\frac{\partial u_i}{\partial x_i} = 0

となる。uiu_i(u-アイ)は速度の第 ii 成分、/xi\partial / \partial x_i(ラウンド x-アイ)は座標 xix_i に沿って変化する割合を表す。繰り返し添字 ii は総和規約により i=1,2,3i=1,2,3 について加算される。つまりこの一行は u1/x1+u2/x2+u3/x3=0\partial u_1/\partial x_1 + \partial u_2/\partial x_2 + \partial u_3/\partial x_3 = 0 の省略表記だ。

物理的には「どんな小さな箱を切り出しても、その箱に入る流れと出る流れの差は0」という意味になる。非圧縮の仮定が崩れる場合(例:超音速の空気)は右辺に密度変化率の項が加わるが、本書の大部分は水や低速空気を扱うので、右辺は最後まで0のままだ。

本論 2 — 運動量方程式:ニュートンの F=ma を流体に移す

運動量方程式は一行のニュートン第二法則である。非圧縮・一定粘性の流体に対して

uit+ujuixj=1ρpxi+ν2uixjxj\frac{\partial u_i}{\partial t} + u_j \frac{\partial u_i}{\partial x_j} = -\frac{1}{\rho}\frac{\partial p}{\partial x_i} + \nu \frac{\partial^2 u_i}{\partial x_j \partial x_j}

左辺は「流体粒子の加速度」、右辺は「その粒子に働く単位質量あたりの力」。一項ずつ読むと:

  • ui/t\partial u_i / \partial t時間微分項。空間上の一点にカメラを固定した時、その点での速度が時間とともに変わる割合。「定点加速度」と呼ぼう。
  • ujui/xju_j \, \partial u_i / \partial x_j移流項(convective term)uju_j 方向に流れる流体が、その方向の速度変化を運んできて一点に加速度として現れる効果。総和規約で j=1,2,3j=1,2,3 について加算される。非線形性の源であり — 次章以降ずっと — 乱流を難しくしているたった一つの項だ。
  • (1/ρ)p/xi-(1/\rho)\, \partial p / \partial x_i圧力勾配の力ρ\rho(ロー、密度、kg/m³)は単位体積あたりの質量、pp(ピー、圧力、Pa)は単位面積あたりの力。圧力の高い側から低い側へ流体を押し出す。負号がその向きを決める。
  • ν2ui/xjxj\nu \, \partial^2 u_i / \partial x_j \partial x_j粘性拡散項ν\nu(ニュー、動粘性係数、m²/s)はすぐ後で定義する。同じ添字 jj が二度現れるので総和規約で 2ui/x12+2ui/x22+2ui/x32\partial^2 u_i / \partial x_1^2 + \partial^2 u_i / \partial x_2^2 + \partial^2 u_i / \partial x_3^2。隣り合う流体層の摩擦で速度差をなめらかに削り取る効果。第1章で挙げた「拡散性」はまさにこの項から出てくる。

四つの項を日本語で覚え直すと — 定点加速度 + 移流 = -圧力 + 粘性拡散。この一行を手で一度書いてみると、これから出てくる RANS や LES の方程式が全部この式の変奏であることが見えてくる。

本論 3 — 粘性係数とレイノルズ数を読み直す

粘性には二つの表記がある。μ\mu(ミュー、粘性係数、Pa·s)はせん断応力と速度勾配の比例定数だ。しかし運動量方程式で粘性が加速度として入ってくる時、必ず μ/ρ\mu / \rho という形でしか現れないので、この比に名前を付けたのが動粘性係数

ν=μρ\nu = \frac{\mu}{\rho}

である。単位は (Pa·s) / (kg/m³) = m²/s。平たく言えば「動粘性係数が大きいほど、その流体はせん断(層と層の滑り)に強く抵抗する」。水(ν106m2/s\nu \approx 10^{-6}\,\mathrm{m^2/s})より空気(ν1.5×105m2/s\nu \approx 1.5\times 10^{-5}\,\mathrm{m^2/s})の方が大きいのは、分子が軽く同じ粘性効果を出すのにより少ない密度で済むためだ。

移流項と粘性項が同じ式に並んでいるということは、どちらが大きいかを一つの数で計れるということでもある。代表速度 UU(m/s)と代表長さ LL(m)を取ると、移流項の大きさは U2/LU^2/L、粘性項の大きさは νU/L2\nu U / L^2。両者の比が

Re=ULν\mathrm{Re} = \frac{U L}{\nu}

すなわちレイノルズ数である。単位を追うと (m/s)m/(m2/s)=1(\mathrm{m/s}) \cdot \mathrm{m} / (\mathrm{m^2/s}) = 1、確かに無次元だ。第1章で見た結論を改めて書いておくと — Re が大きければ移流項が粘性項を圧倒し、粘性が生んでいた「なめらかに削る効果」が移流の作る非線形なかき混ぜに追いつけなくなる。その結果が乱流である。

Pythonで確かめる

水20°Cの粘性データから ν\nu を求め、直径2 cmの管に毎秒1 mで水が流れるときの Re を計算し、第1章の表で分類してみる。

import numpy as np

# 水20°Cの粘性係数 mu と密度 rho
mu = 1.002e-3        # Pa·s
rho = 998.0          # kg/m^3

# 動粘性係数 nu = mu / rho
nu = mu / rho        # m^2/s

# 代表速度 U と代表長さ L (直径2 cmの管、歩く程度の速さ)
U = 1.0              # m/s
L = 0.02             # m
Re = U * L / nu      # 無次元

# 第1章の境界値で流れの状態を分類
if Re < 2300:
    state = "層流 (laminar)"
elif Re < 4000:
    state = "遷移 (transition)"
else:
    state = "乱流 (turbulent)"

print(f"mu  = {mu:.3e} Pa·s")
print(f"rho = {rho:.1f} kg/m^3")
print(f"nu  = {nu:.3e} m^2/s")
print(f"Re  = {Re:.0f}  ->  {state}")

出力は Re ≈ 19,900 — 4,000をはるかに超える明らかな乱流域だ。水道の蛇口から落ちる水流があるところで断面全体に散らばっていく様子こそが、この数の視覚的な証拠である。

次章へ

第4章: 渦度と渦度方程式では運動量方程式の回転(curl)を取って圧力項を落とし、第1章で「乱流の核心」と呼んだ渦度伸長が式の中でどう生き返るかを見る。本章の移流項 ujui/xju_j \, \partial u_i / \partial x_j が、そこで再び主役になる。