Vol. 24

ナビエ・ストークス方程式
会員

堀城之

序Dv
1
1＝ F − gradp ＋ νgradθ＋νΔv
Dtρ
3
これが，ナビエ・ストークス方程式である。主に流体力学において用いられる。土木工学，建築工学，機械工学，航空工学，船舶工学，物理学を履修した方は学んだかもしれない。流体に係る発明も少なくないであろう。しかし，大学院の流体研究室或いは，当該数学研究室に入らない限り圧縮項（右辺第 3 項）までは導かないと思われる。そこで，本稿ではナビエ・ストークス方程式について圧縮項まで含めて説明する。なお，特別な境界条件の場合を除いて未だに誰も解けていない方程式である。

１．ナビエ・ストークス方程式の生みの親

ここで，問題です。感覚的に分かってもらえばよい

土木技術者であったアンリ・ナビエが粘性を考慮しので，細かい条件等は抜きにして大雑把に書きます。

た方程式を記載した論文を 1822 年にフランス学士院

Q：一様に定常的な¢流れるプール£にあなたが飛び込んでプールの底に立つとする。あなたは水平方向に提出した。にどのような力を受けますか？

しかし，技術者の論文は無視されるのはいつの時代でも世の常で，提出から 20 年以上経過した 1845 年にあなたはこう答えるでしょう。

物理数学者のストークスがナビエとは別個に粘性流体

A：下流に流される力を受ける。流れが強いと下流に流される。

方程式を導いた。

答えは×です。

そこで，二人の名前をつけてナビエ・ストークス方
程式とした。ナビエとストークスとの間に¢・£が付

オイラーの運動方程式に支配される流れでは，あないているのはその為である。第一発見者の名前をとったは下流方向への力を受けない。無風状態で打たれたてナビエ方程式の方が正しいような気がする。因みにゴルフボールは，重力以外に力を受けず空気中を飛んナビエは，その後，エリート養成学校であるエコーでいくのである。
非常に奇妙な話である。あなたも流れるプールや川ル・ポリテクニークの教授になっている。

で流されたかもしれない。流体中に存在するものは上流から下流に流されるように力を受けることを感覚的に理解している。ゴルフをやる人なら流体抵抗が無かったらどんなに楽なことかと思うかもしれない。
アンリ・ナビエ

しかし，オイラーの運動方程式で扱う流体では下流

from Wikipedia

２．ナビエ・ストークス方程式が導かれた理由

方向に流される力を受けないのである。オイラーの運動方程式に支配される流体を完全流体という。

―ダランベールの背理―
右辺第 3 項（圧縮項）及び第 4 項（粘性項）を除い

そして，このパラドックスをダランベールの背理という。

た式は，オイラーの運動方程式である。この方程式

当該パラドックスを解くために考え出されたのが，は，ナビエ・ストークス方程式の約 100 年前に，あのナビエ・ストークス方程式なのである。ナビエは，土大数学者オイラーによって，ニュートンの第 2 法則か木橋梁技術者であり，橋脚に受ける力を研究した結ら導かれた。

果，ナビエ・ストークス方程式を導き出した。

Vol. 66

No. 11−1−

パテント 2013

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

w(x, y, z, t)＋δwナビエ・ストークス方程式自体を導くことはそれほど難しくはない。高校の数学，物理の知識があれば導くことが出来る。
（１） オイラーの運動方程式
まず，ナビエ・ストークス方程式の圧縮項（筆者が勝手に名付けた。学部・学科で名称が異なるようである）及び粘性項（右辺第 4 項）が無い方程式を導く。
当該方程式をオイラーの運動方程式という。
``xDv
1
Dt ＝ F − ρ gradp~

grad ＝``y``z
オイラーの運動方程式は，右辺は単位質量当たりの力，左辺は質量力と圧力である。

2 法則F=Mcから導くことが出来る。
x ＋δx，y ＋δy，z ＋δz
図１の速度成分は，場所と時間の関数，
u ＝ u（x, y, z, t）
v＝v（x, y, z, t）
w＝w（x, y, z, t）

dv ∂v
∂v
∂v
＝dt ∂x u ＋ ∂y v ＋ ∂z w
dv ∂w
∂w
∂w
dt ＝ ∂x u ＋ ∂y v ＋ ∂z w ⑦
左辺を見れば分かるように，所謂，加速度である。
∂u∂tを時間的加速度，∂u∂u∂u∂x u ＋ ∂y v ＋ ∂z w
を場所的加速度，という。

①

代表される。
δt 時間経過後の移動距離は，
δx ＝ u(x, y, z, t)δt
δy ＝ v(x, y, z, t)δt
δz ＝ w(x, y, z, t)δt⑤を多変数テイラー展開（※）すると，∂∂
u(x, y, z, t)＋δu ＝ u(x, y, z, t)＋( δx ＋
δy ＋
∂x
∂y∂
∂
∂z δz ＋ ∂t δt)u(x, y, z, t）∂∂∂∂＋( δx ＋δy ＋
δz ＋ δt)2 u(x, y, z, t)＋ 3
∂x
∂y∂z
∂t
乗＋…，…＋ n 乗∂∂v(x, y, z, t)＋δv ＝ v(x, y, z, t)＋( δx ＋
δy ＋
∂x
∂y∂
∂
δz ＋ δt)u(x, y, z, t）
∂z
∂t
∂
∂
∂
∂
＋( δx ＋
δy ＋ δz ＋ δt)2 u(x, y, z, t）＋ 3
∂x
∂y
∂z
∂t
乗＋…，…＋ n 乗∂∂
w(x, y, z, t)＋δw ＝ w(x, y, z, t)＋( δx ＋
δy ＋
∂x
∂y∂
∂
δz ＋ ∂t δt)u(x, y, z, t）
∂z
∂
∂
∂
∂
＋( δx ＋
δy ＋
δz ＋ δt)2 u(x, y, z, t)＋ 3
∂x
∂y∂z
∂t
乗＋…＋ n 乗
影響が小さい高次項を無視し，②を代入すると，∂u∂uu(x, y, z, t)＋δu ＝ u(x, y, z, t)＋( ∂x uδt ＋ ∂y vδt ＋∂u∂u∂z wδt ＋ ∂t δt)∂v∂vv(x, y, z, t)＋δu ＝ v(x, y, z, t)＋(uδt ＋vδt ＋∂x∂y∂u∂v∂z wδt ＋ ∂t δt)∂w∂w
w(x, y, z, t)＋δw ＝ w(x, y, z, t)＋( ∂x uδt ＋ ∂y vδt∂w∂w＋ ∂z wδt ＋ ∂t δt)⑥
δulim
δt → 0 δtとして⑥を整理すると，
du ∂u∂u∂u＝u＋v＋w
dt∂x∂y∂zゆえに，オイラーの運動方程式は，ニュートンの第図 1 に示す，流体中における微少流体

⑤②

座標は，
x ＋δx

前者は理解できると思うが，後者がなぜ加速度と思y ＋δyうかもしれない。つまり場所で微分してなぜ加速度な③のかと思うかもしれない。場所的加速度は，ある場所③を①に代入すると，から他の場所に移動したときの速度の変化率（加速z ＋δzu ＝ u(x ＋δx, y ＋δy, z ＋δz, t ＋δt)度）だと考えて頂ければよい。以上で加速度が求めらv ＝ v(x ＋δx, y ＋δy, z ＋δz, t ＋δt)

れた。つまり，オイラー方程式の左辺である。

w ＝ w(x ＋δx, y ＋δy, z ＋δz, t ＋δt)

次に，右辺，すなわち，力について考える。

④

完全流体に加わる力は，質量力と圧力のみである。

微少時間経過後の速度の変化量をδu として，④を以下のように書き換える。

質量力は，流体工学の本に出てくる用語なのだが，筆u(x, y, z, t)＋δu

者が訳すと慣性力となる。
流体中に微少な立方体を考える。

v(x, y, z, t)＋δv

パテント 2013−2−

Vol. 66

No. 11①

粘性係数μ剪断係数τとの関係は，∂v ⎞τ＝μ⎛⎝ ∂x ⎠で表される。これをニュートンの摩擦法則という（ニュートンの粘性法則ともいう）。
②応力テンソル

図２
図 2 に示す立方体の質量は，ρδxδyδz。である。
Fxv

各面に加わる単位重量当たりの力を F ＝ FyvFxFzとすると，質量力＝ Fρδxδyδz ＝ Fy ρδxδyδz⑧Fzとなる。

次に，微少流体の中心を（x，y，z）とすると，微少流体に作用する圧力は，例えば x 軸に垂直な面に作用図３する圧力差であるから，∂P∂Pdx{(x ＋ dx , y, z)−(x −, y, z)}δyδz∂x2∂x2∂P＝− ∂x δxδyδz ⑨

σは着目面に垂直に作用する応力（材料力学で言えば面外応力），τは着目面に作用する剪断応力（同面内他の 2 面に作用する圧力も同様に求めると，∂P∂P−∂y δxδyδz，− ∂z δxδyδz ⑩，⑪⑦〜⑨をニュートンの第 2 法則に代入し，質量で除応力）である。最初の添え字は着目面を，2 番目のそせば
∂u
∂u
∂u
∂u
1 ∂p
∂x u ＋ ∂y v ＋ ∂z w ＋ ∂t ＝ Fx− ρ ∂x∂v
1 ∂p∂v∂v∂v−u＋∂x∂y v ＋ ∂z w ＋ ∂t ＝ Fy ρ ∂y∂w
∂w
∂w
∂w
1 ∂p
∂x u ＋ ∂y v ＋ ∂z w ＋ ∂t ＝ Fz− ρ ∂zD∂∂∂∂
微分演算子 Dt (＝ ∂x u ＋ ∂y v ＋ ∂z w ＋ ∂t )を用いて書き換えると，
Du
1 ∂p
Dt ＝ Fx− ρ ∂x

この図をしっかり覚えて頂ければ終わったようなもれは作用方向を示す。
以上を行列で表したものが応力テンソル（略）である。のです。
着面に作用する剪断応力は，係る 2 つの剪断応力の和で表されるので，ニュートンの粘性法則から，∂v ∂u ⎞τxy＝τyx＝μ⎛＋⎝ ∂x ∂y ⎠∂w ∂v ⎞τyz＝τzy＝μ⎛⎝ ∂y ＋ ∂z ⎠∂u ∂w ⎞τxy＝τyz＝μ⎛⑫⎝ ∂z ＋ ∂x ⎠また，着面に作用する面外応力は，圧力と，係る 2Dv
1 ∂p
Dt ＝ Fy− ρ ∂yつの剪断応力の和と，圧力による歪みの項との和で表Dw
1 ∂p
Dt ＝ Fz− ρ ∂zgrad を用いて記載し直せば，
Dv
1
Dt ＝ F − ρ gradp
これで，オイラーの運動方程式が導けた。

されるので，∂u
σxx＝−p ＋ 2μ∂x ＋λθ∂v
σyy＝−p ＋ 2μ∂y ＋λθ∂w
σzz＝−p ＋ 2μ ∂z ＋λθ（２） 粘性流体における基礎方程式
1 μ―圧縮項 3 ρ gradθ，粘性項νΔv のオイラー運動方程式への導入―

μを第 1 粘性係数，λを第 2 粘性係数という。
∂u ∂v ∂w ⎞θ＝ div v ＝⎛⎝ ∂x ＋ ∂y ＋ ∂z ⎠非圧縮性流体では歪み 0 故，θ＝ 0 となる。

ダランベールの背理が生じるのは粘性を考慮していないからである。また，圧縮性も考慮されていない。
因みに，水は非圧縮性流体である。

Vol. 66

No. 11⑬

次いで，式⑫，⑬を使って，図 3 に示す立方体に働く力を考える。

−3−

パテント 2013

すると，x 軸方向の¢力£は対面の差で表されるから∞
f(x, y, z)＝Σn=0

x 面では
∂σxx
δx・δyδz∂xy 面では∂τyx
δy・δxδz∂yZ 面では∂τzx
δz・δxδy∂zであり，これらの総和を単位質量で表すと，
1 ⎛ ∂σxx ∂τyx ∂τzx ⎞⑭ρ ⎝ ∂x ＋ ∂y ＋ ∂z ⎠となる。

＝ f(a, b, c)
1∂∂∂⎞
1＋ ⎛hf(a, b, c)
1! ⎝ ∂x ＋ k ∂y ＋ j ∂z ⎠1⎛ ∂∂∂⎞2＋h＋k＋jf(a, b, c)
2! ⎝ ∂x∂y∂z ⎠1⎛ ∂∂∂⎞3＋h＋k＋jf(a, b, c)
3! ⎝ ∂x∂y∂z ⎠＋…
1∂∂∂⎞n＋ ⎛h＋k＋jn! ⎝ ∂x∂y∂z ⎠ f(a, b, c)＋ Rn＋11∂∂⎞⎛ ∂n＋1
Rn＋1 ＝h＋k＋jf(a ＋θh,(n−1)! ⎝ ∂x∂y∂z ⎠b ＋θk, c ＋θj)
0 ＜θ＜ 1⑭を 3．（1）で求めたオイラーの運動方程式に代入する
Du1 ∂σxx ∂τyx ∂τzx ⎞＝ Fx− ⎛Dtρ ⎝ ∂x ＋ ∂y ＋ ∂z ⎠となる。

⑮：多変数ラグランジュ剰余項となる。

圧力項が消えていると思われるかもしれないが，式ここで，⑬にあるように，σの中に入っている。
式⑮に式⑫，⑬の当該式を代入してやれば，
Du
1 ∂p1 ∂θ ⎛ ∂2u ∂2u ∂2u ⎞＝
F＋＋＋
νx−Dt
ρ ∂x
3 ∂x ＋ν⎝ ∂x2 ∂y2 ∂z2 ⎠μ∵ ν＝ ρ
ν：動粘性係数1⎛ ∂∂∂⎞nh＋k＋jn! ⎝ ∂x∂y∂z ⎠ f(a, b, c)

x ＝ a ＋ ht，y ＝ b ＋ kt，z ＝ c ＋ jt
つまり，流線上の 1 点 x, y, z から微少時間δt 経過（微小変化）した場合（故に 1 次）の予測値（座標）を表している（f(x＋δx, y＋δy, z＋δz）＝ f(x, y, z)＋δf(x, y, z)）。

y 面，z 面も同様に求めると，
Dv
1 ∂p 1 ∂θ ⎛ ∂2u ∂2u ∂2u ⎞＝ Fx−Dtρ ∂y ＋ 3ν∂y ＋ν⎝ ∂x2 ＋ ∂y2 ＋ ∂z2 ⎠Dw
1 ∂p 1 ∂θ ⎛ ∂2u ∂2u ∂2u ⎞＋＋＝ Fx−＋ ν ＋ν
⎝ ∂x2 ∂y2 ∂z2 ⎠Dtρ ∂z3 ∂z

―導出―
関数 f(x, y, z)の座標 f(a, b, c)における x, y, z 軸方向となる。
これらの式を grad，Δを用いて表すと
Du11
Dt ＝ F − ρ gradp ＋ 3νgradθ＋νΔv
Δ＝∇2＝∂2∂2∂22＋2＋∂x∂y∂z2'∂∂∂∇＝ ∂x ＋ ∂y ＋ ∂zの傾きに単位時間当たりの各軸方向における変化量をそれぞれ掛けてやれば変化量δf(x, y, z)が出てくるから∂∂∂
δf(x, y, z)≡ dW ＝⎛h＋k＋j ⎞dt⎝ ∂x∂y∂z⎠∂∂∴ f(x, y, z)＋δf(x, y, z)＝ f(a, b, c)＋⎛⎝h ∂x ＋k ∂y ＋∂j ∂z⎞⎠f(a, b, c)

非圧縮性流体では，Du
1＝ F − gradp ＋νΔvDtρとなる。
※多変数テイラー展開
1f(x m)＝ 6 [(x m−Y m)` m] nf(Y m)from Wikipedia
n=0 n
イギリスの数学者ブルック・テイラーが導入した級数式である。ある関数の 1 点から少し離れたところの値を級数を使って真値に近づけるものである。高校の授業で 1 変数テイラー展開（マクローリン展開（a ＝0））を使って円周率や三角関数を単純な数の和として⎛ h ∂ ＋ k ∂ ＋ j ∂ ⎞ の偏微分係数も変化するの⎝ ∂x∂y∂z ⎠で，それを考慮するとδt 時間経過後の偏微分係数の変化量は，
222⎛∂f(x, y, z)⎞ ⎛∂ f(x, y, z)⎞ ⎛∂ f(x, y, z)⎞δt＋⎛∂ f(x, y, z)⎞δtδ⎝ ∂x ⎠＝⎝δt＋2⎠ ⎝ ∂x∂y ⎠ ⎝ ∂x∂z ⎠∂x

222⎛∂f(x, y, z)⎞ ⎛∂ f(x, y, z)⎞ ⎛∂ f(x, y, z)⎞δt＋⎛∂ f(x, y, z)⎞δtδ⎝ ∂y ⎠＝⎝δt＋2⎠ ⎝ ∂x∂y ⎠ ⎝ ∂y∂z ⎠∂y222f(x,y,z)∂f(x,y,z)∂⎞⎛⎞δt＋⎛∂ f(x, y, z)⎞δt＋⎛∂ f(x, y, z)⎞δtδ⎛＝⎝ ∂z ⎠ ⎝ ∂z2 ⎠ ⎝ ∂z∂y ⎠ ⎝ ∂x∂z ⎠求めた方もいるかもしれない。

δt
δt を二つに分割すると，前半の 2 後の変化量は，工学，物理で用いられる便利な式である。重要なので 3 次元の場合を導く。

1⎛ ∂∂
∂⎞
δt
11! ⎝h ∂x ＋ k ∂y ＋ j ∂z ⎠ f(a, b, c) 2

ナビエ・ストークス方程式は 3 次元なので 3 変数場

δt
後半の 2 後の変化量は，合は，関数 f(x, y, z)を座標 a, b, c 周りでテイラー展開

パテント 2013−4−

Vol. 66

No. 11⎡ 1 ⎛h ∂ ＋k ∂ ＋j ∂⎞1 ⎛∂f(x, y, z)⎞∂yz ⎠＋δ⎝ ∂x⎣ 1!⎝ ∂x⎠

つ非圧縮性であれば粘性項及び圧縮項を無視できるので 1 階線形微分方程式となり解くことは可能である。

∂f(x, y, z)⎞ ⎛∂f(x, y, z)⎞⎤ δt＋δ⎛＋δ
f(a, b, c)＝⎝ ∂y⎠ ⎝ ∂z⎠⎦ 2
∂ ∂⎞1 ⎛∂2f(x, y, z)⎞⎡1⎛ ∂h＋k＋j ⎠＋⎝⎣ 1! ⎝ ∂x
⎠
∂yz
∂x2
2
∂ f(x, y, z)⎞ ⎛∂2f(x, y, z)⎞ ⎛∂2f(x, y, z)⎞⎛＋＋＋⎝ ∂x∂y ⎠ ⎝ ∂x∂z ⎠ ⎝ ∂y2⎠∂2f(x, y, z)⎞ ⎛∂2f(x, y, z)⎞ ⎛∂2f(x, y, z)⎞⎛＋＋＋⎝ ∂x∂y ⎠ ⎝ ∂y∂z ⎠ ⎝ ∂z2⎠
22f(x,y,z)∂f(x,y,z)∂⎛⎞ ⎛⎞⎤ δt＋⎝ ∂z∂y ⎠＋⎝ ∂x∂z ⎠⎦ 2 f(a, b, c)

しかし，上記の如き特定の場合を除き，数学的には，厳密解が求められていない。求めることができるのかも分かっていない。
そこで，CMI（クレイ数学研究所）は賞金を懸けている。
http://www.claymath.org/millennium/
CMI は，解けないという証明をしても賞金をもらえδt 時間経過後の変化量は，前半と後半との和で表されるから，

1
∂∂
∂ 1
δt
δtδf(x, y, z)＝ 1!⎛h ∂x ＋k ∂y ＋j ∂z⎞f(a, b, c)⎝⎠
21⎛ ∂∂∂⎞
1 δtf(a, b, c)＋ ⎝h＋k＋j
1! ∂x∂y ∂z⎠ 222∂ f(x, y, z)⎞ ⎛∂ f(x, y, z)⎞ ⎛∂2f(x, y, z)⎞＋⎡⎛＋＋⎣⎝ ∂x2⎠ ⎝ ∂y2⎠ ⎝ ∂z2⎠22⎛∂ f(x, y, z)⎞ ⎛∂ f(x, y, z)⎞＋ 2⎝ ∂x∂y ⎠＋ 2⎝ ∂y∂z ⎠2⎛∂ f(x, y, z)⎞⎤ δt＋ 2⎝ ∂x∂z ⎠⎦ 2 f(a, b, c)
1
∂∂
∂ ⎞1
δtδf(x, y, z)＝ ⎛f(a, b, c)δth1!⎝ ∂x ＋k ∂y ＋j ∂z⎠⎡⎛∂f(x, y, z)⎞ ⎛∂f(x, y, z)⎞＋⎣⎝ ∂x⎠＋⎝ ∂y⎠∂f(x, y, z)⎞⎤δt＋⎛f(a,b,c)⎝ ∂z⎠⎦ 21⎛ ∂∂∂ ⎞1∴δf(x, y, z)＝ ⎝hf(a, b, c)
1! ∂x ＋k ∂y ＋j ∂z⎠1∂∂∂ ⎞2＋ 2!⎛⎝h ∂x ＋k ∂y ＋j ∂z⎠ f(a, b, c)∂⎛ ∂∴ f(x, y, z)＋δf(x, y, z)＝ f(a, b, c)＋⎝h ∂x ＋k ∂y ＋∂1⎛ ∂∂∂ ⎞2j ∂z⎞
⎠ f(a, b, c)＋ 2!⎝h ∂x ＋k ∂y ＋j ∂z⎠ f(a, b, c)

るそうである。
ところで，CMI は，2006 年にポアンカレ予想を証明できたとして，グリゴリー・ヤコヴレヴィチ・ペレル
マンに賞金を出した（本人は受け取り拒否）。先日
NHK で「数学者はキノコ狩りの夢を見る」という題名で放送していたので見た方もいるかもしれない。なお，見逃した人は，オンデマンドで見ることが出来る。
http://www.nhk-on δ eman δ .jp/
ポアンカレ予想とは，位相幾何学の問題で，紐の基端を手で持ち先端をロケットに結びつけて飛ばし，宇宙の隅々まで飛行し地球に帰還した後，紐の基端と先端とをたぐり寄せると紐の全てを回収できるという予想である（星はツルツルで紐は星に引っかからない。
紐は燃えない，切れない）。ペレルマンの論文は，例えば，Cornel University Library でダウンロードすることが出来る。
http://arxiv.org/pδf/math/0211159v1.pδfn ＝ 0, 1, 2 についての導出終わり。高次項は無視さhttp://arxiv.org/pδf/math/0303109v1.pδfれるのでこれで打ち止め。

http://arxiv.org/pδf/math/0307245v1.pδf

なお，ラグランジュ剰余項を含めて，ロルの定理等

ペレルマンは，ウィリアム・サーストンの幾何化予でテイラー級数を導入することが出来る。
Fin

想を，微分幾何学に物理工学的アプローチも加えて解４．賞金＄100 万ドル

き，ポアンカレ予想を解明した。サーストンの幾何化予想とは，宇宙では最大 8 つの図形に分類されるとい

ナビエ・ストークス方程式は，上記の通り，2 階非線うものである。

形偏微分方程式である。
以前，パテントに掲載した振動方程式は 2 階非線形

ナビエ・ストークス方程式も純粋に数学的アプロー偏微分方程式であり（線材は 4 階），ナビエ・ストーク

チではなく，他のアプローチにより解けるかもしれな

ス方程式と同じである。断面変化等の場合を除き，時い。たとえば，生みの親であるナビエが土木技術者だ

間微分のみである。ナビエ・ストークス方程式では時から土木工学的アプローチも面白そうである。土木工

間微分と場所微分とが混在している。よって，振動方学的に，流体を類型化し，それぞれ圧縮項，粘性項が程式よりも解きにくい。

0 とならない境界条件を定めてプロットしていくと，解法のヒントが見つかるかもしれない。

特定の条件下で境界条件を設定して解くことは可能以上である。振動方程式では減衰を無視すれば 2 階線形偏参考書：日野幹雄著「流体力学」株式会社朝倉書店

微分方程式となり解くことが出来る。ナビエ・ストー本間

クス方程式では，例えばマッハのような高速流体で且Vol. 66

No. 11−5−

仁著「標準水理学」

丸善株式会社

パテント 2013