流体力学II教材讲解

流体力学II（Viscous Fluid and Gas Dynamics）讲义第一章、粘性不可压缩流体运动基本方程组（学时数：6）1-1．绪论
流体力学是力学的一个重要分支，主要研究流体介质（液体、气体、等离子体）的特性、状态，在各种力的作用下发生的对流、扩散、旋涡、波动现象和质量、动量、能量传输，以及同化学、生物等其他运动形式之间的相互作用。

它既是一门经典学科，又是一门现代学科，对自然科学和工程技术具有先导作用。

历史上，力学包括流体力学，曾经经历基于直观实践经验的古代力学、基于严密数学理论的经典力学、基于物理洞察能力的近代力学三个阶段。

在人类早期的生产活动过程中，力学即与数学、天文学一起发展。

17世纪，Newton基于前人的天文观测和力学实验，发明了微积分，并总结出机械运动三大定律和万有引力定律，发表了著名的《自然哲学的数学原理》一书。

由于原理是普适自然与工程领域的规律，从而使力学成为自然科学的先导。

从17世纪开始，人们逐步建立了流体力学的基本理论体系，从Pascal定律、Newton粘性定律、Pitot 管测速，到Euler方程和Bernoulli方程，标志着流体动力学正式成为力学的一个分支学科。

18世纪，人们着重发展无粘流体的位势理论。

到了19世纪，为了解决工程实际问题，开始注重粘性的影响，Navier-Stokes方程的建立为流体力学的进一步发展奠定了完整的理论基础，但该方程解的存在性与光滑性的证明至今仍是一大难题。

20世纪初，Prandtl凭借出色的物理洞察能力，提出边界层理论，从而开创了流体力学的近代发展阶段，使力学成为人类实现“飞天”梦想的重要理论先导。

60年代以来，由于超级计算机、先进测试技术的发展和应用，力学进一步凸显宏微观结合和学科交叉的特征，进入现代力学发展新阶段。

刚刚过去的2011年，人类遭遇了一系列极端事件：日本海底地震导致海啸和福岛核电站泄露事故；澳大利亚飓风；我国干旱洪水灾害等异常气候问题。

这些极端事件的预测预警都是流体力学的前沿问题。

同一年，美国航天飞机历经30年130多次飞行之后终于宣布全面退役，其中一个重要原因是存在防热系统不可靠的安全隐患，这也是流体力学工作者亟待解决的一个重要课题。

因此，现代流体力学不仅是一门重要的基础学科，而且在航空航天、海洋海岸、环境能源、生物医学、材料信息等诸多与国家经济、社会发展密切相关的工程技术领域里，仍然具有不可或缺的先导作用。

本课程的教学目的是：在流体力学I的基础上，针对各种与科学技术发展和人类生活、经济活动紧密相关的流体力学问题，主要是粘性流体和气体动力学问题，建立相应的数学模型，综合运用数学、力学和数值模拟方法对模型进行求解，并根据求解结果解释实验现象，进而认识客观事物的规律性。

粘性流体和气体动力学知识的应用面广，涉及问题非常复杂，是力学专业高年级学生具备解决实际问题能力的重要一环，并且有利于形成新的学科交叉型思维。

主要内容：粘性流体和气体动力学的各种模型方程，Navier-Stokes方程的精确解和近似解，量纲分析方法，层流边界层理论，激波现象，湍流初步。

基本要求：掌握粘性流体和气体动力学的基本理论知识和主要物理概念，掌握综合运用数学、力学和数值模拟方法对模型进行求解并对实验现象进行解释的基本方法，能运用所学知识分析和解决一些实际问题。

主要参考书：
1.周光炯等著《流体力学》（第二版），高等教育出版社，2000年6月
2.叶敬棠等著《流体力学》，复旦大学出版社，1989年5月
3.[美]W.F.休斯等著《流体动力学》，科学出版社，2002年3月
4.[德]H.欧特尔等著《普朗特流体力学基础》（第11版），科学出版社，2008年6月
5.林建忠等著《流体力学》，清华大学出版社，2005年9月
6.[苏]谢多夫著《力学中的相似方法与量纲理论》，科学出版社，1982年12月
7.[美]F.M.怀特著《粘性流体动力学》，机械工业出版社，1982年12月
1-2．Navier-Stokes 方程组的导出
在单相单组分连续介质、各向同性Newton 流体、热力学过程为准静态过程的前提下，流体运动的基本方程组可写成
0d V dt ρ
ρ+∇= （1-1） 1
[]dV F P dt ρ
=+∇ （1-2）
211
()([])()2R d V e F V P V k T q dt ρρ
+=+∇+∇∇+ （1-3）
2
[]2[]()[][][]3
P S p V I p I μμτ=-+∇=-+ （1-4）
p R T ρ=， v e C T = （1-5）
1[]()()2j
i ij j i
u u S s x x ∂∂==+∂∂ （1-6）
(),[]()12(),3j
i j i ij i i u u x x u V x μττμ∂⎧∂+⎪
∂∂⎪
==⎨∂⎪-∇⎪∂⎩
i j i j ≠= （1-7）
3
200110()110
T T T T μμ+=+ （1-8） 上列方程组中，自变量是t 和x （123{,,}T x x x x =），F 和R q 是单位质量流体上的体积力和辐射热，k 是热传导系数，μ是动力粘性系数，R 是气体普适常数，v C 是等容比热，[]P 是应力张量，[]S 是应变率张量，[]τ是粘性应力张量，[]I 是单位张量，未知量有密度ρ、速度V （123{,,}T V u u u =）、内能e 、温度T 和压强p ，微分算子123
{,,}T
x x x ∂∂∂∇=∂∂∂，方程组是封闭的。

如果假定流体为不可压缩，则上述方程组可以简化为
0V ∇= （1-9） 21dV F p V dt μ
ρρ=-∇+∇
（1-10） 2R dT C k T q dt
ρρ=Φ+∇+ （1-11）
方程组（1-9），（1-10）就称为Navier-Stokes 方程组，其中只包含未知量V 和p ，可以与（1-11）解耦。

能量方程（1-11）是温度T 的控制方程，其中Φ称为粘性耗散函数，它与粘性应力
张量[]τ以及速度场V 之间的关系为
([])([])V V ττΦ=∇-∇ （1-12） []()T
V V τμ=∇+∇ （1-13）
这里V ∇是并矢张量，T
V ∇是它的转置，即
()(
)j ij i
u V V x ∂∇=∇=∂，
()(
)T T i
ij
j
u V V x ∂∇=∇=∂ （1-14）
直角坐标系中Navier-Stokes 方程组的分量形式为
2222222
2
2
222222222
01()1()1()x y z u v w x y z
u u u u p u u u u v w F t x y z x x y z v v v v p v v v u v w F t x y z y x y z w w w w p w w w u v w F t x y z z x y z μρρμρρμρρ∂∂∂++=∂∂∂∂∂∂∂∂∂∂∂+++=-+++∂∂∂∂∂∂∂∂∂∂∂∂∂∂∂∂+++=-+++∂∂∂∂∂∂∂∂∂∂∂∂∂∂∂∂+++=-+++∂∂∂∂∂∂∂∂ （1-15）
进一步假定为定常流动，体积力为重力，则沿着流线有变形的Bernoulli 方程
2
2222122111
1122p p V gz V gz V dl μ
ρρρ++=+++∇⎰ （1-16）
1-3．Navier-Stokes 方程组的定解条件提法
（1）初始条件：0(,,,0)(,,)V x y z V x y z =，0(,,,0)(,,)p x y z p x y z =；
（2）自由面条件：一般仍用理想流体自由面上的运动学条件和动力学条件，也有假定
自由面形状为已知而要增加切向应力平衡的情况； （3）无穷远条件：V
V ∞∞
=，p p ∞∞=；
（4）固壁边界条件：W W
V V =，固壁静止时0W
V
=，比理想流体多了一个切向速
度限制条件，所以Euler 方程或Laplace 方程都不能适定求解，无旋性也无法保持。

在固壁边界处经常要计算无量纲的壁面摩擦系数f C 和热量吸收系数h ，分别定义为
212
W
f C V τρ∞=
， W s f s f
T
k
q n h T T T T ∂-∂==-- （1-17） 其中W τ和W q 分别为单位面积上的壁面摩擦阻力和热流通量，s T 和f T 分别为壁面附近固壁和流体的温度，n 是从流体指向固体的固壁法线方向。

1-4．柱面坐标系r z O θ和球面坐标系R O θλ中Navier-Stokes 方程组的分量形式
（1）速度场的散度
1()r r z
r z V V V V V r r r z θθθ∂∂∂∇=+++∂∂∂ （1-18） 11()2cot sin R R R V V V V V V R R R R R θθλθλ
θθθλ
∂∂∂∇=++++∂∂∂ （1-19）
（2）迁移加速度
2
(())(())r r z r V V V V V e r θθ∇=∇-+
(())(())r z z V V
V V e V V e r
θθθ+∇++∇ （1-20）
其中r z V V V V r r z
θθ∂∂∂∇=++∂∂∂。

222cot (())(())(())R R R R V V V V V V V V V e V V e R R θλθλθθλθθ
+-∇=∇-+∇++
cot (())R V V V V V V e r
λθλλλθ
++∇+ （1-21）
其中sin R V V V V R R R θλθθλ
∂∂∂
∇=++∂∂∂。

（3）压强梯度
1()(,,)r z p p p p r r z θθ∂∂∂∇=∂∂∂ （1-22） 11()(,,)sin R p p p p R R R θλ
θθλ
∂∂∂∇=∂∂∂ （1-23） （4）粘性项
22222()()r r r z r V V V V e r r θ
θθ∂∇=∇-
-+∂ 2
2222()r z z V V V e V e r r θθθθ
∂+∇-++∇∂ （1-24） 其中2222
222
211r r r r z
θ∂∂∂∂∇=+++∂∂∂∂。

22221()((cot ))sin R R R R V V V V V V e R θλ
θλθθθθλ
∂∂∇=∇-++++∂∂
22
2
22((sin cos ))sin 2R V V V V e R θλθθθθθθλ
∂∂+∇+--+∂∂ 2
22
2((sin cos ))sin 2
R V V V V e R θλλλθθθλλ∂∂+∇++-∂∂ （1-25） 其中2222
2222222
21cot 1sin R R R R R R θθθθλ∂∂∂∂∂∇=++++∂∂∂∂∂。

（5）速度场的梯度
111()r
z r r
z r z
r z V V V r r r V V V V V V r r r r r V V V z z
z θ
θθθθθθθ∂∂∂⎡⎤⎢⎥
∂∂∂⎢
⎥∂∂∂⎢⎥
∇=-+⎢⎥∂∂∂⎢
⎥∂∂∂⎢⎥⎢⎥∂∂∂⎣
⎦
（1-26）
111()sin cos sin cos sin sin sin R R R
R R R V V V R R R V V V V V V R R R R R V V V V V V V R R R θ
λθ
θλθλ
θ
λ
λ
λθθθθ
θ
θθθλλ
λ
θθ
θ⎡⎤
∂∂⎢⎥
∂⎢⎥∂∂∂⎢
⎥
∂∂∂⎢⎥∇=-+⎢⎥∂∂∂⎢
⎥
∂∂∂⎢⎥
--++⎢⎥∂∂∂⎢⎥⎣
⎦
（1-27）
（6）粘性耗散函数
22
21[2(
)2()2()r r z r z V V V V r r r z
θθμθ∂∂∂Φ=++++∂∂∂ 222
11()()()]r z r z V V V V V V V r r r r z z r θθθθθ∂∂∂∂∂∂+-+++++∂∂∂∂∂∂ （1-28）
22211[2()2()2(cot )sin R R R R V V V V V V R R R R R R θλθ
θλμθθθλ∂∂∂Φ=++++++∂∂∂
22
111()(cot )sin R V V V V V V R R R R R R θθλλθθθθθλ∂∂∂∂+-+
+-++∂∂∂∂ 2
1()]sin R V V V R R R λλθλ∂∂++-∂∂ （1-29）
（7）热流通量
(,,)r z T k T T
q k k r r z θθ∂∂∂=---∂∂∂ （1-30） (,,)sin R T k T k T q k R R R θλ
θθλ
∂∂∂=---∂∂∂ （1-31）
1-5．涡输运方程
定义速度场的涡量V Ω=∇⨯，将Navier-Stokes 方程组写成Lamb 形式
21()2V V V F p t μ
ρρ
∂+∇+Ω⨯=-∇-∇⨯Ω∂ （1-32）
在体积力有势的条件下取旋度得
2()d V dt μ
ρ
Ω=Ω∇+∇Ω （1-33）
由于方程（1-33）相比于理想流体的Helmholdz 方程右端多了一个粘性扩散项，旋涡将
在流体中逐渐扩散，所以一般粘性流体运动都是有旋的。

不过，在一些特殊的条件下，Navier-Stokes 方程组也可以存在无旋解，例如： （1）整个流场以匀速运动，速度到处相同；
（2）一个半径为a 的无限长直圆柱绕着对称轴以角速度ω旋转，引起周围流体的旋转运动，就相当于一根环量为2
2a πωΓ=的无限长直涡线诱导的理想不可压缩流动（满足同样的边界条件和无穷远条件），所以是无旋的。

如果这个圆柱突然停止旋转，则在其邻近由于粘性而产生速度梯度和旋涡，并且将逐渐向流场中扩散，下面我们来求解这个涡量场。

采用柱面坐标系，根据问题特征有
V V e θθ=， z z e Ω=Ω， r e r
∂
∇=
∂， 0V ⇒⋅∇=， 0Ω⋅∇= 涡输运方程（1-33）简化为
221()z z z
t r r r
μρ∂Ω∂Ω∂Ω=+∂∂∂ （1-34）
定解条件则有
(,0)0z r Ω=，
(,)0z t Ω∞= （1-35） 解析解为
24(,)r t
z A r t e
t
ρμ-
Ω=， 2
2a A ρωμ
= （1-36）
图1.1 涡输运方程解析解的示意图
从图1.1中可以看到，流场中任一确定的r 位置处，涡量都从初始时刻的零值开始，先有一个比较快速的增加过程，达到最大值后则开始逐渐减小，经过很长时间重新趋向于零。

随着r 增大（离开圆柱较远），涡量的最大值变小，并且到达最大值的时间也变长。

这个例子描述了一个典型的在流体粘性作用下旋涡扩散过程。

1-6．方程组的守恒形式
对一般可压缩流动方程组有
()0V t ρ
ρ∂+∇=∂ （1-37） ()
([][])V VV p I F t ρρτρ∂+∇+-=∂ （1-38） (()[])s
s R E V E p V k T F V q t τρρ∂+∇+--∇=+∂ （1-39） 其中，2
()2
s V E e ρ=+，VV 是并矢张量，定义为()()ij i j VV VV u u ==。

对不可压缩流动则可简化为
0V ∇= （1-40） 1
([][])V p VV I F t τρρ
∂+∇+-=∂ （1-41）
在二维直角坐标系中分量形式可写成
()()
()u F u G u H u t x y
∂∂∂++=∂∂∂ （1-42）
其中
{0,,}T u u v = （1-43）
{,2,()}T p u v u
F u uu vu x x y μμρρρ∂∂∂=+--+∂∂∂ （1-44）
{,(),2}T
u v p v G v uv vv y x y μμρρρ∂∂∂=-++-∂∂∂ （1-45）
{0,,}T x y H F F = （1-46）
1-7．量纲分析方法与∏定理
（1）用长度[]L 、时间[]T 、质量[]M 和温度[]θ四个基本量纲表示的常用物理量量纲列出如下：速度1[]LT -，密度3[]ML -，体积力2[]LT -，应力（压强、壁面剪应力）12[]ML T --，运动粘度21[]L T -，动力粘度11[]ML T --，应变率（角速度）1[]T -，
内能22[]L T -，气体常数（比热）221[]L T θ--，热传导系数31[]MLT θ--，耗散函数13[]ML T --，辐射热23[]L T -，壁面热流通量3[]MT -，热量吸收系数31[]MT θ--。

（2）∏定理：独立的无量纲参数个数=物理量个数-基本量纲个数。

（3）在单相单组分连续介质、各向同性Newton 流体、热力学过程为准静态过程的前提下，流体运动的基本方程组中特征物理量有ρ∞、U ∞、p ∞、T ∞、L 、α、0t 、g 、μ、
k 、R 、p C 、v C 共13个，所以应有9个独立的无量纲参数。

无量纲参数的取法并不是唯
一的，一种常见的形式是
0t L
S U t ∞=， α， 2
r U F gL
∞=， 2
u p E U ρ∞∞∞=， e U L R ρμ∞∞=， p r C P k μ=， 2
c p U E C T ∞∞=， p v C C γ=，
a M = （1-47）
其中γ和a M 在可压缩流动问题中使用，r P 和c E 则在考虑热力相似时使用。

另外还有两个无量纲参数也比较常见
p e e r C U L
P R P k
ρ∞∞==
， 3
2
(1)u e a
U L
N P M kT ργ∞∞∞
=-= （1-48）
它们分别称为Peclet 数和Nussel 数。

1-8．方程组的无量纲化（不考虑热辐射）
对于粘性不可压缩流体方程组，设特征物理量为U ∞，L ，ρ，T ∞，g ，μ，k ，C ，
0L t U ∞=和2
p U ρ∞∞=，引入无量纲参数2
r U F gL
∞=， e U L R ρμ∞=， r C P k μ=，
2
c U E CT ∞
∞
=以后得到 *
*
0V ∇= （1-49）
*
*******
*2*
11()r e V V V F p V t F R ∂+∇=-∇+∇∂ （1-50） ****
*2***
11()()c e r
T V T T E t R P ∂+∇=∇+Φ∂ （1-51） 其中*
x x L =
， *0t t t =， *V V U ∞
=， *F F g =， *p p p ∞=， *T T T ∞=，*22/U L μ∞ΦΦ=为无量纲变量。

对于粘性可压缩流体方程组，在不考虑体积力的条件下，又可设特征物理量为U ∞，L ，
ρ∞，T ∞，μ，p C ，v C ，R ，k ，0L
t U ∞
=
，2p U ρ∞∞∞=和2e U ∞∞=，引入无量纲参数e U L R ρμ∞∞=， 3
u U L
N kT ρ∞∞∞
=以后得到
()0V t ρ
ρ∂+∇=∂ （1-52）
()1
()[]e V VV p t R ρρτ∂+∇=-∇+∇∂ （1-53）
211
(())([])s s u e
E V E p T V t N R τ∂+∇+=∇+∇∂ （1-54）
其中变量已经全部无量纲化，为了书写简单起见把“*”号省略了。

第二章、Navier-Stokes 方程的精确解（学时数：6）
2-1．Navier-Stokes 方程的精确解、近似解和数值解概述
粘性流体运动的基本方程，即Navier-Stokes 方程的建立为流体力学的发展和实际应用奠定了完整的理论基础，但该方程解的存在性与光滑性的证明至今仍是一大难题。

不仅如此，从1844年Stokes 导出Navier-Stokes 方程算起，已经过去了170多年，但已知的该方程的精确解平均还不到一年一个，说明求解该方程存在着数学上的极大困难。

20世纪60年代以后，伴随着超级计算机的发展和应用，计算流体力学（Computational Fluid Dynamics ，CFD ）迅速崛起，采用各种离散逼近方法求得Navier-Stokes 方程的数值解占了流体力学研究的半壁江山。

除了精确解、数值解之外，还有一种方法叫近似解，也就是凭借物理的洞察能力，求解前先对Navier-Stokes 方程进行各种简化，然后求得简化后方程的解析解。

有些文献中对精确解和近似解不加区分，只要采用的是解析推导方法而不是数值离散方法得到的解都称为精确解或解析解。

本课程则对这两者提出明显的区别：精确解的表达式能够满足完整的Navier-Stokes 方程，近似解则不可能。

这一章将讨论Navier-Stokes 方程能够找到精确解的几种典型情况：Couette 类型的流动，Poiseuille 类型的流动，不定常的Stokes 流动和旋转坐标系中的Ekman 流动。

第三章和第四章将分别讨论小Reynolds 数和大Reynolds 数条件下Navier-Stokes 方程的近似解。

至于Navier-Stokes 方程的数值解法，将在《计算力学》等课程中再介绍，本课程中只给出几个简单的例子作为说明。

2-2．Couette 类型的定常流动
对于理想流体来说，如果固体壁面在其自身所占空间范围里切向移动或旋转，都不会引起流体的运动。

粘性流体则不一样，由于粘性作用，紧贴着固壁的那层流体会跟着固壁移动或旋转，然后它又带动相邻的第二、第三…层流体移动或旋转，并且由于存在切向的速度差，流体微团会发生旋转，产生的涡量从壁面逐步向流体内部扩散。

这种由于固体壁面切向运动引起的粘性流动一般就称为Couette 类型的流动，它又可以分为平面固壁平移、圆柱体的旋转和圆柱体的平移等几种情况。

（1）平面问题：假设位于0y =处的无限大平板以常速度1u 沿着x 轴正向平移，同时位于y h =处的无限大平板则以常速度2u 沿着x 轴正向平移，两平板间充满动力粘度为μ的不可压缩流体，要求流体速度分布。

根据问题的物理特征，可以判断（注意：这种判断是求Navier-Stokes 方程精确解的关键，实质是寻找一种特殊形式的精确解）0v w ==，以及
0t z
∂∂==∂∂，未知函数是(,)u u x y =和(,)p p x y =。

又根据连续性方程
0u
x
∂=∂ （2-1） 和y 方向动量方程
0p
y
∂=∂ （2-2） 可得()u y 的控制方程为
22120
1,y y h d u dp
const
dy dx u u u u μ==⎧==⎪⎨
⎪==⎩ （2-3） 121()()(1)2dp y y
u y y y h u u dx h h
μ⇒=-+-+
（2-4） 不失一般性令10u =， 222h dp
P u dx
μ=-，得
2(1)y y
u u P P h h
=+- （2-5）
图2.1 平面Couette 流动速度剖面
在板面恒温且不考虑压力梯度的条件下又有
22
212
212
0(),y h y d T u u dy kh T T T T μ==⎧=--⎪⎨
⎪==⎩
（2-6） 2
21122
()()()(1)2y y T y u u y y h T T kh h h
μ⇒=---+-+ （2-7） 相对静止或无粘流体而言，由于粘性作用，流体中的温度升高，并且升高的幅度随着无量纲参数2
2
r c u P E kh μ∆=
的增大而增大。

（2）旋转问题：假设半径1r r =的无限长圆柱面以常角速度1ω绕着对称轴旋转，同时半径2r r =的无限长圆柱面则以常角速度2ω绕着对称轴旋转，两圆柱同心，其间充满动力粘度为μ的不可压缩流体，要求流体速度分布。

根据问题的物理特征，取圆柱坐标系r z O θ，可以判断0r z V V ==，以及
0t z
θ∂∂∂
===∂∂∂，未知函数是()V V r θθ=和()p p r =。

此时连续性方程自动满足，()V r θ的控制方程为
12
22
2112210,r r r r d V dV V dr r dr r V r V r θθθ
θθωω==⎧+-=⎪⎨⎪==⎩ （2-8） 1()V r Ar B r θ⇒=+， 22221122
21
r r A r r ωω-=-， 22
12122221()r r B r r ωω-=- （2-9） 进一步可得()p r 的控制方程为
21V dp
r dr
θρ=
（2-10） 2
2202()(4ln )2B p r A r AB r p r
ρ⇒=-++ （2-11）
（3）圆柱粘度计原理：将（2-9）代入切向应力公式
()r V d r
dr r
θ
θτμ= （2-12） 2211222
21
2W r r r ωωτμ-⇒=-， 221
2122212W r r r ωωτμ-=-（方向相反！） （2-13） 所以在不考虑自身转动惯量的条件下，维持这样两个圆柱旋转所需的外力矩
22
2212111
21220021
4()L
W r r M r d dz L r r π
τθπμωω==-=-⎰
⎰
222220
L
W r d dz M π
τθ=
=⎰⎰
（2-14）
其中L 为圆柱长度，两个力矩大小相等方向相反，在12ωω=时外力矩为零。

若已知L 、1r 、2r 、ω，则M 就是μ的线性函数，测出M 就能推算出粘度μ，这就是圆柱粘度计测量原理。

（4）旋转流体稳定性分析初步：考虑某一0r z V V ==的旋转流体（不一定是旋转圆柱之间流场），设想由于某种扰动使原来位于1r 处的流体质点改变到了2r （不失一般性
21r r >）处，角动量守恒，所以速度从1V θ改变到了
1
12
r V r θ。

在2r 处当地压力梯度是2222
r r V dp
dr r θρ
==，如果该压力梯度大于为了平衡离心力所需的压力梯度，那么流体质点趋向回到原来的位置1r 处，此时流体是稳定的。

因此可得流体的稳定条件是
221122
()r r r dp V dr r r θρ=> ⇔ 22
2211()()rV rV θθ> （2-15）
定义流体的单位环量C V r θ=，则稳定性判据也可写成0d C
dr
>。

将此判据应用于旋转圆柱之间流场，因为2
C Ar B =+，
2dC
Ar dr
=，在内圆柱1r 处可任取10ω>和0C >，外圆柱2r 处2ω可正可负，则稳定性判据为
222211r r ωω> （2-16）
反之就是不稳定的。

这个结果是Rayleign 首先得到的，Taylor 进一步给出了考虑速度分量
有轴对称的类波扰动时更加精确的解。

在不稳定范围里，流动将会在子午面上发展成具有一系列圆环形涡胞的另一种定常状态，该状态取决于12/ωω值。

（5）其他典型Couette 类型流动有：两个同心圆柱相对平移而引起的流动（习题8-17）；有固定分界面的两层流体的流动（习题8-11）；有固定自由面或斜面的Couette 流（习题8-13）；有壁面渗透作用的Couette 流（习题8-12）。

2-3．Poiseuille 类型的定常流动
Couette 类型的流动产生的原因是固体壁面切向运动，压强梯度则不重要，可以考虑也可以不考虑。

Poiseuille 类型的流动不一样，它就是在压强梯度作用下管道内的流动，最简单的例子就是打开水龙头时管道内水的流动。

因为这类问题最终归结为线性的Poisson 方程的齐次边值问题，所以可以求得精确解的情况（管道形状）相对较多，构成了Navier-Stokes 方程精确解的很大一部分。

本节主要介绍平面管、圆管以及另外几种可以写成简单初等函数解的管道，通过其中无限长直圆管的解可以得到著名的Hagen-Poiseuille 公式，正是这一公式奠定了Navier-Stokes 方程合理性的基础。

（1）平面问题：与平面Couette 流动基本相同，只是两块平板均为静止不动。

所以同样可以判断0v w ==，得到()u y 的控制方程为
220
10,0y y h d u dp
const
dy dx u u μ==⎧==⎪⎨
⎪==⎩ （2-17） 1()()2dp
u y y y h dx
μ⇒=- （2-18）
()2xy u dp h
y y dx τμ∂==-∂， 002W xy y h dp dx ττ===->(0)u > （2-19） 3012h
v h dp Q udy dx μ==-⎰， 212v Q h dp
u h dx
μ==-
， max 32u u = （2-20） 同样假定板面恒温，又有
222
22120()()4,y h y dp d T h dx y hy dy k T T T T μ==⎧
⎪⎪=--+⎨⎪
==⎪⎩
（2-21） 2
322312()()(1)(243)24dp y y
dx T y T T y y hy h y h h h k μ
⇒=-+--+- （2-22）
这时温度分布是y 的四次函数。

若又设12T T =，则
2
30()024W y dp T
dx q k h y μ
=∂=-=-<∂ （2-23）
说明流体内部温度总是升高，向壁面放出热量。

温度升高幅度与r c P E 有关，在同样的r c P E 下单纯Poiseuille 流动的升高幅度总比单纯Couette 流动大。

对水或空气而言，这些温度升
高都还比较小，可以忽略不计，粘性系数μ可视为常数，油则不然。

（2）无限长等截面直圆管内流动问题：将圆管的对称轴取为z 轴，则可以判断在直角坐标系中有0u v ==，在圆柱坐标系中有0r V V θ==。

（2-1）截面为(,)0s x y =时控制方程的一般形式是
2222(,)010s x y w w dp
const
x y dz w μ=⎧∂∂+==⎪∂∂⎨
⎪=⎩
（2-24） 或者
22222
(,)001110,z z z z z s r r V V V dp
const r
r r r dz V V θθμ==⎧∂∂∂++==⎪∂∂∂⎨⎪=<∞⎩
（2-25） （2-2）对于半径为a 的等截面圆管来说，则有
220110,z z z z r a r V V dp
const r
r r dz V V μ==⎧∂∂+
==⎪∂∂⎨⎪=<∞⎩
（2-26） 22
1()()4z dp V r a r dz
μ⇒=-
- （2-27）
4028a
v z a dp
Q V rdr dz
ππμ=⋅=-
⎰ （2-28） 以及
28z a dp V dz μ=-， 2max 24z a dp
V V dz μ=-=
（与平面问题作比较） （2-29） 2z dV r dp dr dz τμ==， 02W rz r a a dp dz ττ==-=->(0)z V > （2-30）
定义阻力系数212W f z
C V τρ=， Reynolds 数z ed V d
R ρμ=（2d a =），又有
16
f ed
C R = （2-31） 最后，同样假定板面恒温，又有温度场方程为
2222
01()()()4,z r a w r V T r dp T r r r r k k r k dz T T T μμ==⎧∂∂∂Φ-∇==-=-=⎪∂∂∂⎨
⎪<∞=⎩
（2-32） 442()()64W a r dp T r T k dz μ-⇒=+， 32
()016W r a T a dp q k r dz
μ=∂==-<∂ （2-33）
同样说明粘性流动使流体温度升高，向壁面放热。

（2-3）Hagen-Poiseuille 公式：根据（2-28）式可知管道中的体积流量、压差、流体动力粘度和管道半径之间成立公式
44088L
v p p a dp a Q dz L
ππμμ-=-=
（2-34） 这里L 为管长，0p 和L p 分别为其入口端和出口端压强。

Hagen （1839）和Poiseuille （1841）分别在实验中发现流量与管径的四次方成正比的规律，Stokes （1845）和Hagenbach （1858）
给出了上述理论解，与实验结果比较后成为N-S 方程的重要理论基础。

（2-34）称为Hagen-Poiseuille 公式，它说明，在同样的压强梯度下，直径大一倍的粗管与四根面积相同但直径小一半的细管相比，前者的流量是后者的四倍。

利用Hagen-Poiseuille 公式可得修正Bernoulli 方程
22000122L L
L V p V p gz gz h ρρ
++=+++ （2-35）
中沿程阻力引起的压头损失
201420
88L
L
z v z
p p L L
h V dz Q V a a νννρ
π-=-∇⋅=
=
=⎰ （2-36）
毛细管粘度计原理：如图2.2所示，毛细管上接一个大容器，下端开口连通大气，大容器内液体的速度可忽略不计。

因此，毛细管两端压差0()L p p g H L ρ-=+，代入Hagen-Poiseuille 公式得
4
(1)8v
a H
g Q L
πμρ=
+
（2-37）
图2.2 毛细管粘度计示意图
（3）其他几种典型截面直管内的速度和流量公式
（3-1）椭圆管：22
22(,)10x y s x y a b =+-=
2222
22
221(,)(1)2dp a b x y w x y dz a b a b
μ=--+ （2-37） 33
22
4v dp a b Q dz a b
πμ=-+ （2-38） （3-2）圆环管：(,)()()0s r r a r b θ=--=（a b <，习题8-18）
222ln
ln
1()()4ln ln z r b
dp a r
V r a b r b b dz a a
μ=-+- （2-39）
222
44()[]8ln v dp b a Q b a b dz a
πμ-=--- （2-40）
（3-3
）正三角形管：2
2
(,)(3)()0
s x y y x y =--
= 22(,)3)()w x y y x y =-
（2-41）
4320v dp
Q dz
μ=-
（2-42）
2-4．非定常流动--Stokes 问题
（1）物理模型：在平面Couette 流动中，位于0y =处的无限大平板以常速度沿着x 轴正向平移，带动其上方的粘性不可压缩流体作平行流动。

现在进一步考虑平板从静止开始的起动过程，流体也从静止开始起动，由于在有限的时间内只有一定厚度的流体被带动，所以可假定平板上方是半无限大流场，这样一个非定常问题就称为Stokes 第一问题。

与此相类似，如果假定位于0y =处的无限大平板以固定频率沿着x 轴来回震荡，带动其上方半无限大流场中的粘性不可压缩流体也作平行震荡流动，这就称为Stokes 第二问题。

对于这两种问题，都可以判断0v w ==，
0z
∂
=∂，又根据连续性方程知u 与x 无关。

（2）Stokes 第一问题的控制方程为（压力梯度仅为时间的函数，这里不考虑）
220
00,t y y u u
t y u u u u ν=→∞=⎧∂∂=⎪∂∂⎨
⎪===⎩ （2-43） 设22000(,)()()'yu u y y f f f u u u t t
ηνν====（相似变换），可得 22
''
20'1,'0d u du d d u u ηηηηη=→∞⎧+=⎪⎨
⎪==⎩
（2-44） 2
00
(,)(1)u y t u d ττ-⇒=
（2-45）
2
4y t z u
y ν-∂Ω=-=∂
（2-46）
这是一个与图1.1所示相同的涡扩散过程，由
00.01u
u = ⇒
2
=
⇒涡扩散距离 4δ====
（2-47） （3）Stokes 第二问题的控制方程为
22000,cos y y u u
t y u u u t νω→∞
=⎧∂∂=⎪∂∂⎨
⎪==⎩ （2-48）
设
120
()cos ()sin u
f y t f y t u ωω=+，可得 122
11122''0,''0(0)1,()(0)()0
f f f f f f f f ωωνν⎧-=+=⎪⎨
⎪=∞==∞=⎩ （2-49）
0(,)cos()u y t u e
t ω⇒= （2-50）
同样由00.01u
u =
，可以看到这时涡扩散距离为~
δ 下面介绍几种直接数值求解（2-48）的差分格式，计算中的初值由（2-50）给出。

1110002
000(2),11,01cos ,0cos 2,00,0n n n n n j j j j j j y j n n J u u s u u u j J n N u u e j y j J u u ns y n N u n N ++--∆⎧=+-+≤≤-≤≤-⎪=∆≤≤⎪⎨=∆≤≤⎪
⎪=≤≤⎩
（2-51） （2-51）称为FTCS 格式，其中(,)(,)n n
j j u u y t u j y n t =∆∆， 2
t s y ν∆=∆，
0y y ∆=。

我们采用100J =， 5000N =， 01u =， 00.05y ∆=， 113
,,424
s =进行数值试验，
所得结果如图2.3所示。

图2.3 Stokes 第二问题的数值解
FTCS 格式是显式推进的格式，知道了前一时刻的函数值n
j u 就能直接计算后一时刻函数值
1
n j
u +，不必求解代数方程组，但其数值稳定性不易得到保证，从图2.3可以看到在11,42
s =
时能够稳定，但在3
4s =
时就已经失稳。

（2-51）中的内点推进格式可以改成 111111(2),11
,01n n n n n
j j j j j u s u u u u j J n N +++++---+=≤≤-≤≤- （2-52） 或者
11111111(2)(2),11,0122
n n n n n n n n
j j j j j j j j s s u u u u u u u u j J n N +++++-+---+=+-+≤≤-≤≤-（2-53）
（2-52）和（2-53）分别称为BTCS 格式和Crank-Nicolson （C-N ）格式，它们都是隐式差
分格式，从n j u 到1
n j u +的推进过程需要求解联立代数方程组，但能无条件数值稳定。

2-5．Ekman 流动--随地球旋转的运动坐标系
（1）Ekman 发现，在北冰洋极地附近，伴随着海面上大冰块的平移，其下方的海水将会出现一种旋转流动，不同深度处的水流方向是不同的，速度大小则随着水深衰减。

这种流动被称为Ekman 流动，其主要成因是流体的粘性和地球自转，所以能够用旋转坐标系中Navier-Stokes 方程的精确解来描述。

旋转坐标系中Navier-Stokes 方程形式为
2'1(')2''r r r d V f r V p V dt μ
ωωωρρ
=-⨯⨯-⨯-∇+∇
（2-53） 其中{,,}T r V u v w =是相对速度分布，它仍应满足连续性方程（1-9）。

（2）物理模型：不考虑相对加速度、体积力和压强梯度，又由于2'r V r ω⨯而忽略
牵连加速度。

根据Ekman 流动特征，选取直角坐标系（以下均省略'记号），其中z 轴垂直
海平面向上，海平面位于0z =处。

可以判断：0w =，
0x y
∂∂
==∂∂ ⇒ ()u u z =， ()v v z =。

（3）控制方程形式为
00''20,''20(0),(0),()()0
u v v u u u v v u v ωωνν⎧
+=-=⎪
⎨
⎪==-∞=-∞=⎩ （2-54） 令W u iv =+， 000W u iv =+，可得
22
020(0),()0
d W i W dz W W W ω
ν⎧-=⎪⎨
⎪=-∞=⎩
（2-55）
00()()u z u v ⇒=- （2-56）
00()(sin )v z u v =+ （2-57）
（4）对Ekman 现象的解释：410.710sec ω--≈⨯， 621
10sec m ν--≈，速度方向旋
转一个周期约等于海面向下0.7m 水深，所以能够明显观测得到。

第三章、Navier-Stokes 方程的近似解（学时数：6）
3-1．小球在无界粘性不可压缩流体中匀速缓慢运动问题的Stokes 解法
（1）Reynolds 数e U L
R ρμ
∞∞=
是粘性流动中最重要的无量纲参数，根据它的大小不同可以形成对Navier-Stokes 方程作近似处理的几种思路：当1e R 时，忽略方程中非线性
的对流项，简化为类似热传导方程的Stokes 方程求解，必要时引入部分对流项进行修正；当1e R 时，忽略方程中的粘性扩散项，简化为理想流体的Euler 方程，并采用Prandtl 边界层理论进行求解。

本章讨论1e R 条件下Navier-Stokes 方程的近似解。

所谓1e
R 条件下的流动，应该满足1U ∞（特征速度很小）
，或1L （流场特征尺度很小）。

显然，由于小球在无界粘性不可压缩流体中作匀速缓慢直线运动而引起的流体
运动，能够满足小e R 数流动的条件。

在小e R 数条件下，Navier-Stokes 方程可简化为Stokes 方程
0V ∇= （3-1）
21V F p V t μ
ρρ
∂=-∇+∇∂ （3-2）
流动定常并且不考虑体积力时，（3-2）式进一步简化为
21
V p μ
∇=
∇ （3-3）
在球面坐标系R O θλ中，方程（3-1）和（3-3）的分量形式为
112cot 0sin R R V V V V V R R R R R θθλ
θθθλ
∂∂∂++++=∂∂∂ （3-4） 22211(cot )sin R R V V p
V V V R R
θλθθθθλμ∂∂∂∇-+++=∂∂∂ （3-5）
222221(sin cos )sin 2R V V V p
V R R θλθθθθθλμθ
∂∂∂∇+--=∂∂∂ （3-6）
22221(sin cos )sin 2sin R V V V p
V R R θλλθθθλλμθλ
∂∂∂∇++-=∂∂∂ （3-7）
其中222
2
2222222
21cot 1sin R R R R R R θθθθλ∂∂∂∂∂∇=++++∂∂∂∂∂。

（2）控制方程与定解条件提法：取R O θλ为原点在球心处、跟随小球一起匀速直线运动的球面坐标系，对称轴（0θ=）的正向与小球运动方向一致，流体相对该坐标系的速
度设为{,,}T R V V V V θλ=。

根据问题的对称性可得0V λ=， 0λ
∂
=∂，所以方程组（3-4）~（3-7）又可简化为
12cot 0R R V V V V R R R R
θθ
θθ∂∂+++=∂∂ （3-8） 2222222121cot 2
(cot )R R R R R
V V V V V p V V R R R R R R R θθθθμθθθ
∂∂∂∂∂∂=+++-++∂∂∂∂∂∂ （3-9） 222222222121cot 2sin R
V V V V V V p R R R R R R R R θθθθθθμθθθθθ
∂∂∂∂∂∂=++++-
∂∂∂∂∂∂ （3-10）
这就是未知函数(,)R V R θ、(,)V R θθ和(,)p R θ的控制方程。

定解条件如下
0R R a R a
V V θ
====，
0cos R R V u θ→∞=-， 0sin R V u θθ→∞
= （3-11）
其中a 为小球半径，0u 为小球的速度。

（3）相对速度与压强分布求解：设()cos R V f R θ=⋅， ()sin V g R θθ=-⋅，
0()cos p h R p μθ=⋅⋅+，代入控制方程可得
2'()0f f g R +-= '2
R g f f =+
224''''()0h f f f g R R --+-= ⇒ 4
''''
h f f R
=+ 2122
'''()0h g g f g R R R
----= ''2''h Rg g f =+- 32''''8'''8''8'0R f R f Rf f ⇒++-=
21234311
()f R C C R C C R R ⇒=+++
23412311
()222C C g R C C R R R =++-
2321
()10h R C R C R
=+
利用定解条件（3-11）可知
10C u =-， 20C =， 3032C au =
， 3401
2
C a u =- 3
031[1()]cos 22R a a V u R R θ⇒=--
+⋅ （3-12） 3
031[1()]sin 44a a V u R R θθ=--⋅ （3-13）
02
3cos 2a p u p R μθ∞=⋅+ （3-14）
相对坐标系下绝对速度则是
3
00313'cos ()cos 2R R a a V V u u R R θθ=+=-⋅ （3-15）
3
00313'sin ()sin 4a a V V u u R R
θθθθ=-=-+⋅ （3-16）
（4）相对坐标系和绝对坐标系中的流动图案：在一个子午面内，仍以小球运动方向为
x 轴正向建立静止坐标系xy O ，则有2220()R x u t y =-+， 1
0tan y
x u t
θ-=-。

所以，绝
对速度分量可以写为
33
033933'cos 'sin [()()cos 2]8x R u a a a a V V V R R R R θθθθ=-=-+- （3-17）
3
0333'sin 'cos ()sin 28y R u a a V V V R R
θθθθ=+=- （3-18）
据此可以作出瞬时流线图并与理想流体有势流动比较（图3.1，上半部分相对静止坐标系，下半部分相对运动坐标系），其中无粘流体的相对速度分布为
30[1()]cos R a
V u R θ=--⋅ （3-19）
301[1()]sin 2a
V u R
θθ=+⋅ （3-20）
绝对速度分布则为
3
03(13cos 2)4x u a V R θ=+ （3-21）
3
033sin 24y a V u R θ= （3-22）
所以，303
2
02x u a V R π
θ==-<，而根据Stokes 流则有3032
3()04x u a a V R R πθ==+>，静止坐标系中两者方向恰好相反！在运动坐标系中看，主要差别则在于流线疏密（速度大小）不同。

Stokes 流动 有势流动
图3.1 小球运动的Stokes 流动与势流瞬时流线图比较
（5）Stokes 球阻公式及其应用：因为
2022932()cos 2R RR V a a
p p p u R R R μμθ∞∂=-+=---⋅∂
3
04
13()sin 2R R V V V a p u R R R R
θθθμμθθ∂∂=+-=⋅∂∂ 所以积分可得整个小球匀速缓慢运动中所受流体阻力为
(cos sin )RR
R R a F p
p ds θθθ==
-⎰⎰
2000
32(cos )sin 62p u a d u a a
π
μ
π
θθθπμ∞=-+
=-⎰
（3-23） 这就是著名的Stokes 球阻公式。

如果固体小球改为粘度为0μ的球形液滴，并假设运动过程中不发生变形（习题8-26），则有
00
21361F u a
μμπμμ
μ+
=-+ （3-24）
特别，气泡在液体中运动，0μμ，所以，04F u a πμ=-。

Stokes 球阻公式对物体形状不太敏感，如一个半径为a 的薄圆盘正面运动时阻力016F u a μ=-，侧面运动时则032
3
F u a μ=-。

根据Stokes 球阻公式可以得到液滴在空气中下落时速度满足的方程
331144
()()6323
du a a g au dt ρπρπρρπμ+=-- （3-25） 其中1ρ和ρ分别为液滴和空气密度，速度u 以向下为正。

假设初始速度为零，则有
2
1
92()2
()[1]t
a t u t u e
μρρ-
+=⋅- （3-26）
式中终极速度
212()9t t g
u u a ρρμ
→∞==
- （3-27）
在一个大气压和20C ︒条件下，31998/kg m ρ=， 31.2/kg m ρ=， 51.810a P s μ-=⨯⋅，
10a m μ=的雾滴0.012/t u m s =，在0.005
7t s =时0.99t u u =，在0.011t s =时0.9999t u u =。

所以，液滴下落过程中不定常阶段只占总时间中的极小比例（10a m μ=的雾滴降落一公里需要将近一天），通常可以忽略不计。

又由于下降速度与a 成正比，所以大
雾滴能够追上小雾滴并形成更大的雾滴，最后发展为降雨。

若已知a 、1ρ、ρ，则根据终极速度公式（3-27）可知μ就是t u 的反比例函数，测出
t u 就能推算出粘度μ，这就是应用于测量石油粘度的落球式粘度计的工作原理。

（6）其他小e R 数假定下Stokes 方程的典型解还有：小球匀速缓慢旋转而引起的流动
（习题8-23）；两块半无限大非平行平板间辐射流动（习题8-24）；球壳匀速缓慢旋转而引起内部流体的流动；两同心球壳以不同角速度缓慢旋转而引起的流动。

3-2．小球在无界粘性不可压缩流体中匀速缓慢运动问题的Ossen 修正
（1）Stokes 解法的局限性：Stokes 流动中有μ，a ，0u ，F 四个特征物理量，组成唯一的独立无量纲参数
0F
u a
μ，Stokes 公式则刻画了这一无量纲参数的相似律。

如果考虑无限长细圆柱匀速缓慢运动时所受流体阻力这样一个二维问题，则特征物理量仍是这四个，但由于此时F 的量纲单位是2[]MT -而不是2[]MLT -，所以相似律变成
F
u μ为常数，也就是说受力与圆柱大小无关！这称为Stokes 悖论，说明Stokes 解误差随着距离增大而增大。

（2）Ossen 修正思想：把相对坐标下的速度看作均匀流与小球引起的扰动的叠加（注意并不是在Stokes 解法的基础上叠加扰动，否则无法求解，称为Whitehead 悖论），保留一阶小量得到控制方程组，定解条件还是和Stokes 解法时一样。

在Ossen 修正下控制方程为
2
0'01''()'V V p V V ρμμ⎧∇⋅=⎪
⎨∇=∇+⋅∇⎪⎩
（3-28） 其中'{,,0}T R V V V θ=为扰动速度，0'V V V =+。

（3-28）式的分量方程可以写成
12cot 0R R V V V V R R R R
θθ
θθ∂∂+++=∂∂ （3-29） 222222221cot 2
(cot )R R R R R
V V V V V V V R R R R R R θθθθθθθ
∂∂∂∂∂+++-++∂∂∂∂∂ 001'[cos sin ()]R R u V V p u V R R R θρ
θθμμθ
∂∂∂=+-+-∂∂∂ （3-30） 22222222221cot 2sin R V V V V V V R R R R R R R θθθθθθθθθθ
∂∂∂∂∂++++-∂∂∂∂∂。