第四章 粘性流体动力学基础
4.1 运算符号
\(\nabla \cdot \overrightarrow{V}\) \(\overrightarrow{V}\) 的散度=\(\frac{\partial {{V}_{x}}}{\partial x}+\frac{\partial {{V}_{y}}}{\partial y}+\frac{\partial {{V}_{z}}}{\partial z}\)
\(\nabla \times \overrightarrow{V}\) \(\overrightarrow{V}\)的旋度=\(\overrightarrow{i}(\frac{\partial {{V}_{z}}}{\partial y}-\frac{\partial {{V}_{y}}}{\partial z})+\overrightarrow{j}(\frac{\partial {{V}_{x}}}{\partial z}-\frac{\partial {{V}_{z}}}{\partial x})+\overrightarrow{k}(\frac{\partial {{V}_{y}}}{\partial x}-\frac{\partial {{V}_{x}}}{\partial y})\)
\(\frac{D()}{Dt}=\frac{\partial ()}{Dt}+(\overrightarrow{V}\cdot \nabla )\)，()实质导数，表示跟中流体微团的运动
\(\frac{\partial }{\partial ()}\) 对于()的偏微分
\((\nabla \cdot \overrightarrow{V})\overrightarrow{V}=\nabla \left( \frac{{{V}^{2}}}{2} \right)-\overrightarrow{V}\times (\nabla \times \overrightarrow{V})\)
\(\nabla (\rho \overrightarrow{V})=(\nabla \cdot \overrightarrow{V})\rho +\rho \nabla \cdot \overrightarrow{V}\)
\({{\sigma }_{ij}}={{\tau }_{ij}}-p{{\delta }_{ij}}\)
x，y，z 笛卡尔坐标
4.2 引言
第三章讲述了可压缩流体定常多维绝热无粘性流动的一般特征，并讨论了描述多维流动的重要概念和重要定理。
粘性是流体的重要属性，在自然界和工程领域的真实流动都是具有粘性的流动。前面几章所讨论的无粘流动只是某种情况下真实流动的近似模型，而在本章将开展粘性流动的讨论。首先，针对层流流动，引入流体粘性的概念，通过对流体微元表面受力分析得到应力张量。然后，在此基础之上，建立流体的本构方程，即应力张量和速度变形张量间的关系。最后，建立粘性流体的控制方程。无论是无粘流动还是粘性流动，根据质量守恒建立的连续性方程都是一致的，粘性主要影响动量方程和能量方程。因此，本章将着重介绍根据牛顿第二定律和热力学第一定律推导粘性流动的动量方程和能量方程。
4.3 流体的输运性质
流体由非平衡状态转向平衡状态时，流体物理量的传递性质称为流体的输运性质。当两流体层之间存在速度差时，通过动量传递，速度差就会减小而使速度趋向均匀；若流体各部分之间存在温度差，则通过能量传递过程，使得温度差减小而使温度趋向均匀。流体内部通过这种自发过程达到新的平衡状态。
流体的输运性质主要指动量输运、能量输运和质量输运。从宏观上看，它们分别表现为流体的粘性、导热性和扩散性。
4.3.1 流体的粘性
在引入流体粘性的概念前，有必要回顾材料力学中关于弹性固体变形的相关概念，以便加深对流体粘性的理解。针对弹性固体，在材料力学中，我们知道弹性固体的变形阻抗是弹性模量，其中对于剪切变形，存在剪切模量为
（4.1）
由式（4.1）知，弹性固体的剪切模量可以表示为固体剪应力与剪应变之比。对于流体而言，也存在着相似的关系式，可以表示为层流中剪应力同流体某一性质的关系，即
（4.2）
如果流体中各流层的流速不相等，那么在相邻的两流体层之间的接触面上，就会形成一对等值而反向的摩擦力（或粘性阻力）来阻碍两流体层做相对运动。流体质点具有抵抗其质点做相对运动的性质，称为流体的粘性。因此，粘性是流体的一个物性参数，它与流体的温度、组分、压力有关，而与应变率无关。
图4.1表示一块平板安装在风洞的试验段中，气流沿平板板面流动，测出沿板面法线方向的气流速度分布。由图可以看出，沿平板的法向方向，在离开平板上方较远的地方，其流速与外部气流速度V0基本相同，而越靠近平板，速度越小。到板面上，流体粘附在板面上流速为零。上述流速分布说明，每一运动较慢的流体层，都是在运动较快的流体层的带动下运动的，同时运动快的流体层也受到运动较慢的流体层的阻碍，使其流速减小。这样一层一层影响下去，就有相当多层的流体，在粘性的作用下受到阻滞而减小了流速。结果就形成了图4-1所示的速度分布。平板上气流速度出现上述分布情况，正是由于气体粘性作用的结果。
图4.1 绕平板的平面流动
粘性阻力产生的物理原因是由于分子间的吸引力和分子不规则运动的动量交换引起的。一方面，当相邻流体层具有相对速度时，快层分子产生的吸引力拖动慢层，而慢层分子产生的吸引力阻滞快层，因此产生了两流体之间吸引力所形成的阻力。另一方面，由于分子不规则运动时，各流体层之间互有分子迁移和掺混，快层分子进入慢层时，给慢层交换动量，使慢层加速，反之，慢层分子迁移到快层时，给快层传递动量而使快层减速。因此产生了分子不规则运动的动量交换所形成的阻力。综上所述，流体粘性现象就是动量输运的结果。
粘性阻力即为流体的内摩擦力，可根据牛顿内摩擦定律确定。该定律由实验得出，其数学表达式为
\(F=\mu \frac{dV}{dy}A\) （4.3.1）
式中，F是做相对运动两流体层之间接触面上的内摩擦力（N）；A是接触面积；dV/dy 是沿接触面外法线方向的速度梯度； \(\mu\)动力粘度（\(Pa\cdot s\) ）。
对于单位面积上的内摩擦力
\(\tau =\frac{F}{A}=\mu \frac{dV}{dy}\) (4.3.2)
在流体力学中，粘度\(\mu\)经常与流体密度\(\rho \)结合在一起，以\({\mu }/{\rho }\;\)的形式出现。所以将该比值定义为运动粘度，并用v来表示，其单位为m/s2，即
\(v=\frac{\mu }{\rho }\) (4.3.3）
动力粘性系数 \(\mu\)随温度和压力而变化，但是压力的影响甚微。对于液体，它随温度升高而升高，对于气体则随温度升高而增大。
对于各种气体的粘性系数可以近似采用下列公式来计算
\(\frac{\mu }{{{\mu }_{0}}}\approx {{(\frac{T}{{{T}_{0}}})}^{n}}\) （4.3.4）
其中， T0=273.16K；\({{\mu }_{0}}\)为一个大气压下0℃时气体的动力粘性系数；n是温度系数。
萨瑟兰（Sutherland）给出了更准确的公式
\(\frac{\mu }{{{\mu }_{0}}}\approx {{(\frac{T}{{{T}_{0}}})}^{1.5}}(\frac{{{T}_{0}}+{{T}_{s}}}{T+{{T}_{s}}})\) （4.3.5）
式中，TS为萨瑟兰常数，与气体性质有关。
4.3.2 流体的导热性
当流体中沿着某个方向存在着温度梯度时，热量就会从温度高的地方传向温度低的地方，这种热量传递的性质称为流体的导热性。温度梯度与热流量的关系遵循傅里叶定律
\(q=-\lambda \nabla T\) （4.3.6）
式中，q是热流矢量，单位为\(W/{{m}^{2}}\)；\(\nabla T\) 是温度梯度；\(\lambda \)为热传导系数，国际制单位为\(W/m\cdot K\)。负号表示热量的传递方向与温度梯度方向相反。
气体的热传导与内摩擦有关，从微观来看，二者都是源自于分子的热运动。由于分子的热运动，两流体层之间有动量交换产生内摩擦力；两流体层之间有能量传递而产生热传导。故气体的热传导系数\(\lambda \) 和粘度\(\mu \)具有内在的联系。
4.3.3 流体的扩散性
当流体的密度分布不均匀时，流体的质量就会从高密度区迁移到低密度区，流体的这种现象称为扩散性。在单组份流体中，由于自身的密度差所引起的扩散称为自扩散。对于两种组分的混合介质，由于各组分各自密度之间的差在组分中所引起的扩散称为交互扩散。
当流体分子进行动量、能量（热能）交换且伴随有质量交换时，质量输运的机理与动量、热能的输运机理相同。对于由双组份A,B所组成的混合物系统，各组分均由其各自的高密度区向低密度区扩散。假设仅考虑组分A在组分B中的扩散，则组分A的定常扩散率与其密度梯度和截面积成正比，单位时间单位面积的质量流量与密度梯度成正比，即
\({{w}_{AB}}=-D\frac{d{{\rho }_{A}}}{dy}\) （4.3.7）
式中，\({{w}_{AB}}\) 为每单位面积质量流率；D是扩散系数，单位为\({{m}^{2}}/s\) ，它的大小取决于压强、温度和混合物的成分。方程（4.3.7）就是著名的一维定常菲克第一扩散定律。
4.4 流体表面应力
为了建立流体动力学方程，需要对流体微元体进行受力分析。微元体受到两种不同性质的力：质量力和表面力。质量力是某种力场作用在全部流体质点上的力，其大小和流体的质量成正比，比如重力、离心力、电磁力等等。表面力是由于控制面上应力作用而产生的力，比如压强p和粘性应力\(\tau \)。
1） 流体微团表面受力分析
图4.2 微元体表面力分析示意图
在上一节中，我们从一维层流出发，对粘性阻力进行了分析。现在我们考虑多维情况，在空间中取一个微元体开展分析。
首先，对微元体表面力进行受力分析。微元体垂直于x轴的面计为x平面，并以x平面上的受力分析为例进行说明，如图4.2所示。x平面上的表面力包括了粘性力\({{\tau }_{x}}\)和压力P。其中，粘性力 \({{\tau }_{x}}\)是矢量，可以将它沿三个坐标轴方向分解为三个矢量，记为\({{\tau }_{xx}},{{\tau }_{xy}},{{\tau }_{xz}},\)其中角标的第一个字母表示粘性力所在的平面，第二个字母表示分量的方向。除了粘性力\({{\tau }_{x}}\)以外，x平面上还受到压力的作用。我们知道，压力方向垂直于x轴平面，且指向流体，因此x平面上的压力为沿x轴负向。下面分别确定粘性力\({{\tau }_{x}}\)和压力P。
2）表面应力确定
图4.3 微元体的粘性应力分析示意图
我们从流体中取一个正六面微元体作为模型来分析粘性应力，如图4.3所示。在垂直于x轴的两个外表面上，分别有合力\({{\tau }_{x}}\)、\({{\tau }_{x}}+\frac{\partial {{\tau }_{x}}}{\partial x}dx\)（下标x表示应力向量作用在与x轴垂直的微元面上）。由此可得到作用在垂直于x轴微元面上的表面力的合力为\(\frac{\partial {{\tau }_{x}}}{\partial x}dxdydz\)。同理可得作用在垂直于y轴和z轴的微元面上的表面力的合力为\(\frac{\partial {{\tau }_{y}}}{\partial x}dxdydz\)、\(\frac{\partial {{\tau }_{z}}}{\partial z}dxdydz\)。于是，单位容积的表面力为
\(\tau =(\frac{\partial {{\tau }_{x}}}{\partial x}dxdydz+\frac{\partial {{\tau }_{y}}}{\partial y}dxdydz+\frac{\partial {{\tau }_{z}}}{\partial z}dxdydz)/(dxdydz)=\frac{\partial {{\tau }_{x}}}{\partial x}+\frac{\partial {{\tau }_{y}}}{\partial y}+\frac{\partial {{\tau }_{z}}}{\partial z}\) (4.4.1)
式中\({{\tau }_{x}},{{\tau }_{y}},{{\tau }_{z}}\)是向量，还可以把它们沿三个坐标方向分解，即分解成法向应力和切应力。图4.3中作用于与x轴垂直的微元面的应力\({{\tau }_{x}}\)可分解为
\({{\tau }_{x}}=\overrightarrow{{{\tau }_{x}}}={{\tau }_{xx}}\overrightarrow{i}+{{\tau }_{xy}}\overrightarrow{j}+{{\tau }_{xz}}\overrightarrow{k}\) （4.4.2a）
同理有
\({{\tau }_{y}}=\overrightarrow{{{\tau }_{y}}}={{\tau }_{yx}}\overrightarrow{i}+{{\tau }_{yy}}\overrightarrow{j}+{{\tau }_{yz}}\overrightarrow{k}\) （4.4.2b）
\({{\tau }_{z}}=\overrightarrow{{{\tau }_{z}}}={{\tau }_{zx}}\overrightarrow{i}+{{\tau }_{zy}}\overrightarrow{j}+{{\tau }_{zz}}\overrightarrow{k}\) （4.4.2c）
对（4.4.2）整理可得
用矩阵表示
上述表达式左边为矢量，右边第一项为系数矩阵，右边第二项为单位向量，各项表达式如下：
由上述关系可以看出，表面力\(\tau \)含有三个方向上的分量，即\({{\tau }_{x}},{{\tau }_{y}},{{\tau }_{z}}\)，每个分量都是矢量。因此，确定给定表面上表面力的本质是求这个面上三个力的合力，而求合力的关键就是确定上述关系式中的系数矩阵。为了表述方便，我们用\({{\tau }_{ij}}\)来代替系数矩阵，\({{\tau }_{ij}}\)又称为应力张量，即
粘性应力张量是二阶对称张量。
3）表面压力确定
表面力除了粘性应力\(\tau \)以外，还受到压力作用。与应力张量类似，压力也有在垂直于三个坐标轴的面上的分力，在每个面上的分力可以分解为三个分量，因此，也能写成压力的应力张量
可以证明，该应力张量可写成下列对角线形式
这里不给出具体证明过程，该形式的物理意义为在与主轴方向垂直的平面上，只有法向应力，切向应力等于零。
对于理想流体，同一点的各个不同方向上，法向应力应是相等的。因此，在理想流体中，只要用一个标量函数即压力p便完全地刻画了任一点上的应力状态，此时应力张量变为
即
其中，取-p的原因是强调压力与作用面的法向方向恰好相反。
4）表面力确定
我们用\(\sigma \)表示单位微元体的应力，综合上面的讨论结果，单位容积微元体的应力为\({{\sigma }_{ij}}={{\tau }_{ij}}-p{{\delta }_{ij}}\)，即
4.5 本构方程
在上一节中我们讨论了流体表面应力，那么如何获得流体表面应力就成为我们最关心的问题。我们知道，物体的应力与运动学参数之间存在着一定的关系。对于大多数流体，如空气和水，应力与应变变化率成正比，服从这种关系的流体称为牛顿流体。我们把应力张量 和变形速率张量\({{\varepsilon }_{ij}}\)联系起来的方程称为本构方程。粘性流体作直线层状运动时，切应力与层间速度梯度成正比，即
\({{\tau }_{ij}}=\mu \frac{d{{V}_{x}}}{dv}\) （4.5.1）
在这种层状运动中，\(\frac{\partial {{V}_{y}}}{\partial x}=0\)，所以此式右端的速度梯度实为应变变化率张量分量 的二倍，所以上式又可以写为
\({{\tau }_{yx}}=\mu \frac{d{{V}_{x}}}{dy}=2\mu {{\varepsilon }_{yx}}\)
（4.5.2）
斯托克斯提出了在牛顿流体中应力张量和变形速率张量之间一般关系的三项假定：
（1）流体是连续的，它的应力张量是应变变化率张量的线性函数；
（2）流体是各向同性的，也就是说流体的物理性质与方向无关，只是坐标位置的函数。
（3）所建立的方程不仅应适合运动的情况，也适合于静止的情况。
在静止状态时，粘性切应力\({{\tau }_{ij}}\)应等于零；而由流体静力学可知，这时流体内任一点受到的压力与方向无关，此时流体的应力张量为
\({{\sigma }_{ij}}=-p{{\delta }_{ij}}\) （4.5.3）
根据第一条假设，并考虑式（4.5.2），最简单的想法便是假设存在如下应力-应变关系
\({{\sigma }_{ij}}=2\mu {{\varepsilon }_{ij}}\)
这一关系式对于切应力的表示是符合牛顿粘性应力公式的，但违背了第三条假设的要求。当静止时，主对角线上的三个应变变化率\({{\varepsilon }_{x}},{{\varepsilon }_{y}},{{\varepsilon }_{z}}\)均为零，因而按上式推知，\({{\sigma }_{xx}},{{\sigma }_{yy}},{{\sigma }_{zz}}\)也都应为零，这显然与式（4.5.3）矛盾。所以我们假设应力-应变的关系为
\({{\sigma }_{ij}}=2\mu {{\varepsilon }_{ij}}+b{{\delta }_{ij}}\) （4.5.4）
根据第二条假设，所寻求的关系与坐标选择无关，由张量理论可知，二阶张量中主对角线上三个分量之和是不变量，故，假设
\(b={{b}_{1}}({{\sigma }_{xx}}+{{\sigma }_{yy}}+{{\sigma }_{zz}})+{{b}_{2}}({{\varepsilon }_{x}}+{{\varepsilon }_{y}}+{{\varepsilon }_{z}})+{{b}_{3}}\)
（4.5.5）
将此式代入式（4.5.4），并取等式两边主对角线上三个分量之和，然后归并同类项后可得
\((1-3{{b}_{1}})({{\sigma }_{xx}}+{{\sigma }_{yy}}+{{\sigma }_{zz}})=(2\mu +3{{b}_{2}})({{\varepsilon }_{x}}+{{\varepsilon }_{y}}+{{\varepsilon }_{z}})+{{b}_{3}}\) （4.5.6）
静止时
\({{\sigma }_{xx}}={{\sigma }_{yy}}={{\sigma }_{zz}}=-p,{{\varepsilon }_{x}}={{\varepsilon }_{y}}={{\varepsilon }_{z}}=0\)
于是上式变为
常数b1，b3有两种选择方法可使此式成立且不再包含任何待定系数
第一种：\({{b}_{1}}=0,{{b}_{3}}=-p\) （4.5.8）
第二种：\({{b}_{1}}=\frac{1}{3},{{b}_{3}}=0\) （4.5.9）
在后面将会看出，第一种选择可包含更广的情况，在进一步的假设条件下可得出第二种选择的结果。若定义
\(\overline{p}=-\frac{{{\sigma }_{xx}}+{{\sigma }_{yy}}+{{\sigma }_{zz}}}{3}\) （4.5.10）
表示运动状态下某点的平均压强，将第一种假设代入式（4.5.6）中，并用 代替 ，则得到
\(\overline{p}=p-(\lambda +\frac{2}{3}\mu )div\upsilon \) （4.5.11）
如果将p理解为热力学参数的压强，则由此式可看出，一般情况下，运动时的平均压强\(\overline{p}\)并不等于热力学压强p。
系数\(\lambda \)和\(\mu \)由相同的量纲，且与体积膨胀率\(div\upsilon \)有关，所以称为体积粘性系数或者第二粘性系数。
为了消除\(\overline{p}\ne p\)的情况，斯托克斯假设
\(\lambda =-\frac{2}{3}\mu \) （4.5.12）
将式（4.5.6）和式（4.5.8）、（4.5.12）代入式（4.5.4），得到
\({{\sigma }_{ij}}=2\mu {{\varepsilon }_{ij}}-\frac{2}{3}\mu div\upsilon {{\delta }_{ij}}-p{{\delta }_{ij}}\) （4.5.13）
公式（4.5.13）可展开为
已知
\({{\sigma }_{xx}}={{\tau }_{xx}}-p\)
\({{\sigma }_{yy}}={{\tau }_{yy}}-p\) （4.5.15）
\({{\sigma }_{zz}}={{\tau }_{zz}}-p\)
所以，结合公式（4.5.15）和式（4.5.14）得到
\({{\tau }_{xx}}=2\mu \frac{\partial {{V}_{x}}}{\partial x}-\frac{2}{3}\mu (\nabla \cdot V)\)
\({{\tau }_{yy}}=2\mu \frac{\partial {{V}_{y}}}{\partial y}-\frac{2}{3}\mu (\nabla \cdot V)\)
\({{\tau }_{zz}}=2\mu \frac{\partial {{V}_{z}}}{\partial z}-\frac{2}{3}\mu (\nabla \cdot V)\) （4.5.16）
\({{\tau }_{xy}}={{\tau }_{yx}}=\mu (\frac{\partial {{V}_{y}}}{\partial x}+\frac{\partial {{V}_{x}}}{\partial y})\)
\({{\tau }_{yz}}={{\tau }_{yz}}=\mu (\frac{\partial {{V}_{z}}}{\partial y}+\frac{\partial {{V}_{y}}}{\partial z})\)
\({{\tau }_{zx}}={{\tau }_{xz}}=\mu (\frac{\partial {{V}_{x}}}{\partial z}+\frac{\partial {{V}_{z}}}{\partial x})\)
4.6 粘性流动动量方程的推导
动量方程是动量守恒定律对于黏性流体运动规律的数学描述，可以由牛顿第二定律推导出。
牛顿第二定律：
\(\overrightarrow{F}=m\cdot \overrightarrow{a}\)
将牛顿第二定律应用在图4.4所示的运动流体微团，则有作用于微团上力的总和等于微团的质量与微团运动时加速度的乘积。这是一个矢量关系式，可以沿\(x,y,z\(轴分解成三个标量的关系式。下面，以x方向动量方程为例，建立动量方程。
图4.4 运动的无穷小微团模型
(图中只画出了x方向的力，用于推导x方向动量方程)
先仅考虑x方向的动量分量，牛顿第二定律具有如下形式：
\({{F}_{x}}=m{{a}_{x}}\) （4.6.1）
其中，\({{F}_{x}}\) 和\({{a}_{x}}\) 分别是微团所受力（分力）和加速度在x方向的分量。
1) 力分析
首先，考虑方程(4.6.1)的左边，是运动的流体微团受到的x方向的力，这个力有以下两个来源：
（1） 体积力
直接作用于流体微团整个体积微元上的力，而且作用是超距离的，比如重力、电场力、磁场力。
（2） 表面力
直接作用于流体微团的表面上的力，只考虑以下两种形式的力：
①由包围在流体微团周围的流体所施加的，作用于微团表面的压力分布；
②由于外部流体推拉微团而产生的，以摩擦的方式作用于表面的切应力和正应力分布。
体积力表达式：
将作用在单位质量流体微团上的体积力记做f，其x方向分量记做\({{f}_{x}}\) 。流体微团的体积为dxdydz，因此
表面力表达式：
流体微团的切应力和正应力与流体微团变形的时间变化率相关联，图4-4给出了xy平面内的情形。
切应力
在图4-4a中用\({{\tau }_{xy}}\)表示，与流体微团剪切变形的时间变化率有关；
正应力
在图4-4b中用\({{\tau }_{xx}}\)表示，与流体微团体积的时间变化率有关。
不论是切应力还是正应力，都依赖于流动的速度梯度，后面将对它们进行分析。在大多数粘性流动中，正应力(例如\({{\tau }_{xx}}\))要比切应力小得多，很多情形下可以忽略。然而，当法向速度梯度很大时(例如，在激波内部)，正应力(x方向就是\({{\tau }_{xx}}\)）就变得重要。
图4.5 正应力与切应力的示意图
施加在流体微团x方向的全部表面力如图4.4所示，我们约定用\({{\tau }_{ij}}\)表示j方向的应力作用在垂直于i轴的平面上。
2) 微元体受力分析（参见图4-3）
（1）面abcd上受到的力：
仅存在由切应力引起的x方向分力\({{\tau }_{yx}}dxdz\)，其中\({{\tau }_{yx}}\)是向左的(与x轴方向相反)。
（2）面efgh上受到的力
面dfgh与面abcd的距离为dy，所以efgh面上x方向的切应力为\(\left( {{\tau }_{yx}}+(\partial {{\tau }_{yx}}/\partial y)dy \right)dxdz\)。其中力\({{\tau }_{yx}}+(\partial {{\tau }_{yx}}/\partial y)dy\)是向右的(与x轴方向相同)。
（3）面abfe上的力
仅存在沿x轴正向的切应力\(({{\tau }_{zx}}+(\partial {{\tau }_{zx}}/\partial z)dz)dxdy\)
（4）面dcgh上的力
仅存在沿x轴负方向的切应力\({{\tau }_{zx}}dxdy\)。
（5）面adhe上的力
存在沿x轴正向的压力pdxdz，同时存在沿x轴负向的应力\({{\tau }_{xx}}dydz\)。
（6）面bcgf上的力
存在压力\((p+(\partial p/\partial x)dx)dydz\)，方向沿x轴负方向。存在正应力\(({{\tau }_{xx}}+(\partial {{\tau }_{xx}}/\partial x)dx)dydz\)，沿x轴正向。
（7） 合力
综上所述，对运动的流体微团，有
3) x方向总的表面力
（4.6.3）
x方向总的力 ，可以由式（4.6.2）和式（4.6.3）相加得到，消去相同的项，得
\({{F}_{x}}=(-\frac{\partial p}{\partial x}+\frac{\partial {{\tau }_{xx}}}{\partial x}+\frac{\partial {{\tau }_{yx}}}{\partial y}+\frac{\partial {{\tau }_{zx}}}{\partial z})dxdydz+\rho {{f}_{x}}dxdydz\) （4.6.4）
式(4.6.4)给出了式(4.6.1)的左边。
下面考虑方程(4.6.1)的右边。运动的流体微团，其质量是固定不变的，等于
\(m=\rho dxdydz\) （4.6.5）
流体微团的加速度是速度的时间变化率，则加速度的x方向分量 ，应等于 的时间变化率。但由于我们考虑运动的流体微团，因此这个时间变化率是由物质导数给出的，即
\({{a}_{x}}=\frac{Du}{Dt}\) （4.6.6）
将式(4.6.1)与式(4.6.4)、式(4.6.6)综合起来，我们得到
\(\rho \frac{Du}{Dt}=-\frac{\partial p}{\partial x}+\frac{\partial {{\tau }_{xx}}}{\partial x}+\frac{\partial {{\tau }_{yx}}}{\partial y}+\frac{\partial {{\tau }_{zx}}}{\partial z}+\rho {{f}_{x}}\) （4.6.7a）
这就是粘性流x方向的动量方程。
用同样的办法，可得到y方向和z方向的动量方程
\(\rho \frac{Dv}{Dt}=-\frac{\partial p}{\partial y}+\frac{\partial {{\tau }_{xy}}}{\partial x}+\frac{\partial {{\tau }_{yy}}}{\partial y}+\frac{\partial {{\tau }_{zy}}}{\partial z}+\rho {{f}_{y}}\) (4.6.7b)
\(\rho \frac{Dw}{Dt}=-\frac{\partial p}{\partial z}+\frac{\partial {{\tau }_{xz}}}{\partial x}+\frac{\partial {{\tau }_{yz}}}{\partial y}+\frac{\partial {{\tau }_{zz}}}{\partial z}+\rho {{f}_{z}}\) (4.6.7c)
方程(4.6.7a、b、c)分别是x、y、z方向的动量方程。注意到，它们都是偏微分方程，是通过将基本的物理学原理应用于无穷小流体微团直接得到的。同时，由于流体微团是运动的，所以方程(4.6.7a、b、c)是非守恒形式的。它们都是标量方程，统称为纳维-斯托克斯方程。
按照下面的方法，可以得到纳维-斯托克斯方程的守恒形式。
根据物质导数的定义，可将方程(4.6.7a)的左边写成
\(\rho \frac{Du}{Dt}=\rho \frac{\partial u}{\partial t}+\rho V\cdot \nabla u\) （4.6.8）
另外，展开下面的导数
\(\frac{\partial (\rho u)}{\partial t}=\rho \frac{\partial u}{\partial t}+u\frac{\partial \rho }{\partial t}\)
整理，得
\(\rho \frac{\partial u}{\partial t}=\frac{\partial (\rho u)}{\partial t}-u\frac{\partial \rho }{\partial t}\) （4.6.9）
利用标量与向量乘积的散度的向量恒等式，有
\(\nabla \cdot (\rho uV)=u\nabla \cdot (\rho V)+(\rho V)\cdot \nabla u\)
或改写成
\((\rho V)\cdot \nabla u=\nabla \cdot (\rho uV)-u\nabla \cdot (\rho V)\) （4.6.10）
将式(4.6.9)和式(4.6.10)代人式(4.6.8)，得
\(\rho \frac{Du}{Dt}=\frac{\partial \rho u}{\partial t}-u\frac{\partial \rho }{\partial t}-u\nabla \cdot (\rho V)+\nabla \cdot (\rho uV)=\frac{\partial (\rho u)}{\partial t}-u(\frac{\partial \rho }{\partial t}-\nabla \cdot (\rho V))+\nabla \cdot (\rho uV)\) （4.6.11)
此式右边方括号里的表达式就是连续性方程(2.7.7)的左边，所以方括号中的项等于零，于是(4.6.11)可以简化为
\(\rho \frac{Du}{Dt}=\frac{\partial (\rho u)}{\partial t}+\nabla \cdot (\rho uV)\) (4.6.12)
再将式(4.6.12)代人式(4.6.7a)，得
\(\frac{\partial (\rho u)}{\partial t}+\nabla \cdot (\rho uV)=-\frac{\partial p}{\partial x}+\frac{\partial {{\tau }_{xx}}}{\partial x}+\frac{\partial {{\tau }_{yx}}}{\partial y}+\frac{\partial {{\tau }_{zx}}}{\partial z}+\rho {{f}_{x}}\) （4.6.13a）
同样，方程(4.6.7b、c)可以写成
\(\frac{\partial (\rho v)}{\partial t}+\nabla \cdot (\rho vV)=-\frac{\partial p}{\partial y}+\frac{\partial {{\tau }_{xy}}}{\partial x}+\frac{\partial {{\tau }_{yy}}}{\partial y}+\frac{\partial {{\tau }_{zx}}}{\partial z}+\rho {{f}_{y}}\) （4.6.13b）
\(\frac{\partial (\rho w)}{\partial t}+\nabla \cdot (\rho wV)=-\frac{\partial p}{\partial z}+\frac{\partial {{\tau }_{xz}}}{\partial x}+\frac{\partial {{\tau }_{yz}}}{\partial y}+\frac{\partial {{\tau }_{zz}}}{\partial z}+\rho {{f}_{z}}\) （4.6.13c）
方程(4.6.13a)到方程(4.6.13c)就是纳维-斯托克斯方程的守恒形式。
流体的切应力与应变的时间变化率，即速度梯度，是成正比的。这样的流体被称为牛顿流体(切应力 与速度梯度不成正比的流体称为非牛顿流体，例如，血液)。在空气动力学的所有实际问题中，流体都可以被看成是牛顿流体。对于这样的流体则有，
\({{\tau }_{xx}}=\lambda (\nabla \cdot V)+2\mu \frac{\partial u}{\partial x}\) (4.6.14a)
\({{\tau }_{yy}}=\lambda (\nabla \cdot V)+2\mu \frac{\partial v}{\partial y}\) (4.6.14b)
\({{\tau }_{zz}}=\lambda (\nabla \cdot V)+2\mu \frac{\partial w}{\partial z}\) (4.6.14c)
\({{\tau }_{xy}}={{\tau }_{yx}}=\mu (\frac{\partial v}{\partial x}+\frac{\partial u}{\partial y})\) (4.6.14d)
\({{\tau }_{xz}}={{\tau }_{zx}}=\mu (\frac{\partial u}{\partial z}+\frac{\partial w}{\partial x})\) (4.6.14e)
\({{\tau }_{yzy}}={{\tau }_{zy}}=\mu (\frac{\partial w}{\partial y}+\frac{\partial v}{\partial z})\) (4.6.14f)
其中\(\mu \)是分子粘性系数，\(\lambda \)是第二粘性系数。斯托克斯提出假设，认为
\(\lambda =-\frac{2}{3}\mu \)
将式(4.6.15)各分式代人方程(4.6.13)各式，得到完整的纳维-斯托克斯方程守恒形式
4.7 粘性流动能量方程的推导
在本节，应用到了第三个物理学原理：能量守恒。为了与4.6节纳维-斯托克斯方程(动量方程)的推导保持一致，采用随流体运动的无穷小微团的流动模型。上述物理学原理其实就是热力学第一定律。对于和流体一起运动流体微团模型而言，这个定律表述如下：
下文用A、B、C代表与文字相对应的各项。
1) 力所做功分析
首先来计算等号左边第二项C，即得到体积力和表面力对运动着的流体微团做功的功率表达式。可以证明，作用在一个运动物体上的力对物体做功的功率等于这个力的大小乘以速度在此力作用方向上的分量。所以，作用于速度为V的流体微团上的体积力，做功的功率为
\(\rho f\cdot V(dxdydz)\)
至于表面力(压力加上切应力和正应力)，只考虑作用x方向上的力，如图4-5所示。在图4-5中，x方向上压力和切应力对流体微团做功的功率，就等于速度的分量u乘以力(比如，在面abcd上为\({{\tau }_{xy}}dxdy\)，即\(u{{\tau }_{xy}}dxdy\)。在其他面上也有类似的表达式。
对比图4-5中作用在面adhe和bcgf上的压力，则压力在x方向上做功的功率为
\(\left( up-\left( up+\frac{\partial (up)}{\partial x}dx \right) \right)dydz=-\frac{\partial (up)}{\partial x}dxdydz\)
类似地，在面abcd和面efgh上，切应力在x方向上做功的功率是
\(\left( \left( u{{\tau }_{yx}}+\frac{\partial (u{{\tau }_{yx}})}{\partial y}dy \right) \right)dxdz=-\frac{\partial (u{{\tau }_{yx}})}{\partial y}dxdydz\)
图4-5中所有表面力对运动流体微团做功的功率为
\(\left( -\frac{\partial (up)}{\partial x}+\frac{\partial (u{{\tau }_{xx}})}{\partial x}+\frac{\partial (u{{\tau }_{yx}})}{\partial y}+\frac{\partial (u{{\tau }_{zx}})}{\partial z} \right)dxdydz\)
上式仅考虑了x方向上的表面力。再考虑y和z方向上的表面力，也能得到类似的表达式。加在一起，对运动流体微团做功的功率是x、y和z方向上表面力贡献的总和，在式(4.7.1)中记作C，即
\(C=\left( \left( \frac{\partial (up)}{\partial x}+\frac{\partial (up)}{\partial x}+\frac{\partial (up)}{\partial x} \right)+\frac{\partial (u{{\tau }_{xx}})}{\partial x}+\frac{\partial (u{{\tau }_{yx}})}{\partial y}+\frac{\partial (u{{\tau }_{zx}})}{\partial z}+\frac{\partial (v{{\tau }_{xy}})}{\partial x}+\frac{\partial (v{{\tau }_{yy}})}{\partial y}+\frac{\partial (v{{\tau }_{zy}})}{\partial z}\text{+}\frac{\partial (w{{\tau }_{xz}})}{\partial x}+\frac{\partial (w{{\tau }_{yz}})}{\partial y}+\frac{\partial (w{{\tau }_{zz}})}{\partial z} \right)dxdydz+\rho f\cdot Vdxdydz\) （4.7.2）
注意，式(4.7.2)右边的前三项就是\(\nabla \cdot \left( pV \right)\)。
图4.6 运动无穷小流体微团的能量通量
(用于推导能量方程。为简单起见，图中只画出了x方向的通量)
2) 运动无穷小流体微团的能量通量
再考虑式(4.7.1)中的B项，即进人微团内的总热流量。这一热流来自于体积加热，如吸收或释放的辐射热；由温度梯度导致的跨过表面的热输运，即热传导。定义\(\dot{q}\)为单位质量的体积加热率。在图4-5中，运动流体微团的质量为\(\rho dxdydz\)，我们由此得到
在图4.6中，热传导从面adhe输运给微团内的热量是\({{\dot{q}}_{x}}dydz\)，其中\({{\dot{q}}_{x}}\)是热传导在单位时间内通过单位面积在x方向上输运的热量。给定方向上的热传导，若以单位时间内通过垂直于该方向的单位面积的能量来表述，称作该方向上的热流。这里\({{\dot{q}}_{x}}\)就是x方向上的热流。经过面bcgf输运到微团外的热量是
\(\left( {{{\dot{q}}}_{x}}-\left( {{{\dot{q}}}_{x}}+\frac{\partial {{{\dot{q}}}_{x}}}{\partial x}dx \right) \right)dydz=-\frac{\partial {{{\dot{q}}}_{x}}}{\partial x}dxdydz\)
再加上图4-5中通过其他面在y和z方向上的热的输运量，我们可以得到
式(4.7.1)中的B项是式(4.7.3)和式(4.7.4)之和，即
\(B=\left( p\dot{q}-\left( \frac{\partial {{{\dot{q}}}_{x}}}{\partial x}+\frac{\partial {{{\dot{q}}}_{y}}}{\partial y}+\frac{\partial {{{\dot{q}}}_{z}}}{\partial z} \right) \right)dxdydz\) （4.7.5）
根据傅里叶热传导定律，热传导产生的热流与当地的温度梯度成正比
\({{\dot{q}}_{x}}=-k\frac{\partial T}{\partial x},{{\dot{q}}_{y}}=-k\frac{\partial T}{\partial y},{{\dot{q}}_{z}}=-k\frac{\partial T}{\partial z}\)
其中k为热导率。所以，式(4.7.5)可写成
\(B=\left( p\dot{q}+\frac{\partial }{\partial x}\left( k\frac{\partial T}{\partial x} \right)+\frac{\partial }{\partial y}\left( k\frac{\partial T}{\partial y} \right)+\frac{\partial }{\partial z}\left( k\frac{\partial T}{\partial z} \right) \right)dxdydz\) （4.7.6）
对于式(4.7.1)中的A项，即运动流体微团的能量，有两个来源：
(1)由于分子随机运动而产生的(单位质量)内能e。
(2)流体微团平动时具有的动能，单位质量的动能为\({{V}^{2}}/2\) 。
因此，运动着的流体微团既有动能又有内能，两者之和就是总能量。在式(4.7.1)中，A项表示的能量是总能量，即内能与动能之和。这一总能量为\(e+{{V}^{2}}/2\)。由于跟随着一个运动的流体微团，单位质量的总能量变化的时间变化率由物质导数给出。流体微团的质量为\(\rho dxdydz\)，所以有
\(A=\rho \frac{D}{Dt}(e+{{V}^{2}}/2)dxdydz\) （4.7.7）
3) 能量方程
将式(4.7.2)、式(4.7.6)和式(4.7.7)代人式(4.7.1)，得到能量方程的最终形式为
这是非守恒形式的能量方程，并且是用总能量\(e+{{V}^{2}}/2\)表示的。再次强调，对一个运动的流体微团运用基本的物理学原理，得到的方程是非守恒形式。
方程(4.7.8)左侧包含了总能量的物质导数，\(D(e+{{V}^{2}}/2)/Dt\)，这只是能量方程许多不同形式中的一种，它是对运动流体微团直接运用能量守恒原理所得到的形式。这个方程很容易从以下两个方面进行改动。
(1)方程左边可以只用内能e、只用焓h或者只用总焓\({{h}_{0}}=h+{{V}^{2}}/2\)来表示，方程的右边也随之变动。例如，我们在下一段会将方程(4.7.8)转化为关于De/Dt的方程，并给出所需的演算。
(2)能量方程，对上述每一种不同形式，都有守恒形式和非守恒形式。这两种形式之间转换的演算也将在下面讨论。
从方程(4.7.8)出发，先将它改写成只用 的形式，将方程(4.6.7a~c)分别乘u、v、w得
\(\rho \frac{D}{Dt}(\frac{{{u}^{2}}}{2})=-u\frac{\partial p}{\partial x}+u\frac{\partial {{\tau }_{xx}}}{\partial x}+u\frac{\partial {{\tau }_{yx}}}{\partial y}+u\frac{\partial {{\tau }_{zx}}}{\partial z}+\rho u{{f}_{x}}\) （4.7.9）
\(\rho \frac{D}{Dt}(\frac{{{v}^{2}}}{2})=-v\frac{\partial p}{\partial x}+v\frac{\partial {{\tau }_{xx}}}{\partial x}+v\frac{\partial {{\tau }_{yx}}}{\partial y}+v\frac{\partial {{\tau }_{zx}}}{\partial z}+\rho v{{f}_{y}}\) （4.7.10）
\(\rho \frac{D}{Dt}(\frac{{{w}^{2}}}{2})=-w\frac{\partial p}{\partial x}+w\frac{\partial {{\tau }_{xx}}}{\partial x}+w\frac{\partial {{\tau }_{yx}}}{\partial y}+w\frac{\partial {{\tau }_{zx}}}{\partial z}+\rho w{{f}_{z}}\) （4.7.11）
将式(4.7.9)~式(4.7.11)各式加在一起，并注意\({{u}^{2}}+{{v}^{2}}+{{w}^{2}}={{V}^{2}}\)，可得
\(\rho \frac{D}{Dt}(\frac{{{V}^{2}}}{2})=-u\frac{\partial p}{\partial x}-v\frac{\partial p}{\partial x}-w\frac{\partial p}{\partial x}+u\left( \frac{\partial {{\tau }_{xx}}}{\partial x}+\frac{\partial {{\tau }_{yx}}}{\partial y}+\frac{\partial {{\tau }_{zx}}}{\partial z} \right)+v\left( \frac{\partial {{\tau }_{xy}}}{\partial x}+\frac{\partial {{\tau }_{yy}}}{\partial y}+\frac{\partial {{\tau }_{zy}}}{\partial z} \right)+w\left( \frac{\partial {{\tau }_{xz}}}{\partial x}+\frac{\partial {{\tau }_{yz}}}{\partial y}+\frac{\partial {{\tau }_{zz}}}{\partial z} \right)+\rho (u{{f}_{x}}+u{{f}_{y}}+u{{f}_{z}})\)
(4.7.12)
从方程(4.7.8)中减去式(4.7.12)，注意可\(\rho f\cdot V=\rho (u{{f}_{x}}+u{{f}_{v}}+u{{f}_{z}})\)，我们有
方程(4.7.13)这种形式的能量方程，其左边只包含了内能的物质导数，动能的物质导数和右边的体积力己经去掉。要注意到方程(4.7.8)中正应力和切应力是与速度相乘，一起出现在x、y、z的导数内。与之相比，方程(4.7.13)中粘性应力单独出现，直接与速度梯度相乘。
在式(4.6.14d~f)中有\({{\tau }_{xy}}={{\tau }_{yx}},{{\tau }_{xz}}={{\tau }_{zx}},{{\tau }_{yz}}={{\tau }_{zy}}\)当流体微团的体积缩成一点的时候，切应力的这种对称性可以避免流体微团的角速度，它与作用在流体微团上的力矩有关，趋于无穷大)，因此可以合并方程(4.7.13)中的一些项，得到
为了用速度梯度表示粘性应力，再次利用式(4.6.14a~f)各式，方程(4.7.14)又可以写成
方程(4.7.15)是完全用流场变量表示的能量方程。利用方程(4.6.14a~f)，对方程(4.7.8)也可以进行类似的变换，得到关于流场变量的能量方程。
再次强调在方程(4.7.15)左边只出现了内能。能量方程的左边可以用不同的能量形式表示。例如，方程(4.7.8)用总能量，方程(4.7.15)用内能。之前曾经讲过，用焓h和总焓h+V2/2表示的形式也可以通过类似的变换得到。
能量方程的左边可以用能量的不同形式表示，而能量方程的右边也有相应的不同形式，这只是能量方程的一个方面。现在我们描述能量方程的另一方面，也就是与连续性方程和动量方程相同的方面：能量方程也可以表达为守恒形式。方程(4.7.8)、式(4.7.13)、式(4.7.15)所给出的能量方程，左边都出现物质导数，因而都是非守恒形式。它们直接由运动流体微团模型得出。考虑方程(4.7.15)的左边，由物质导数的定义
\(\rho \frac{De}{Dt}=\rho \frac{\partial e}{\partial t}+\rho V\cdot \nabla e\) (4.7.16)
但
\(\rho \frac{\partial (\rho e)}{\partial t}=\rho \frac{\partial e}{\partial t}+e\frac{\partial \rho }{\partial t}\)
或
\(\rho \frac{\partial e}{\partial t}=\rho \frac{\partial (\rho e)}{\partial t}-e\frac{\partial \rho }{\partial t}\) (4.7.17)
另一方面，对于标量与向量乘积的散度，有向量恒等式
\(\nabla \cdot (\rho eV)=e\nabla \cdot (\rho V)+\rho V\cdot \nabla e\)
或写成
\(\rho V\cdot \nabla e=\nabla \cdot (\rho eV)-e\nabla \cdot (\rho V)\) (4.7.18)
将式(4.7.17)和式(4.7.18)代人式(4.7.16)，得
\(\rho \frac{De}{Dt}=\rho \frac{\partial (\rho e)}{\partial t}-e\left( \frac{\partial \rho }{\partial t}+\nabla \cdot \left( \rho V \right) \right)\nabla \cdot \left( \rho eV \right)\) (4.7.19)
由连续性方程(4.7.15)可知，式(4.7.19)右边方括号内的式子等于零，于是式(4.7.19)就变成
\(\rho \frac{De}{Dt}=\rho \frac{\partial (\rho e)}{\partial t}+\nabla \cdot (\rho Ve)\) (4.7.20)
将式(4.7.20)代入方程(4.7.15)，我们有
这是用内能表示的守恒形式的能量方程。
将内能e改为总能量e+V2/2，重复由式(4.7.16)到式(4.7.20)的推导过程，可以得到
\(\rho \frac{D}{Dt}\left( e+\frac{{{V}^{2}}}{2} \right)=\frac{\partial }{\partial t}\left( \rho \left( e+\frac{{{V}^{2}}}{2} \right) \right)+\nabla \cdot \left( \rho \left( e+\frac{{{V}^{2}}}{2} \right)V \right)\) (4.7.22)
将式(4.7.22)代人方程(4.7.8)的左边，我们得到
方程(4.7.23)是用总能量e+V2/2表示的守恒形式的能量方程。
要将方程的非守恒形式转化为守恒形式，只需要改变方程的左边就可以了，方程的右边保持不变。例如，对比方程(4.7.15)与方程( 4.7.20)，两者都是用内能表示的，方程(4.7.15)是非守恒形式的，而方程(4.7.21)是守恒形式的。它们只是左边不同，右边则是相同的。比较一下方程(4.7.8)和方程(4.7.23)，也是一样。
4.8 小结
在这一章中，通过建立流体本构方程，根据动量守恒定律和能量守恒定律我们推导了粘性流动的动量方程和能量方程。为参考方便，这些方程分别总结在表4.1中。
表4.1 粘性流动动量方程和能量方程
1. 试画出下述条件下流体微元的变形图：
（a）\(\partial {{v}_{x}}/\partial y\gg \partial {{v}_{v}}/\partial x.\);
（b）\(\partial {{v}_{v}}/\partial x\gg \partial {{v}_{x}}/\partial y\)。
2. 考虑两个平行板之间的粘性不可压缩流体的运动。设两板为无限平面，间距为h，上板不动，下板以常速U沿板向运动。设板向压力梯度为常数，运动定常，流体所受外力不计。
（1）研究流体的运动规律，即求速度分布、流量、平均速度、最大速度、内摩擦力分布及作用在运动板上的摩擦力。
（2）若沿板向没有压力梯度，流体的速度分布如何？
（3）若沿板向的压力梯度为常数，但两板均不动，流体的速度分布又怎样？
3. 已知单位体积的剪切功率为\(\tau v\)，试求一个圆管内抛物线的速度分布，并求出距管壁多远处剪切功最大？
4. 对于一个速度\({{v}_{x}}={{v}_{x}}(y)\(的二维、不可压缩流动，试画出一个三维流体微元，并标明每个应力分量的大小、方向及作用面。
5. 已知在一个\({{v}_{x}}={{v}_{x}}(x)\)的一维流动中，轴的应变率为\(\partial {{v}_{x}}/\partial x\)。试问其体积变化率是多少？若推广到三维微元体，其体积变化率又是多少？
6. 应用圆柱微元体证明，斯托克斯的粘度关系可以导出下述的剪应力分量
\({{\tau }_{r\theta }}={{\tau }_{\theta r}}=\mu \left[ r\frac{\partial }{\partial r}(\frac{{{v}_{\theta }}}{r})+\frac{1}{r}\times \frac{\partial {{v}_{r}}}{\partial \theta } \right]\)
\({{\tau }_{\theta z}}={{\tau }_{z\theta }}=\mu \left[ \frac{\partial {{v}_{\theta }}}{\partial z}+\frac{1}{r}\times \frac{\partial {{v}_{z}}}{\partial \theta } \right]\)
\({{\tau }_{zr}}={{\tau }_{r\theta }}=\mu \left[ \frac{\partial {{v}_{z}}}{\partial r}+\frac{\partial {{v}_{r}}}{\partial z} \right]\)
7. 直径为36.02cm的柱塞在直径为36.04cm的圆筒中滑动，构成了汽车的行程，中间的环状区域充满了油，油的运动粘性系数为0.00037m2/s，相对密度为0.85。如果柱塞运行的速度是0.15m/s，试估算3.14m的柱塞在圆筒中运动所形成的摩擦阻力。
8. 上题中如果柱塞和汽车滑轨的总质量为680kg，当仅有重力和粘性摩擦作用时，试估算柱塞和滑轨下沉速度的最大值。假设柱塞作用的长度为2.44m。
