1.1 流体力学基础

流体力学是连续介质力学的一个分支，是研究流体（包含气体及液体）现象以及相关力学行为的科学。

1.1.1 流体力学概述

1738年，伯努利在他的专著中首次采用了水动力学这个名词并作为书名；1880年前后出现了空气动力学这个名词；1935年以后，人们概括了这两方面的知识，建立了统一的体系，统称为流体力学。

在人们的生活和生产活动中随时随地都可遇到流体，因此流体力学是与人类日常生活和生产密切相关的。大气和水是最常见的两种流体，大气包围着整个地球，地球表面的70%是水面。大气运动、海水运动（包括波浪、潮汐、中尺度涡旋、环流等）乃至地球深处熔浆的流动都是流体力学的研究内容。

20世纪初，世界上第一架飞机出现以后，飞机和其他各种飞行器得到迅速发展。20世纪50年代开始的航天飞行，使人类的活动范围扩展到其他星球和银河系。航空航天事业的蓬勃发展是同流体力学的分支学科——空气动力学和气体动力学的发展紧密相联的。这些学科是流体力学中最活跃、最富有成果的领域。

石油和天然气的开采、地下水的开发利用，要求人们了解流体在多孔或缝隙介质中的运动，这是流体力学分支之一——渗流力学研究的主要对象。渗流力学还涉及土壤盐碱化的防治，化工中的浓缩、分离和多孔过滤，燃烧室的冷却等技术问题。

燃烧离不开气体，这是有化学反应和热能变化的流体力学问题，是物理化学流体动力学的内容之一。爆炸是猛烈的瞬间能量变化和传递过程，涉及气体动力学，从而形成了爆炸力学。

沙漠迁移、河流泥沙运动、管道中的煤粉输送、化工中气体催化剂的运动等，都涉及流体中带有固体颗粒或液体中带有气泡等问题，这类问题是多相流体力学研究的范围。

等离子体是自由电子、带等量正电荷的离子以及中性粒子的集合体。等离子体在磁场作用下有特殊的运动规律。研究等离子体运动规律的学科称为等离子体动力学和电磁流体力学，它们在受控热核反应、磁流体发电、宇宙气体运动等方面有广泛的应用。

风对建筑物、桥梁、电缆等的作用使它们承受载荷和激发振动，废气和废水的排放造成环境污染，河床冲刷迁移和海岸遭受侵蚀，研究这些流体本身的运动及其同人类、动植物间的相互作用的学科称为环境流体力学（其中包括环境空气动力学、建筑空气动力学）。这是一门涉及经典流体力学、气象学、海洋学和水力学、结构动力学等学科的新兴边缘学科。

生物流变学研究人体或其他动植物中有关的流体力学问题。例如血液在血管中的流动，心、肺、肾中的生理流体运动和植物中营养液的输送。此外，还研究鸟类在空中的飞翔，动物在水中的游动等。

因此，流体力学既包含自然科学的基础理论，又涉及工程技术科学方面的应用。

目前，研究流体力学问题的方法有理论分析研究、实验模拟研究和数值模拟方法研究3种。

流体力学理论分析的一般过程是：建立力学模型，用物理学基本定律推导流体力学数学方程，用数学方法求解方程，然后检验和解释求解结果。理论分析结果能揭示流动的内在规律，物理概念清晰，物理规律能公式化，具有普遍适用性，但分析范围有限，只能分析简单的流动。而且，线性问题能得到结果，非线性问题分析非常困难。

实验研究的一般过程是：在相似理论的指导下建立模拟实验系统，用流体测量技术测量流动参数，处理和分析实验数据。

典型的流体力学实验有风洞实验、水洞实验、水池实验等。测量技术有热线、激光测速，粒子图像、迹线测速，高速摄影，全息照相，压力密度测量等。现代测量技术在计算机、光学和图像技术配合下，在提高空间分辨率和实时测量方面已取得长足进步。

实验结果能反映工程中的实际流动规律，发现新现象，检验理论结果等，现象直观，测试结果可靠。但流体的实验研究对测试设备要求较高，设计制造周期长，且调试复杂。实验研究的方法只能得到有限的实验数据，真实模拟物理问题比较困难。

数值研究的一般过程是：对流体力学数学方程进行简化和数值离散化，编制程序进行数值计算，将计算结果与实验结果比较。

常用的数值模拟方法有有限差分法、有限元法、有限体积法、边界元法、谱分析法等。计算的内容包括飞机、汽车、河道、桥梁、涡轮机等流场的计算，湍流、流动稳定性、非线性流动等的数值模拟。大型工程计算软件已成为研究工程流动问题的有力武器。数值模拟方法的优点是能计算理论分析方法无法求解的数学方程（适用于线性和非线性问题），能处理各种复杂流动问题，比实验方法省时省钱。但毕竟是一种近似解方法，适用范围受数学模型的正确性和计算机的性能所限制。

流体力学的3种研究方法各有优缺点，在实际研究流体力学问题时，应结合实际问题，取长补短，互为补充和印证。

1.1.2 连续介质模型

如固体一样，流体也是由大量的分子组成的，而分子间都存在比分子本身尺度大得多的间隙，同时，由于每个分子都在不停地运动，因此，从微观的角度看，流体的物理量在空间分布上是不连续的，且随时间不断变化。

在流体力学中仅限于研究流体的宏观运动，其特征尺寸（如日常见到的是m、cm、mm那样的量级）比分子自由程大得多。描述宏观运动的物理参数，是大量分子的统计平均值，而不是个别分子的值。在这种情形下，流体可近似用连续介质模型处理。

连续介质模型认为，物质连续地分布于其所占有的整个空间，物质宏观运动的物理参数是空间及时间的可微连续函数。

根据连续介质模型假设，可以把流体介质的一切物理属性，如密度、速度、压强等都看作是空间的连续函数。因而，对于连续介质模型，微积分等现代数学工具可以加以应用。

连续介质模型假设成立的条件是建立在流体平均自由程远远小于物体特征尺寸的基础上的，即

在式1-1中，L为求解问题中物体或空间的特征尺寸；λ为流体分子的平均自由程。

在某些情况下，例如，在120km的高空，如果空气分子的平均自由程和飞行器的特征尺寸在同一数量级，连续介质模型假设就不再成立。这时，必须把空气看成是不连续的介质，这个范围属于稀薄空气动力学范畴。

1.1.3 流体的基本概念及性质

1．密度

流体的密度定义为单位体积所含物质的多少，以ρ表示。密度是流体的一种固有物理属性，国际单位为kg m3。对于均质流体，设其体积为V，质量为m，则其密度为

对于非均质流体，不同位置的密度不同。若取包含某点的流体微团，其体积为ΔV，质量为Δm，则该点的密度定义为

液体和气体的密度可以查询相关手册，这里为了方便读者，表1-1给出部分常见液体和气体的密度。

2．质量力和表面力

作用在流通微团上的力可分为质量力与表面力。

与流体微团质量大小有关并且集中作用在微团质量中心上的力称为质量力，如在重力场中的重力mg、直线运动的惯性力ma等。质量力是一个矢量，一般用单位质量所具有的质量力来表示，其形式为

式1-4中，fx、fy、fz为单位质量力在x、y、z轴上的投影，或简称为单位质量分力。

大小与表面面积有关而且分布作用在流体表面上的力称为表面力。表面力按其作用方向可分为两种：一种是沿表面内法线方向的压力，称为正压力；另一种是沿表面切向的摩擦力，称为切向力。

作用在静止流体上的表面力只有表面内法线方向的正压力，单位面积上所受到的表面力称为这一点处的静压强。静压强有两个特征：

● 静压强的方向垂直指向作用面。

● 流场内一点处静压强的大小与方向无关。

对于理想流体流动，流体质点只受到正压力，没有切向力。

对于黏性流体流动，流体质点所受到的表面力既有正压力，也有切向力。单位面积上所受到的切向力称为切应力。对于一元流动，切应力由牛顿内摩擦定律给出；对于多元流动，切应力可由广义牛顿内摩擦定律求得。

3．绝对压强、相对压强和真空度

单位面积上受到的压力叫作压强。在流体静力学中压强的定义是：作用在浸没于流体中物体表面上单位面积上的法向正应力。

压强在国际单位制中的单位是N/m2或者Pa（Pascal或帕斯卡）。另外压强也常用液柱高（汞柱、水柱等）、标准大气压（atm）和bar等单位进行度量，常用转换关系如下：

一个标准大气压的压强是760mmHg，相当于101325 Pa，通常用patm表示。

若压强大于标准大气压，则以标准大气压为计算基准得到的压强称为相对压强，也称为表压强（Gauge Pressure），通常用pr表示。

在FLUENT软件中求解器运算过程中实际使用的压强值都是表压强。FLUENT中默认的参考压力为一个标准大气压，用户可以指定参考压力。

若压强小于标准大气压，则压强低于大气压的值就称为真空度，通常用pv表示。

绝对压强ps、相对压强pr和真空度pv之间的关系为

所以在FLUENT中，如果用户将参考压力设置为0，表压强就等于绝对压强，有时候这么做可以方便边界条件的设置。

在流体力学中有如下约定：对于液体来说，压强用相对压强；对于气体，特别是马赫数大于0.3的流动，压强用绝对压强。

4．静压、动压和总压

物体在流体中运动时，在正对流体运动的方向的表面，流体完全受阻，此处的流体速度为0，其动能转变为压力能，压力增大，其压力称为全受阻压力（简称全压或总压），它与未受扰动处的压力（即静压）之差称为动压，即

式1-7中，pt为总压，p为静压，为动压，ρ为流体密度，v为流体速度。也可以表达为：对于不考虑重力的流动，总压就是静压和动压之和。

根据伯努利方程的物理意义可知，对于一条理想流体，在一条流线上流体质点的机械能是守恒的，不可压缩流动的表达式为

式1-8中，称为压强水头，也是压能项，p为静压；为速度水头，也是动能项；第三项z为位置水头，也是重力势能项。这三项之和就是流体质点的总机械能，等式右边的H称为总水头。

将式1-8的左右两边同时乘以ρg，则有

5．黏性

黏性是施加于流体的应力和由此产生的变形速率以一定的关系联系起来的流体的一种宏观属性，表现为流体的内摩擦。由于黏性的耗能作用，在无外界能量补充的情况下，运动的流体将逐渐停止下来。

黏性对物体表面附近的流体运动产生重要作用，使流速逐层减小并在物面上为零，在一定条件下也可使流体脱离物体表面。

黏性又称黏性系数、动力黏度，记为 μ。牛顿黏性定律指出，在纯剪切流动中，相邻两流体层之间的剪应力（或黏性摩擦应力）为式中，τ为剪切应力，为垂直流动方向的法向速度梯度。黏性数值上等于单位速度梯度下流体所受的剪应力。

速度梯度也表示流体运动中的角变形率，故黏性也表示剪应力与角变形率之间比值关系。按国际单位制，黏性的单位为kg/(m·s)。

黏性系数与密度之比称为运动黏性系数，常记作ν，则有

液体的黏性系数随温度的增加而减小，气体的黏性系数随温度的增加而变大。对于气体而言，黏性系数与温度的关系可以用萨瑟兰公式表示：

式1-12中，μ0和T0分别为参考黏性系数和参考温度。

这里为了方便读者，表1-2列出了部分常见液体和气体的动力黏度和运动黏度。

6．传热性

当气体中沿某一方向存在温度梯度时，热量就会由温度高的地方传向温度低的地方，这种性质称为气体的传热性。实验证明，单位时间内所传递的热量与传热面积成正比，与沿热流方向的温度梯度成正比，即

式中，q 为单位时间通过单位面积的热量，负号表示热流量传递的方向永远和温度梯度的方向相反。

流体的导热系数λ随流体介质而不同，同一种流体介质的导热系数随温度变化而略有差异。

7．扩散性

当流体混合物中存在组元的浓度差时，浓度高的地方将向浓度低的地方输送该组元的物质，这种现象称为扩散。

流体的宏观性质，如扩散、黏性和热传导等，是分子输运性质的统计平均。由于分子的不规则运动，在各层流体间交换着质量、动量和能量，使不同流体层内的平均物理量均匀化。这种性质称为分子运动的输运性质。质量输运在宏观上表现为扩散现象，动量输运表现为黏性现象，能量输运则表现为热传导现象。

理想流体忽略了黏性，即忽略了分子运动的动量输运性质，因此，在理想流体中也不应考虑质量和能量输运性质——扩散和热传导，它们具有相同的微观机制。

8．流线和迹线

所谓流线，就是这样一种曲线，在某时刻，曲线上任意一点的切线方向正好与那一时刻该处的流速方向重合。可见，流线是由同一时刻的不同流点组成的曲线，它给出了该时刻不同流体质点的速度方向，是速度场的几何表示。

所谓迹线，就是流体质点在各时刻所行路径的轨迹线（或流体质点在空间运动时所描绘出来的曲线），如喷气式飞机飞过后留下的尾迹、台风的路径、纸船在小河中行走的路径等。

流线具有以下性质。

（1）同一时刻的不同流线，不能相交。

（2）流线不能是折线，而是一条光滑的曲线。

（3）流线族的疏密反映了速度的大小。

（4）实际流场中，除驻点或奇点外，流线不能突然转折。

流线的微分方程为

式1-14中，u(x, y, z, t)、v(x, y, z, t)和w(x, y, z, t)分别表示点(x, y, z)在t时刻的速度在x、y和z方向上的分量。

迹线的微分方程为

式1-15中，t为关于时间的自变量，x(t)、y(t)、z(t)为t时刻此流体质点的空间位置。u、v和w分别为流体质点速度在x、y和z方向上的分量。

可见，对于定常流动，流线的形状不随时间变化，而且流体质点的迹线与流线重合。

9．流量和净通量

（1）流量。单位时间内流过某一控制面的流体体积称为该控制面的流量Q，其单位为m3/s。若单位时间内流过的流体是以质量计算的，则称为质量流体 Qm，不加说明时，“流量”一词概指体积流量。在曲面控制面上有

（2）净流量。在流场中取整个封闭曲面作为控制面A，封闭曲面内的空间称为控制体。流体经一部分控制面流入控制体，同时也有流体经另一部分控制面从控制体中流出。此时流出的流体减去流入的流体，所得出的流量称为流过全部封闭控制面A的净流量（或净通量），计算式为

对于不可压流体来说，流过任意封闭控制面的净通量等于0。

10．流速和音速

当把流体视为可压缩流体时，扰动波在流体中的传播速度是一个特征值，称为音速。音速方程式的微分形式为

音速在气体中的传播过程是一个等熵过程。将等熵方程式p=c ρk代入式1-18，并由理想气体状态方程p=ρRT，得到音速方程为

对于空气来说，k=1.4, R=287 J/(kg·K)，得到空气中的音速为

流速是流体流动的速度，而音速是扰动波的传播速度，两者之间的关系为

式1-21中，Ma称为马赫数。

11．马赫数和马赫锥

（1）马赫数。流体流动速度v与当地音速c之比称为马赫数，用Ma表示为

Ma<1的流动称为亚音速流动，Ma>1的流动称为超音速流动，Ma>3的流动称为高超音速流动。

（2）马赫锥。对于超音速流动，扰动波传播范围只能允许充满在一个锥形的空间内，这就是马赫锥，其半锥角θ称为马赫角，计算式如下：

12．滞止参数和等熵关系

（1）滞止参数。流场中速度为0的点上的各物理量称为滞止参数，用下标“0”表示，如滞止温度T0、滞止压强p0、滞止密度ρ0等。

（2）等熵流动基本关系式。流动参数与滞止参数及马赫数之间的基本关系为

式1-24的使用条件是绝热流动，式1-25和式1-26的使用条件是等熵流动。

13．正激波和斜激波

气流发生参数发生显著、突跃变化的地方，称为激波。激波常在超音速气流的特定条件下产生；激波的厚度非常小，约为10-4mm，因此一般不对激波内部的情况进行研究，所关心的是气流经过激波前后参数的变化。气流经过激波时受到激烈的压缩，其压缩过程是很迅速的，可以看作是绝热的压缩过程。

（1）正激波。

激波面与气流方向垂直，气流经过激波后方向不变，这称为正激波。

假设激波固定不动，激波前的气流速度、压强、温度和密度分别为v1、p1、T1和ρ1，经过激波后突跃地增加到v2、p2、T2和 ρ2。设激波前气流马赫数为Ma1，则激波前后气流应满足如下公式。

连续性方程

动量方程

能量方程（绝热过程）

状态方程

由可将式1-30改写成

在以上几个基本关系式的基础上，可导出以下的重要关系式。

（2）斜激波。

气流经过激波后方向发生改变，这种激波称为斜激波。v1t、v2t和v1n、v2n各表示斜激波前后速度v1和v2的切向分速度和法向分速度，α为气流折转角，β为激波角。

由于沿激波面没有切面压强的变化，所以气流经过斜激波后，沿激波面的分速度没有变化，即有

因为发生变化的只有法向分速度，所以斜激波相当于法向分速度的正激波。

由于斜激波是相当于法向分速度的正激波，所以只要把正激波关系式中各下标“1”、“2”换成“1n”、“2n”，则正激波有关方程式可以应用于斜激波。

对于斜激波前后气流参数之间的关系，有如下的关系式：

激波角β与气流折转角α之间满足如下的关系式：

1.1.4 流体流动分类

1．理想流体与黏性流体

如果忽略流动中流体黏性的影响，则可以近似地把流体看成是无黏的，称为无黏流体（inviscid liquid），也叫作理想流体。这时的流动称为理想流动，理想流体中没有摩擦，也就没有耗散损失。

事实上，真正的理想流体是不存在的。但是在一定的情形下，至少在特定的流动区域中，某些流体的流动非常接近于理想的流动条件，在分析处理中可以当作理想流体。

例如，在空气绕物体的流动（空气动力学）中，除去邻近与物体表面的薄层（称为边界层）之外，在其余的流动区域中，空气动力学中都处理成理想流动，此时所求解的控制方程组是不考虑黏性的欧拉方程组。

2．牛顿流体与非牛顿流体

根据内摩擦剪应力与剪应力变率的关系不同，黏性流体又可分为牛顿流体与非牛顿流体。

如果流体的剪应力与剪应力变率遵守牛顿内摩擦定律，即公式，则这种流体就称为牛顿流体。尽管这个线性的牛顿关系式只是一种近似，但是却很好地适用于一类范围很广的流体。水、空气和气体等绝大多数工业中常用的流体都是牛顿流体。

但是，对于某些物质而言，剪应力不只是速度梯度的函数（速度梯度和剪应变率是相同的），通常还可以是应变的函数，这种物质称为黏-弹性流体。剪应力只依赖于速度梯度的简单黏性流体，也可以不是牛顿流体。

事实上存在这样的流体，其剪应力与应变率之间有着相当复杂的非线性关系。如果流体的应力-应变关系还取决于实际的工况，即应变工况，则称为触变流体（如印刷油墨）。

非牛顿流体具有塑性行为，其特征是有一个表现的屈服应力，在达到表现的屈服应力之前，流体的性态像固体一样，一旦超过这个表现的屈服应力，则和黏性流体一样。

塑性流体的另一个极端情形是：在低应变率时，黏性系数很小，很容易流动，但随着应变率的增加，变得越来越像固体（如流沙）。这种流体称为膨胀流体。在图1-1中用曲线说明了牛顿流体和非牛顿流体的特征。

3．可压流体与不可压流体

流体的压缩性是指在外界条件变化时，其密度和体积发生了变化。这里的条件有两种，一种是外部压强发生了变化，另一种就是流体的温度发生了变化。描述流体的压缩性常用以下两个量。

（1）流体的等温压缩率β。

当质量为m，体积为V的流体外部压强发生Δp的变化时，相应其体积也发生了ΔV的变化，则定义流体的等温压缩率为

这里的负号是考虑到ΔV与Δp总是符号相反的缘故；β 的单位为1/Pa。流体等温压缩率的物理意义是，当温度不变时，每增加单位压强所产生的流体体积的相对变化率。

考虑到压缩前后流体的质量不变，式1-42还有另外一种表示形式，即

将理想气体状态方程代入式1-43，得到理想气体的等温压缩率为

（2）流体的体积膨胀系数α。

当质量为m，体积为V的物体温度发生ΔT的变化时，相应其体积也发生了ΔV的变化，则定义流体的体积膨胀系数为

考虑到膨胀前后流体的质量不变，式1-45还有另外一种表示形式，即

这里的负号是考虑到随着温度的升高，体积必然增大，则密度必然减少；α的单位为1/K。体积膨胀系数的物理意义是，当压强不变时，每增加单位温度所产生的流体体积的相对变化率。

对于物理气体，将气体状态方程代入式1-46，得到

在研究流体流动过程时，若考虑到流体的压缩性，则称为可压缩流动，相应地称流体为可压缩流体，如马赫数较高的气体流动。若不考虑流体的压缩性，则称为不可压缩流动，相应地称流体为不可压缩流体，如水、油等液体的流动。

4．定常流动与非定常流动

根据流体流动过程以及流动过程中的物理参数是否与时间相关，可将流动分为定常流动与非定常流动两种。

流体流动过程中，各物理量均与时间无关的流动称为定常流动。

流体流动过程中，某个或某些物理量与时间有关的流动称为非定常流动。

5．层流流动与湍流流动

流体的流动分为层流流动和湍流流动，从试验的角度来看，层流流动就是流体层与层之间相互没有任何干扰，层与层之间既没有质量的传递也没有动量的传递；而湍流流动中层与层之间相互有干扰，而且干扰的力度还会随着流动的加速而加大，层与层之间既有质量的传递，又有动量的传递。

判断流动是层流还是湍流，是看其雷诺数是否超过临界雷诺数。雷诺数的定义如下

式中，v为截面的平均速度，L为特征长度，ν为流体的运动黏度。

对于圆形管内流动，特征长度L取圆管的直径d，即

一般认为临界雷诺数为2000，当Re<2000时，管内流动是层流，否则为湍流。

对于异型管道内的流动，特征长度L取水力直径dH，则雷诺数的计算式为

异型管道水力直径的定义为

式中，A为过流断面的面积；S为过流断面上流体与固体接触的周长，称为湿周。

1.1.5 流体流动描述的方法

描述流体物理量有两种方法，一种是拉格朗日描述，另一种是欧拉描述。

拉格朗日（Lagrange）描述也称随体描述，它着眼于流体质点，并将流体质点的物理量认为是随流体质点及时间变化的，即把流体质点的物理量表示为拉格朗日坐标及时间的函数。

设拉格朗日坐标为（a, b, c），以此坐标表示的流体质点的物理量，如矢径、速度、压强等在任意时刻t的值，便可以写为a、b、c及t的函数。

若以f表示流体质点的某一物理量，其拉格朗日描述的数学表达式为

例如，设时刻t流体质点的矢径，即t时刻流体质点的位置以r表示，其拉格朗日描述为

同样，质点速度的拉格朗日描述为

欧拉描述也称空间描述，它着眼于空间点，认为流体的物理量随空间点及时间而变化，即把流体物理量表示为欧拉坐标及时间的函数。设欧拉坐标为（q1, q2, q3），用欧拉坐标表示的各空间点上的流体物理量，如速度、压强等，在任意时刻t的值，可写为q1、q2、q3及t的函数。

从数学分析知道，当某时刻一个物理量在空间的分布一旦确定，该物理量在此空间形成一个场。因此，欧拉描述实际上描述了一个个物理量的场。

若以f表示流体的一个物理量，其欧拉描述的数学表达式为（设空间坐标取用直角坐标）

如流体速度的欧拉描述为

拉格朗日描述着眼于流体质点，将物理量视为流体坐标与时间的函数；欧拉描述着眼于空间点，将物理量视为空间坐标与时间的函数。它们可以描述同一物理量，必定互相相关。

设表达式f =f(a, b, c, t)表示流体质点（a, b, c）在t时刻的物理量；表达式f=F(x, y, z, t)表示空间点（x, y, z）在时刻t的同一物理量。如果流体质点（a, b, c）在t时刻恰好运动到空间点（x, y, z）上，则应有

事实上，将式1-57代入式1-58左端，即有

或者反解式1-57，得到

将式1-60代入式1-58的右端，也应有

由此，可以通过拉格朗日描述推出欧拉描述，同样也可以由欧拉描述推出拉格朗日描述。

1.1.6 流体力学基本方程组

流体力学基本方程组包括质量守恒方程、动量守恒方程、组分质量守恒方程、能量守恒方程、本构方程、状态方程及通用形式守恒方程。

下面对各种形式方程进行总结和对比，并分析它们之间的转化关系，以帮助读者理解流体力学基本方程组的数学物理意义，为离散计算这些方程组打下基础。

1．质量守恒方程

质量守恒方程可由4种方法得到，分别在拉格朗日法（L法）下对有限体积和体积元应用质量守恒定律、在欧拉法（E法）下对有限体积应用质量守恒定律及在直角坐标系中直接应用质量守恒定律。

（1）L法有限体积分析。

取体积为τ，质量为m的一定流体质点团，则有

因为速度散度的物理意义是相对体积膨胀率及密度的随体导数，即

代入式1-62，得

运用奥高定理

得

式1-67即是连续性方程的积分形式。

假定被积函数连续，而且体积τ是任意选取的，由此可知被积函数必须等于0，即

或

在直角坐标系中，连续性方程为

或

式1-71表明，密度变化（随时间和位置）等于密度和体积变形的乘积。

（2）L法体积元分析。

考虑质量为dm的体积元dτ，对其用拉格朗日观点，根据质量守恒定律有

两边同除以ρdτ，得

或写成

式1-75表明，要维持质量守恒定律，相对体积变化率必须等于负的相对密度变化率。

（3）E法有限体积分析。

着眼坐标空间，取空间中以S面为界的有限体积τ，则称S面为控制面，τ为控制体。取外法线方向为法线的正方向，n 为外法线方向的单位矢量。考虑该体积内流体质量的变化，该变化主要由以下两方面原因引起。

① 通过表面S有流体流出或流入，单位时间内流出流入变化的总和为

② 由于密度场的不定常性（注意，欧拉观点下空间点是固定的，密度的变化只由场的不定常性刻画），单位时间内体积τ的质量将变化，变化量为

上述两者应相等，即

由于体积τ是任意的，且被积函数连续，则

（4）E法直角坐标系分析。

控制体体积如图1-2所示，单位时间内通过表面EFGH的通量为ρudydz。

通过表面ABCD的通量为

其他两对表面类似。另外，该控制体内质量的变化率为

则

特殊情况下的连续性方程为

定常态

不可压缩流体

2．动量守恒方程

任取一个体积为τ的流体，它的边界为S。根据动量定理，体积τ中流体动量的变化率等于作用在该体积上的质量力和面力（应力）之和。

单位面积上的面力pn=n · P，其中P是二阶对称应力张量，所以pn不是通常指的 p在 n（单位体积面元的法线方向）方向的分量。单位质量上的质量力为F，则作用在该体积上的质量力和面力分别为

及

动量变化率为

上述动量变化率的表达式可采用如下两种处理方法。

（1）求解式1-85右边第二项内对体积元的随体导数，则

（2）对动量变化率表达式1-85右边第二项应用质量守恒定律

由上可得两种积分形式的动量方程，即

或

动量方程的微分形式为

或

微分方程中各项的物理意义为，表示单位体积上惯性力，ρF为单位体积上的质量力，divP为单位体积上应力张量的散度，它是与面力等效的体力分布函数（由奥高公式转化而来）。

在直角坐标系下以应力表示的运动方程可采取下列形式。

或

这两种表达方式的等号左边实际只差了一个连续性方程，由基本微分公式

得

由连续性方程知

所以有

上述运动方程是以应力表示的黏性流体的运动方程，它们对任何黏性流体和任何运动状态都是适用的。但它们没有反映出不同属性的流体受力后的不同表现。

另外，方程数和未知量之数不等，运动方程有3个，加上连续性方程共4个，但未知量却有9个（6个应力张量分量（9个张量分量因对称关系减少为6个）和3个速度分量），所以该方程组不封闭。为使该方程组可解，必须考虑应力张量和变形速度张量之间的关系（将应力张量用速度分量表示出来），补足所需的方程。

3．本构方程

本构方程是表征流体宏观性质的一种微分方程，它用于表达流体黏性定律的应力张量和变形速度张量之间的关系。

最简单的应力与应变之间的关系是对牛顿流体作一维运动，即牛顿剪切定律可写为

要得到普遍意义上的广义牛顿定律需做一定的假设，但首先应理解流体速度分解定理和变形速度张量。

（1）流体速度分解定理。

刚体运动包括平动和转动两部分，一般可表为

其中v0是刚体中选定一点O上的平动速度，ω是刚体绕O点转动的瞬时角速度矢量，r是要确定速度那一点到O点的矢径。转动角速度可用v表示。因为

故

流体运动除平动、转动外，还有变形运动。设微团内M0点的速度为v0，邻域内任一点M的速度为v。将v在M0点泰勒展开并略去二阶无穷小项，得

显然，是一个二阶张量（局部速度梯度张量），由张量分解定理可将该张量分解成对称张量S和反对称张量A之和，于是

所以

式1-104右边第二项、第三项可具体表示为

及

其中

另外

所以

式1-109表明流体运动可分为平动、转动和变形3种形式，S称为变形速度张量，该定理称为亥姆霍兹（Helmholtz）速度分解定理。

另外，流变学中常用应变速率张量．来表示流体的变形和拉伸（或压缩），而用转动张量Ωij表示转动，它们与流体力学中的变形速度张量和转动张量的关系为

（2）变形速度张量的物理意义。

变形速度v3的表达式为

经分析可得

其中、及是角变形速率，亦称剪切应变速率（称拉伸应变速率）。

由上可知，变形速度张量的对角线分量ε1、ε2、ε3的物理意义分别是x、y、z轴线上线段元dx、dy、dz的相对拉伸速度或相对压缩速度。

而非对角线分量θ1、θ2、θ3的物理意义分别是y轴与z轴、z轴与x轴、x轴与y轴之间夹角的剪切速率的负值。

（3）广义牛顿定律及基本假设。

① 运动流体的应力张量在运动停止后应趋于静止流体的应力张量。据此将应力张量P写成各向同性部分-pI和各向异性部分P’是方便的，因此

P ＇是除去-pI后得到的张量，称为偏应力张量。当运动消失时它趋于0。可见，偏应力张量和应力张量一样也是对称张量。

② 偏应力张量的各分量τij是局部速度梯度张量各分量的线性齐次函数。当速度在空间均匀分布时，偏应力张量为0；当速度偏离均匀分布时，在黏性流体中产生了偏应力，它力图使速度恢复到均匀分布情形。

③ 流体是各向同性的，即流体性质不依赖于方向或坐标系的转换。

根据假设②，有

显然cijkl是一个四阶张量，它是表征流体黏性的常数，共34=81个。根据假设③, cijkl是各向同性张量且对i、j对称，故

观察式1-117可知，cijkl对k、l也是对称的，物性常数减少至只有2个，即第二黏度λ和黏度μ。将式1-117代入偏应力表达式1-116（反对称项为零），得

则应力张量为

引进

则

根据式1-121，可知

将式1-122、式1-123和式1-124相加，得

对于可压缩流体，流体的体积在运动过程中发生膨胀或收缩，将引起平均法应力（由奥高公式可证某固定点处所有方向上法应力的平均值等于x、y、z三个方向上法应力的平均值，这是一个不随坐标系改变的不变量）的值发生μ＇divv的改变，μ＇称为第二黏性系数，亦称膨胀黏性系数。应用斯托克斯假定，即μ＇=0，则本构方程为

一般处理的是不可压缩流体，则

在直角坐标系下有

这里pij=σij，有的文献中将应力张量用σij表示。将上述应力张量与变形速度张量的关系式代入运动方程，得

即

写成直角坐标系下的形式为

有些文献中常把式1-136等号右边表示分子黏性作用的3项做如下变化，以第一式为例。

其中

据此，有

及

其中广义源项定义为

当流体黏度不变且不可压缩时（牛顿流体），有

所以运动方程简化为

其中υ是运动黏度，亦是动量扩散系数，单位为m2/s。

本构方程和运动方程是紧密联系在一起的，通过本构方程可将应力张量用变形速度张量表示出来，即应力可用应变速率表示，而应变速率实际由速度分量决定，故使运动方程和连续性方程原则上封闭可解。

4．能量守恒方程

由能量守恒定律知，体积τ内流体的动能和内能的变换率等于单位时间内质量力和表面力所做的功加上单位时间内给予体积τ的热量。

体积τ内流体的动能和内能的总和为

其中U是单位体积内的流体内能。质量力对体积τ内流体所做的功（单位时间内移动距离v，点积求做功）为

表面力对体积τ内流体所做的功（单位时间内移动距离v，点积求做功）为

单位时间内以热传导方式通过表面S传给体积τ的热量为

式1-147中的被积函数实际就是傅里叶热传导定律，即热流密度矢量正比于传热面法向温度梯度。单位时间内由于辐射或其他原因（反应、蒸发等）传入τ内的总热量（q为单位时间内传入单位质量的热量分布函数）为

将式1-144～式1-148进行守恒计算，得

这是积分形式的能量守恒方程。求解体积分的随体导数并运用奥高公式把面积分转化为体积分，可得到微分形式的能量守恒方程，即

由质量守恒定律可得

所以

另外

则能量方程的微分形式为

或

或

式1-157各项的物理意义为：左边第一、第二项代表内能和动能的随体导数，右边第一项是单位体积内的质量力做的功，第二项是单位体积内面力所做的功，第三项是单位体积内热传导输入的热量，最后一项表示由于辐射或其他物理或化学原因的热量贡献。

能量守恒方程的另一种形式为

式1-158的物理意义为：单位体积内由于流体变形，面力所做的功加上热传导及辐射等其他原因传入的热量恰好等于单位体积内的内能在单位时间内的增加。将该式进一步简化，有

设

φ为由于黏性作用机械能转化为热能的部分，称为耗散函数。另外，在考虑液体流体时，比焓h与内能值可看作相等，即h=U+pV=U，压力不做功，则

所以有

其中S h=ρq是单位体积内热源或由于辐射或其他物理或化学原因的热量贡献。p div v一般较小，可以忽略。对液体及固体可以取h=cTp，进一步取cp为常数，并把耗散函数纳入源项ST=Sh+φ，于是

对于不可压缩流体，有

对于可以忽略黏性耗散作用的稳态低速流，能量方程可以简化为

及

取速度为0，则可得到稳态的热传导方程（对流项消失）为

5．状态方程

由连续性方程、运动方程、能量方程确定的未知量有u、v、w、p、T、ρ共6个，但方程数只有5个，为使方程组封闭，需补充一个联系p、ρ的状态方程：

6．组分质量守恒方程

在一个特定的系统中可能存在质的交换，或者存在多种化学组分，每一种组分都需要遵守组分质量守恒定律，即系统内某种化学组分对时间的变化率等于通过系统界面的净扩散流量与由反应产生的生成率之和，可表示为

其中代表单位体积内组分l的质量变化率，ρvml是组分l的对流流量密度。Jl代表扩散流量密度，它由费克定律给出。Sl是单位体积内组分l的生成率。费克定律可表示为

其中Dl为扩散系数。将扩散定律代入守恒方程，得

为便于以后引用，表1-3列出了本节的守恒型控制方程。

1.1.7 湍流模型

描述流体运动（层流）的流体力学基本方程组是封闭的，而描述湍流运动的方程组由于采用了某种平均（时间平均或网格平均等）而不封闭，必须对方程组中出现的新未知量采用模型而使其封闭，这就是CFD中的湍流模型。

湍流模型的主要作用是将新未知量和平均速度梯度联系起来。目前，工程应用中湍流的数值模拟主要分三大类：直接数值模拟（DNS）、大涡模拟（LES）和基于雷诺平均N-S方程组（RANS）的模型。

1．直接数值模拟（DNS）

直接数值模拟（DNS）方法直接求解湍流运动的N-S方程，得到湍流的瞬时流场，即各种尺度的随机运动，可以获得湍流的全部信息。

随着现代计算机的发展和先进数值方法的研究，DNS方法已经成为解决湍流的一种实际的方法。但由于计算机条件的约束，目前只能限于一些低Re数的简单流动，不能用于处理工程中的复杂流动问题。

目前国际上正在做的湍流直接数值模拟还只限于较低的雷诺数（Re约200）和非常简单的流动外形，如平板边界层、完全发展的槽道流以及后台阶流动等。

2．大涡模拟（LES）

大涡模拟（LES）方法即对湍流脉动部分的直接模拟，将N-S方程在一个小空间域内进行平均（或称之为滤波），以便从流场中去掉小尺度涡，导出大涡所满足的方程。

小涡对大涡的影响会出现在大涡方程中，再通过建立模型（亚格子尺度模型）来模拟小涡的影响。

由于湍流的大涡结构强烈地依赖于流场的边界形状和边界条件，难以找出普遍的湍流模型来描述具有不同的边界特征的大涡结构，所以宜做直接模拟。

相反地，由于小尺度涡对边界条件不存在直接依赖关系，而且一般具有各向同性性质，所以亚格子尺度模型具有更大的普适性，比较容易构造，这是它比雷诺平均方法要优越的地方。

LES方法已经成为计算湍流的最强有力的工具之一，应用的方向也在逐步扩展，但是仍然受计算机条件等的限制，使之成为解决大量工程问题的成熟方法仍有很长的路要走。

3．基于雷诺平均N-S方程组（RANS）的模型

目前能够用于工程计算的方法就是模式理论。所谓湍流模式理论，就是依据湍流的理论知识、实验数据或直接数值模拟结果，对雷诺应力做出各种假设，即假设各种经验的和半经验的本构关系，从而使湍流的平均雷诺方程封闭。

从对模式处理的出发点不同，可以将湍流模式理论分成两大类：一类称为二阶矩封闭模式（雷诺应为模式），另一类称为涡黏性封闭模式。

（1）雷诺应力模式。

雷诺应力模式即二阶矩封闭模式，是从雷诺应力满足的方程出发，将方程右端未知的项（生成项、扩散项、耗散项等）用平均流动的物理量和湍流的特征尺度表示出来。

典型的平均流动的变量是平均速度和平均温度的空间导数。这种模式理论由于保留了雷诺应力所满足的方程，如果模拟得好，可以较好地反映雷诺应力随空间和时间的变化规律，因而可以较好地反映湍流运动规律。

因此，二阶矩封闭模式是一种较高级的模式，但是，由于保留了雷诺应力的方程，加上平均运动的方程，整个方程组总计15个方程，是一个庞大的方程组，应用这样一个庞大的方程组来解决实际工程问题，计算量很大，这就极大地限制了二阶矩封闭模式在工程问题中的应用。

（2）涡黏性封闭模式。

在工程湍流问题中得到广泛应用的模式是涡黏性封闭模式。这是由Boussinesq仿照分子黏性的思路提出的，即设雷诺应力为这里是湍动能，νT称为涡黏性系数，这是最早提出的基准涡黏性封闭模式，即假设雷诺应力与平均速度应变率成线性关系。当平均速度应变率确定后，6个雷诺应力只需要通过确定一个涡黏性系数νT就可完全确定，且涡黏性系数各向同性，可以通过附加的湍流量来模化，如湍动能k、耗散率ε、比耗散率ω以及其他湍流量τ k/= ε、、。根据引入的湍流量的不同，可以得到不同的涡黏性模式，如常见的k-ε、k-ω模式，以及后来不断得到发展的k-τ、q-w、k-l等模式，涡黏性系数可以分别表示为

为了使控制方程封闭，引入多少个附加的湍流量，就要同时求解多少个附加的微分方程，根据求解的附加微分方程的数目，一般可将涡黏性模式划分为4类：零方程模式、半方程模型、一方程模式和两方程模式。

① 零方程模式。

所谓零方程模式，就是试图直接用平均流动物理量模化νT，而不引入任何湍流量（如k、ε等）。

例如，Prandtl的混合长理论就是一种零方程模式：

式中，l称为混合长。

在零方程模式的框架下，得到最为广泛应用的是Baldwin-Lomax模式（BL模式）。该模式是对湍流边界层的内层和外层采用不同的混合长假设。

这是因为靠近壁面处，湍流脉动受到很大的抑制，含能涡的尺度减小很多，因此长度尺度减小很多；另一方面，在边界层外缘，湍流呈间歇状，质量、动量和能量的输运能力大大下降，即湍流的扩散能力减小。

这样，应用混合长理论来确定涡黏性系数在这两个不同的区域应该有不同的形式。Baldwin-Lomax模式的具体数学描述如下。

这里yc是(vT)in=(vT)out的离壁面最小距离y值。

对于内层，即y≤yc，有

Ω是涡量，j, l是长度尺度。

式中，k=0.4是Karman常数，A+是模化常数，y+是无量纲法向距离。

其中uτ是摩擦速度，其含义为

此处下标w表示壁面。

对于外层，即y＞yc，有

其中

Fmax是函数F(y)的最大值。

而ymax是F(y)达到最大值的位置。Fkleb是Klebanoff间歇函数。

Udif是平均速度分布中最大值和最小值之差。

几个模化常数的值如下。

A+ =26.0 C=0.02668 Ckleb=0.3 Cwk=1.0 k=0.4

由上述模化关系可以看出，雷诺应力完全由当时当地的平均流参数用代数关系式决定。

平均流场的任何变化立刻为当地的湍流所感知，这表明零方程模式是一个平衡态模式，假定湍流运动永远处于与平均运动的平衡之中。

实际上对于大多数湍流运动而言并非如此，特别是对于平均流空间和时间有剧烈变化的情形，再有因为坐标y显式地出现在湍流模式中，零方程模式不具有张量不变性，所以将它应用到复杂几何外形的流动的数值模拟会带来困难。

当流动发生分离时，Baldwin-Lomax模式会遇到困难，这是因为在分离点和再附点附近，摩擦速度uτ为0，此时要引入一些人为的干涉来消除这些困难。

计算实践表明，只要流动是附体的，零方程模式一般都可以较好地确定压强分布，但是摩阻和传热率的估算不够准确，特别是当流动有分离和再附时。这是因为附体流压强分布对湍流应力不敏感。

总之，对于附体流动，如果只关心压强分布，应用零方程模式通常可以给出满意的结果，而且模式应用起来十分简便。但是对于计算摩阻的需求，零方程模式是不能满足要求的。对于有分离、再附等复杂流动，零方程模式是不适用的。

② 半方程模式。

为了能计算具有较强压强梯度，特别是较强逆压梯度的非平衡湍流边界层，Johnson-King于1985年提出了一个非平衡代数模式（JK模式），该模式仍采用涡黏性假设，把涡黏性的分布与最大剪切应力联系在一起，内层涡黏性与外层涡黏性分布用一个指数函数作为光滑拟合，外层涡黏性系数作为一个自由参数，由描述最大剪切应力沿流向变化的常微分方程来确定，此常微分方程是由湍流动能方程导出的，故此模式又称为半方程模式。

JK模式虽然仍采用涡黏性假设，却包含雷诺应力模式的特点。由于求解常微分方程比一方程、两方程模式中求解偏微分方程要简单、省时得多，故用JK模式的工作量只略高于通常平衡状态的零方程代数模式的工作量。

③ 一方程模式。

Baldwin-Barth（BB）模式是在两方程模型中，将某一导出的应变量作为基本物理量而得到的，应用此一方程模式可避免求解两方程时会遇到的某些数值困难。BB一方程模式所选择的导出应变量为“湍流雷诺数”Ret。

BB模式对计算网格的要求低，壁面的网格可以与采用BL代数模式的相当，而不像两方程k-ε模式那样要求壁面网格很细，这样就避免了在k-ε模式中流场求解的刚性问题。

Spalart-Allmaras（SA）模式与BB模式不同，它不是直接利用k-ε两方程模式加以简化而得，而是从经验和量纲分析出发，由针对简单流动再逐渐补充发展，进而适用于带有层流流动的固壁湍流流动的一方程模式。该模式中选用的应变量是与涡黏性νT相关的量，除在黏性次层外，与νT是相等的。

上述两种一方程模式具有相似的特点，它们不像代数模式那样需要分为内层模式、外层模式或壁面模式、尾流模式，同时也不需要沿法向网格寻找最大值，因此易于应用到非结构网格中；但由于在每个时间步长内，需要对整个流场求解一组偏微分方程，故比BL和JK模式更费机时。

④ 两方程模式。

常用的两方程模式有标准k-ε两方程模式、可实现型k-ε两方程模式、低雷诺数k-ε模式及k-ω两方程模式等。

● 标准k-ε两方程模式。

k-ε 模式是最为人所知和应用最广泛的两方程涡黏性模式，为积分到壁面的不可压缩/可压缩湍流的两方程涡黏性模式。

雷诺应力的涡黏性模式为

其中μ t为涡黏性（eddy viscosity）, Sij为平均速度应变率张量（mean-velocity strain-rate tensor）, ρ为流体密度，k为湍动能，δ ij为克罗内克算子（Kronecker delta）。涡黏性定义为湍动能k和湍流耗散率ε的函数。

基于量纲分析，涡黏性由流体密度ρ、湍流速度尺度k2和长度尺度来标度，衰减函数f由湍流雷诺数Ret=ρk2/(εμ)来模化。

湍流输运方程可表示成湍流能量输运方程1-185和能量耗散输运方程1-186。

其中右端项分别表示生成项（production term）、耗散项（dissipation term）和壁面项（wall term）。

模式中各常数的值如下。

Cμ=0.09 cε1=1. 45 cε2 =1. 92

σk=1. 0 σε=1. 3 Prt=0. 9

近壁衰减函数

壁面项

其中us为平行于壁面的流动速度。

积分到壁面的无滑移边界条件为

k=0 ε=0

● 可实现型k - ε两方程模式。

上述标准k - ε模式，对于高平均切变率流动会出现非物理的结果（如当Sk/ε＞3.7时，其中）。为了保证模式的可实现性，模式函数Cμ不应该是常数，而应当是平均切变率的函数。实验表明，对边界层流动和均匀切变流，Cμ的值是不同的。为此，根据可实现性对模式的约束条件，建议采用以下形式的Cμ。

式1-191中

是在以角速度ωk旋转的旋转坐标系中得到的平均旋转速率。

关系式1-191中唯一未确定的系数是A0。为简单起见，可以设其为常数。对边界层流动，可以取A0= 4.0；对于其他流动，A0的数值可以调节。

● 低雷诺数k - ε模式。

上述几种k - ε模式适用于高雷诺数情形。但是对于近壁区，湍流雷诺数很低，对湍流动力学而言，黏性效应非常重要，此时湍流雷诺数的效应必须加以考虑。研究摩阻的计算关注的恰恰是近壁区，因此低雷诺数k - ε模式的研究是十分重要的。

现将有关结果整理如下。

低雷诺数下的涡黏性和k - ε模式方程为

式中

所有模化常数如下。

其中

此处，fμ和f1、f2称为阻尼函数，是用来反映近壁区低雷诺数效应的一个经验公式，系数ai、＇ai见表1-4。

● k-ω两方程模式。

k-ω模式是积分到壁面的不可压缩/可压缩湍流的两方程涡黏性模式，最主要的文献来自Wilcox。

下面求解湍动能k和它的ω=ε/k, τ=k/ε,,,的对流输运方程。

雷诺应力的涡黏性模型为

这里μt为涡黏性，Sij为平均速度应变率张量，ρ为流体密度，k为湍动能，δij为克罗内克算子。涡黏性定义为湍动能k和比耗散率ω的函数。

k和ω的输运方程为

模式中各常数的值如下。

对于边界层流动，壁面无滑移边界条件为

这里y1为离开壁面第一个点的距离，且。

除以上介绍的FLUENT中常用的湍流模型外，还有众多诸如SST两方程模式、k-τ模型、q-ω 模型、双尺度两方程模式等湍流模型，在此不再一一详细介绍，感兴趣的读者可以查阅相关文献。

在实际求解中，选用什么模型要根据具体问题的特点来决定。选择的一般原则是精度高，应用简单，节省计算时间，同时也具有通用性。

不同软件所包含的湍流模型有略微区别，但常用的湍流模型一般CFD软件都包含。图1-3为常用的湍流模型及其计算量的变化趋势。