第 3 章:物理现象的数学描述(Mathematical Description of Physical Phenomena)
3.1 介绍(Introduction)
本章概述支配流体流动、传热与传质以及其他相关输运现象的守恒原理,并将这些原理所对应的物理定律翻译成数学关系——表达为偏微分方程(PDE)的形式,作为后续数值模拟的载体。作者首先推导连续性方程、动量方程与能量方程(统称为 Navier-Stokes 方程),分别对应质量、动量与总能量的守恒原理;随后建立一个针对任意标量、矢量或张量量的通用守恒方程。本章还考察这些方程中各项的数学性质,并介绍一种常见做法——将守恒方程写成无量纲形式——以及由此产生的、在工程参数研究中极为有用的若干无量纲数。
Navier-Stokes 方程已有超过一百七十年的历史;如今所说的 Navier-Stokes 方程已经不仅是线性动量守恒方程,而被用来统称质量、动量与能量三个守恒方程。这些方程可用来刻画广泛尺度与场景的流体流动——飓风内部、涡轮机械内部、机翼或潜艇周围的绕流、动脉与肺中的流动、泵与压气机内部——Navier-Stokes 方程都能描述这些现象。
3.2 流体的分类(Classification of Fluid Flows)
流体指液体与气体,是在大的应力(单位面积上的力)作用下不会发生永久形变的物质。固体能通过形变来抵抗所施加的剪切或切向应力,而流体不能;施加于流体的剪切应力将使流体进入运动状态。与固体具有明确形状不同,流体没有确定的形状。气体是完全充满所占据区域的流体,而液体则在重力场作用下形成自由表面。
在分析流体流动现象时,关注的是宏观尺度而非微观尺度上发生的过程;同时假设流体是连续介质(continuum),因此其物理性质与流动性质在空间的每一点上都有定义。在这一假设下,流体流动行为可分为牛顿型与非牛顿型两类:牛顿流体的剪切应力与剪切率之间存在线性关系,分子黏度 \(\mu\)(衡量流体在应力下抵抗形变的能力)代表该线性函数的斜率;对非牛顿流体,该关系是非线性的。流体流动还可按多种维度加以分类:一维或多维、单相或多相、稳态或非稳态、真实(黏性)或理想(无黏)、可压缩或不可压缩、湍流或层流、旋转或无旋,等等。分类的目的在于简化流体流动现象的分析与建模过程。
流体流动也可按描述它们的偏微分方程进行数学分类。例如,两个自变量的二阶偏微分方程可分为双曲型、抛物型与椭圆型。方程的信息沿着两条特征线传播,根据方程分别是双曲型、抛物型还是椭圆型,特征线可以是实的且相异的、实的且重合的或复的。方程性质的这种差异要求不同的求解方法,任何用于数值求解这些方程的方法也必须识别这些性质。
本章将看到,流体流动由 Navier-Stokes 方程所支配,这些方程在一般情况下是四自变量的二阶高度非线性偏微分方程(流动通常是非稳态、三维的)。因此,上述分类严格来说并不完全适用于 Navier-Stokes 方程;然而由于它们与二自变量二阶方程分享诸多性质,同样的术语仍被用于其分类。非稳态流与超声速流动属双曲型,边界层流动属抛物型,回流流动属椭圆型。当流动在某些区域为亚声速而在另一些区域为超声速(例如收敛-扩张喷管内的流动),或在壁面附近为黏性支配而在核心区为近乎无黏时,则难以将流动归为单一类型,一般而言是混合型。在数值求解中,这种分类被翻译为:抛物型流动仅受上游影响,椭圆型流动同时受上游与下游影响,双曲型流动则可能在解中出现间断(如激波)。
3.3 守恒律的 Euler 描述与 Lagrange 描述(Eulerian and Lagrangian Description of Conservation Laws)
守恒原理指出:对于孤立系统,某些物理可测量量在局部区域内守恒。这一守恒原理或守恒律是一条无法用数学证明的公理,但可以用数学关系表达。此类定律支配多个物理量,如质量、动量与能量(即 Navier-Stokes 方程)。
涉及流体流动及相关输运现象的守恒律可按 Lagrange(物质体积,material volume,MV)或 Euler(控制体,control volume,CV)方式给出数学表述。在 Lagrange 描述(图 3.1a)中,流体被细分为若干流体微团(fluid parcels),并跟踪每一个微团在空间与时间中的运动。这些微团用一个与时间无关的位矢场 \(\mathbf{x}_0\) 进行标记,通常取该微团在某一初始时刻 \(t_0\) 的质心位置;流动由函数 \(\mathbf{x}(t;\mathbf{x}_0)\) 描述。流体微团在不同时刻所占据位置的集合形成其迹线(path line)。
在 Euler 描述(图 3.1b)中,则关注流动区域中的固定空间位置随时间的变化。因此流动变量是位置 \(\mathbf{x}\) 与时间 \(t\) 的函数,速度为 \(\mathbf{v}(t,\mathbf{x})\)。由于流体微团位矢 \(\mathbf{x}_0\) 对时间的导数代表其速度,两种描述之间由下式关联:
\[
\mathbf{v}\!\left(t,\mathbf{x}(\mathbf{x}_0,t)\right) = \frac{\partial}{\partial t}\mathbf{x}(t,\mathbf{x}_0)
\]
基于上述描述,流动流体性质的变化既可以在空间的固定点上测量(流体粒子穿越该点时,Euler 方式),也可以沿流体微团的迹线进行跟踪(Lagrange 方式)。
3.3.1 物质导数与当地导数(Substantial Versus Local Derivative)
场变量 \(\phi(t,\mathbf{x}(t))\)(可以是标量或矢量,代表密度、速度、温度等)关于空间固定位置的导数称为 Euler 导数 \(\partial \phi/\partial t\);而跟随运动流体微团所取的导数称为 Lagrange 导数、物质导数(substantial 或 material derivative),记为 \(D\phi/Dt\)。对变量 \(\phi\) 的物质导数,可通过链式法则对沿迹线的所有自变量诱导的改变求和得到:
\[
\frac{D\phi}{Dt} = \frac{\partial \phi}{\partial t} + \frac{\partial \phi}{\partial x}\frac{dx}{dt} + \frac{\partial \phi}{\partial y}\frac{dy}{dt} + \frac{\partial \phi}{\partial z}\frac{dz}{dt} = \frac{\partial \phi}{\partial t} + u\frac{\partial \phi}{\partial x} + v\frac{\partial \phi}{\partial y} + w\frac{\partial \phi}{\partial z} = \frac{\partial \phi}{\partial t} + \mathbf{v}\cdot\nabla\phi
\]
其中 \(\mathbf{v}\) 为速度矢量,\(\nabla\) 为梯度算子(已在第 2 章定义)。式 (3.2) 表明:在由 Euler 描述的速度场 \(\mathbf{v}\) 下,\(\phi\) 沿迹线从时间 \(t\)、位置 \(\mathbf{x}\) 移动到时间 \(t+dt\)、位置 \(\mathbf{x}+\mathbf{v}\,dt\) 时的总变化率等于其当地变化率与对流变化率之和。
物质导数的一个重要例子是 \(D\mathbf{v}/Dt\)——即跟随流动所取的速度变化率(加速度矢量):
\[
\frac{D\mathbf{v}}{Dt} = \frac{\partial \mathbf{v}}{\partial t} + (\mathbf{v}\cdot\nabla)\mathbf{v}
\]
本书采用 Euler 形式来描述守恒律,将焦点放在空间中某一指定区域(称为控制体)内的流动。这一选择的依据是:Euler 描述跟踪的是场(系统)而非粒子,避免了对单个粒子进行追踪的繁琐而往往不必要的工作,将注意力集中在固定点(或固定体积)上随粒子流过的变化。Lagrange 描述的一个关键缺陷在于无法控制所关注的区域,因为流体微团会被流体带向任何地方,可能离开所关注的区域;这限制了该方法的实用性,因为在大多数流体流动应用中,人们关心的是某一固定区域内的流体性质(如运动列车表面的剪切应力),而非运动物质体积的性质。尽管如此,Euler 描述将流体流动通过对流项 \(\mathbf{v}\cdot\nabla\phi\)(未知速度场与未知变量场梯度的乘积)的局部输运效应引入守恒方程。这种非线性导致了流体流动中最引人入胜、也最具挑战性的现象。
3.3.2 Reynolds 输运定理(Reynolds Transport Theorem)
前述守恒律适用于流体的运动物质体积,而非固定点或控制体。为了以 Euler 形式表达这些定律,需要知道运动物质体积上积分所对应的 Euler 等价形式——这正是 Reynolds 输运定理所提供的。
转换公式随控制体是固定的、运动的或可变形的而略有差异。设 \(B\) 为流体的任意性质(质量、动量、能量等),\(b = dB/dm\) 为 \(B\) 在流体任意小微元中的强度量(单位质量所含的 \(B\) 量)。对于图 3.1 所示的任意运动且可变形的控制体,物质体积内 \(B\) 的瞬时总变化率等于控制体 \(V\) 内 \(B\) 的瞬时总变化率加上通过控制面 \(S\) 进出控制体的 \(B\) 的净流量。设 \(\rho\) 为流体密度,\(\mathbf{n}\) 为控制体表面的外法向,\(\mathbf{v}(t,\mathbf{x})\) 为流体速度,\(\mathbf{v}_s(t,\mathbf{x})\) 为可变形的控制体表面的速度,
\[
\mathbf{v}_r(t,\mathbf{x}) = \mathbf{v}(t,\mathbf{x}) - \mathbf{v}_s(t,\mathbf{x})
\]
为流体进出控制体的相对速度,则 Reynolds 输运定理给出:
\[
\frac{dB}{dt}_{MV} = \left[\frac{d}{dt}\int_{V(t)}\rho b\, dV\right] + \int_{S(t)}\rho b\,\mathbf{v}_r\cdot\mathbf{n}\, dS
\]
对于固定的控制体,\(\mathbf{v}_s = 0\),且几何形状与时间无关,因此式 (3.4) 右端的时间导数项可由 Leibniz 法则写成:
\[
\frac{d}{dt}\int_{V}\rho b\, dV = \int_{V}\frac{\partial}{\partial t}(\rho b)\, dV
\]
从而式 (3.4) 简化为:
\[
\frac{dB}{dt}_{MV} = \int_{V}\frac{\partial}{\partial t}(\rho b)\, dV + \int_{S}\rho b\,\mathbf{v}\cdot\mathbf{n}\, dS
\]
应用散度定理将面积分转化为体积分,式 (3.6) 化为:
\[
\frac{dB}{dt}_{MV} = \int_{V}\left[\frac{\partial}{\partial t}(\rho b) + \nabla\cdot(\rho\mathbf{v}b)\right] dV
\]
式 (3.7) 的另一种等价形式,可通过对括号内第二项展开并利用物质导数得到:
\[
\frac{dB}{dt}_{MV} = \int_{V}\left[\frac{D}{Dt}(\rho b) + \rho b\,\nabla\cdot\mathbf{v}\right] dV
\]
式 (3.7) 或式 (3.8) 均可用来推导固定区域下守恒律的 Euler 形式。
3.4 质量守恒(连续方程)(Conservation of Mass / Continuity Equation)
质量守恒原理表明:在没有质量源与汇的情况下,区域将在局部守恒其质量。对图 3.3 所示的具有质量 \(m\)、密度 \(\rho\)、速度 \(\mathbf{v}\) 的物质体积,质量守恒在物质(Lagrange)坐标系下可写为:
\[
\frac{dm}{dt}_{MV} = 0
\]
对 \(B = m\),对应的强度量为 \(b = 1\)。基于式 (3.8),质量守恒在 Euler 坐标系下的等价表达式为:
\[
\int_{V}\left[\frac{D\rho}{Dt} + \rho\nabla\cdot\mathbf{v}\right] dV = 0
\]
为使式 (3.10) 中的积分对任意控制体 \(V\) 成立,被积函数必须为零,从而得到质量守恒或连续方程的微分形式:
\[
\frac{D\rho}{Dt} + \rho\nabla\cdot\mathbf{v} = 0
\]
连续方程的通量形式可利用式 (3.7) 推导得到:
\[
\int_{V}\left[\frac{\partial \rho}{\partial t} + \nabla\cdot(\rho\mathbf{v})\right] dV = 0
\]
同样为使式 (3.12) 对任意 \(V\) 成立,被积函数必须为零,从而得到质量守恒(连续方程)的通量形式:
\[
\frac{\partial \rho}{\partial t} + \nabla\cdot(\rho\mathbf{v}) = 0
\]
在不存在显著的压力或温度变化时,可以假设流动为不可压缩——即压力的变化对密度没有显著影响。这一情形在液体中几乎总是成立,对远低于声速的气体也是良好的近似(注:声波是可压缩现象)。在流体动力学中最重要的后果是:质量守恒(连续)方程不再能用来求解密度。
不可压缩条件表示 \(\rho\) 不随流动变化,可写成 \(D\rho/Dt = 0\)。利用式 (3.11),这等价于说不可压缩流动的连续方程为:
\[
\nabla\cdot\mathbf{v} = 0
\]
或其积分形式:
\[
\int_{S}\mathbf{v}\cdot\mathbf{n}\, dS = 0
\]
式 (3.15) 表明:对不可压缩流动,任意控制体上的净流量为零,即"流出量"等于"流入量"。注意,\(D\rho/Dt = 0\) 并不意味着 \(\rho\) 在空间上处处相同(尽管在许多水力学应用中确实如此),而是说 \(\rho\) 沿流线不发生变化。更准确地说,不可压缩近似意味着每个流体微元在运动过程中保持其原始密度。在工程实际中,由于盐浓度差异引起的水的密度差、由温差引起的空气密度差,都会引发显著的浮力。
3.5 线性动量守恒(Conservation of Linear Momentum)
线性动量守恒原理指出:在没有外力作用于物体时,物体保持其总动量(即质量与速度矢量的乘积)。由于动量为矢量,其在任意方向上的分量也将守恒。对物质体积,Newton 第二定律指出:该特定体积的动量只有在作用于其上的合力(可能包括表面力与体积力)下才会改变。因此对图 3.4 所示的具有质量 \(m\)、密度 \(\rho\)、速度 \(\mathbf{v}\) 的物质体积,Newton 定律在 Lagrange 坐标下可写为:
\[
\frac{d}{dt}_{MV}(m\mathbf{v}) = \left[\int_{MV}\mathbf{f}\, dV\right]
\]
其中 \(\mathbf{f}\) 为作用于该物质体积上每单位体积的外力。式 (3.16) 右端是对物质坐标下的体积分,积分域为流体瞬时占据的体积,因此:
\[
\left[\int_{V}\mathbf{f}\, dV\right]_{MV} = \int_{V}\mathbf{f}\, dV
\]
式 (3.16) 在 Euler 坐标下的等价表达式可写成两种形式,即守恒形式(conservative)与非守恒形式(non-conservative)。
注意到此处 \(b = \mathbf{v}\),利用式 (3.8) 进行推导即可得到非守恒形式:
\[
\int_{V}\left[\frac{D}{Dt}(\rho\mathbf{v}) + (\rho\nabla\cdot\mathbf{v})\mathbf{v} - \mathbf{f}\right] dV = 0
\]
为使该积分对任意控制体为零,被积函数必须为零:
\[
\frac{D}{Dt}(\rho\mathbf{v}) + (\rho\nabla\cdot\mathbf{v})\mathbf{v} = \mathbf{f}
\]
对动量项的物质导数进行展开与重组,得到非守恒形式:
\[
\rho\frac{D\mathbf{v}}{Dt} + \mathbf{v}\underbrace{\left[\frac{D\rho}{Dt} + \rho\nabla\cdot\mathbf{v}\right]}_{\text{Continuity}} = \mathbf{f}
\]
应用连续性约束并展开物质导数,动量方程的非守恒形式简化为:
\[
\rho\left[\frac{\partial \mathbf{v}}{\partial t} + (\mathbf{v}\cdot\nabla)\mathbf{v}\right] = \mathbf{f}
\]
守恒(或通量)形式通过对 Reynolds 输运定理应用式 (3.7) 给出:
\[
\int_{V}\left[\frac{\partial}{\partial t}(\rho\mathbf{v}) + \nabla\cdot[\rho\mathbf{v}\mathbf{v}] - \mathbf{f}\right] dV = 0
\]
为使积分对任意体积 \(V\) 为零,令被积函数为零,得到动量方程的守恒形式:
\[
\frac{\partial}{\partial t}(\rho\mathbf{v}) + \nabla\cdot[\rho\mathbf{v}\mathbf{v}] = \mathbf{f}
\]
其中 \(\rho\mathbf{v}\mathbf{v}\) 为第 2 章所述的并矢积(dyadic product),是张量积的一种特例,其散度为矢量。
本书中两种形式都将使用,以便更好地描述离散化概念并给出实际实现细节。后续推导采用守恒形式;非守恒形式可由守恒形式在任意步骤通过调用上述连续性约束得到。动量方程的完整形式,需要在指定作用于控制体上的外部表面力与体积力后给出。力 \(\mathbf{f}\) 分为两部分:\(\mathbf{f}_s\) 表示表面力,\(\mathbf{f}_b\) 表示体积力:
\[
\mathbf{f} = \mathbf{f}_s + \mathbf{f}_b
\]
下面分别给出这两种力的细节。
3.5.3 表面力(Surface Forces)
对图 3.4 所示的任意宏观体积微元,作用于其表面的力来自压力与黏性应力,可由总应力张量 \(\boldsymbol{\mathcal{R}}\) 表达,如图 3.5 所示。在一般情形下,每一点有九个应力分量:在每个坐标面内,一个法向分量与两个(平行于受应力面的)切向分量。因此在笛卡尔坐标下应力张量为:
\[
\boldsymbol{\mathcal{R}} = \begin{bmatrix} \mathcal{R}_{xx} & \mathcal{R}_{xy} & \mathcal{R}_{xz} \\ \mathcal{R}_{yx} & \mathcal{R}_{yy} & \mathcal{R}_{yz} \\ \mathcal{R}_{zx} & \mathcal{R}_{zy} & \mathcal{R}_{zz} \end{bmatrix}
\]
其中 \(\mathcal{R}_{ii}\) 形式的项代表法向应力,\(\mathcal{R}_{ij}\) 代表切向应力。法向应力可以是压应力(若 \(\mathcal{R}_{ii} < 0\))或拉应力(若 \(\mathcal{R}_{ii} > 0\));最重要的压法向应力通常来自压力而非黏性效应。分量 \(\mathcal{R}_{ij}\) 代表作用于 \(i\) 面、沿 \(j\) 方向的应力,\(i\) 面的方向当面的外法向沿正方向时取正。
实际中,应力张量常被分解为两项之和:
\[
\boldsymbol{\mathcal{R}} = \underbrace{\begin{bmatrix} -p & 0 & 0 \\ 0 & -p & 0 \\ 0 & 0 & -p \end{bmatrix}}_{-p\mathbf{I}} + \underbrace{\begin{bmatrix} \mathcal{R}_{xx}+p & \mathcal{R}_{xy} & \mathcal{R}_{xz} \\ \mathcal{R}_{yx} & \mathcal{R}_{yy}+p & \mathcal{R}_{yz} \\ \mathcal{R}_{zx} & \mathcal{R}_{zy} & \mathcal{R}_{zz}+p \end{bmatrix}}_{\boldsymbol{\sigma}}
\]
即
\(\boldsymbol{\mathcal{R}} = -p\mathbf{I} + \boldsymbol{\sigma}\)
,其中 \(\mathbf{I}\) 为 \(3\times 3\) 单位张量,\(p\) 为压力,\(\boldsymbol{\sigma}\) 为偏应力张量(deviatoric 或 viscous stress tensor)。压力定义为法向应力的平均值的负值:
\[
p = -\frac{1}{3}(\mathcal{R}_{xx} + \mathcal{R}_{yy} + \mathcal{R}_{zz})
\]
对面积为 \(dS\)、取向为 \(\mathbf{n}\) 的微分面元(图 3.5),作用在其上的表面力为 \((\boldsymbol{\mathcal{R}}\cdot\mathbf{n})\,dS\)。应用散度定理,作用于控制体的总表面力为:
\[
\int_{V}\mathbf{f}_s\, dV = \int_{S}\boldsymbol{\mathcal{R}}\cdot\mathbf{n}\, dS = \int_{V}\nabla\cdot\boldsymbol{\mathcal{R}}\, dV \;\;\Longrightarrow\;\; \mathbf{f}_s = [\nabla\cdot\boldsymbol{\mathcal{R}}] = -\nabla p + [\nabla\cdot\boldsymbol{\sigma}]
\]
3.5.4 体积力(Body Forces)
体积力以单位体积的力表示,可能由多种效应产生。下面给出占主导地位的几类体积力。
3.5.4.1 重力
重力场中每单位体积物质体积的重量所对应的力称为重力(图 3.6),由下式给出:
\[
\mathbf{f}_b = \rho\mathbf{g}
\]
其中 \(\mathbf{g}\) 为重力加速度矢量。
3.5.4.2 系统旋转
在旋转参考系下求解流体流动问题时,必须计及参考系刚体旋转所产生的力。这些力可视为形如式 (3.30) 的体积力:
\[
\mathbf{f}_b = \underbrace{-2\rho[\boldsymbol{\Omega}\times\mathbf{v}]}_{\text{Coriolis forces}} - \rho\underbrace{[\boldsymbol{\Omega}\times(\boldsymbol{\Omega}\times\mathbf{r})]}_{\text{Centrifugal forces}}
\]
其中 \(\boldsymbol{\Omega}\) 为旋转参考系的角速度,\(\mathbf{r}\) 为位置矢量(图 3.7)。注意,重力与离心力依赖于位置但不依赖于速度,因此可被吸收到修正压力中,从而在不出现于边界条件时不必作为独立项处理。但 Coriolis 力必须显式处理。其他力(如磁力、电力)可根据具体情形添加。由于体积力种类繁多,本书后续方程不采用任何特定类型,而保留一般性的 \(\mathbf{f}_b\) 力。
将式 (3.23) 中的外力 \(\mathbf{f}\) 用其等价表达式代入,得到一般守恒形式的动量方程:
\[
\frac{\partial}{\partial t}(\rho\mathbf{v}) + \nabla\cdot[\rho\mathbf{v}\mathbf{v}] = -\nabla p + [\nabla\cdot\boldsymbol{\sigma}] + \mathbf{f}_b
\]
3.5.5 应力张量与动量方程(Newton 流体的动量方程)(Stress Tensor and the Momentum Equation for Newtonian Fluids)
为使动量方程进一步推进,需要明确流体类型以便将应力张量 \(\boldsymbol{\sigma}\) 与流动变量联系起来。对 Newtonian 流体,应力张量是应变率的线性函数:
\[
\boldsymbol{\sigma} = \mu\left[\nabla\mathbf{v} + (\nabla\mathbf{v})^T\right] + \lambda(\nabla\cdot\mathbf{v})\mathbf{I}
\]
其中 \(\mu\) 为分子黏性系数,\(\lambda\) 为体积黏性系数(通常取 \(\lambda = -(2/3)\mu\)),上标 \(T\) 表示 \(\nabla\mathbf{v}\) 的转置,\(\mathbf{I}\) 为 \(3\times 3\) 单位张量:
\[
\mathbf{I} = \begin{bmatrix} 1 & 0 & 0 \\ 0 & 1 & 0 \\ 0 & 0 & 1 \end{bmatrix}
\]
三维笛卡尔坐标系下应力张量的展开形式为:
\[
\boldsymbol{\sigma} = \begin{bmatrix} 2\mu\frac{\partial u}{\partial x} + \lambda\nabla\cdot\mathbf{v} & \mu\left(\frac{\partial v}{\partial x}+\frac{\partial u}{\partial y}\right) & \mu\left(\frac{\partial u}{\partial z}+\frac{\partial w}{\partial x}\right) \\ \mu\left(\frac{\partial v}{\partial x}+\frac{\partial u}{\partial y}\right) & 2\mu\frac{\partial v}{\partial y} + \lambda\nabla\cdot\mathbf{v} & \mu\left(\frac{\partial w}{\partial y}+\frac{\partial v}{\partial z}\right) \\ \mu\left(\frac{\partial u}{\partial z}+\frac{\partial w}{\partial x}\right) & \mu\left(\frac{\partial w}{\partial y}+\frac{\partial v}{\partial z}\right) & 2\mu\frac{\partial w}{\partial z} + \lambda\nabla\cdot\mathbf{v} \end{bmatrix}
\]
应力张量的散度为一个矢量:
\[
[\nabla\cdot\boldsymbol{\sigma}] = \nabla\cdot\mu\left[\nabla\mathbf{v} + (\nabla\mathbf{v})^T\right] + \nabla(\lambda\nabla\cdot\mathbf{v})
\]
其三个分量形式为:
\[
[\nabla\cdot\boldsymbol{\sigma}]_x = \frac{\partial}{\partial x}\!\left(2\mu\frac{\partial u}{\partial x} + \lambda\nabla\cdot\mathbf{v}\right) + \frac{\partial}{\partial y}\!\left[\mu\left(\frac{\partial v}{\partial x}+\frac{\partial u}{\partial y}\right)\right] + \frac{\partial}{\partial z}\!\left[\mu\left(\frac{\partial u}{\partial z}+\frac{\partial w}{\partial x}\right)\right]
\]
其余两分量类推。
代入式 (3.31),得到 Newtonian 流体动量方程的最终守恒形式:
\[
\frac{\partial}{\partial t}(\rho\mathbf{v}) + \nabla\cdot[\rho\mathbf{v}\mathbf{v}] = -\nabla p + \nabla\cdot\!\left\{\mu\left[\nabla\mathbf{v} + (\nabla\mathbf{v})^T\right] + \lambda(\nabla\cdot\mathbf{v})\right\} + \mathbf{f}_b
\]
为后续引用方便,将动量方程展开为:
\[
\frac{\partial}{\partial t}(\rho\mathbf{v}) + \nabla\cdot[\rho\mathbf{v}\mathbf{v}] = \nabla\cdot[\mu\nabla\mathbf{v}] - \nabla p + \nabla\cdot\mu(\nabla\mathbf{v})^T + \nabla(\lambda\nabla\cdot\mathbf{v}) + \mathbf{f}_b
\]
并改写为:
\[
\frac{\partial}{\partial t}(\rho\mathbf{v}) + \nabla\cdot[\rho\mathbf{v}\mathbf{v}] = \nabla\cdot[\mu\nabla\mathbf{v}] - \nabla p + \mathbf{Q}_v
\]
对不可压缩流动,速度矢量的散度为零,即 \(\nabla\cdot\mathbf{v} = 0\),动量方程简化为:
\[
\frac{\partial}{\partial t}(\rho\mathbf{v}) + \nabla\cdot[\rho\mathbf{v}\mathbf{v}] = -\nabla p + \nabla\cdot\!\left\{\mu\left[\nabla\mathbf{v} + (\nabla\mathbf{v})^T\right]\right\} + \mathbf{f}_b
\]
若黏度为常数,动量方程可进一步简化。取矢量方程的第一分量并设 \(\mu\) 为常数,可写出:
\[
\mu\frac{\partial}{\partial x}\!\left(2\frac{\partial u}{\partial x}\right) + \mu\frac{\partial}{\partial y}\!\left(\frac{\partial v}{\partial x}+\frac{\partial u}{\partial y}\right) + \mu\frac{\partial}{\partial z}\!\left(\frac{\partial u}{\partial z}+\frac{\partial w}{\partial x}\right) = \mu\!\left(\frac{\partial^2 u}{\partial x^2} + \frac{\partial^2 u}{\partial x^2} + \frac{\partial^2 u}{\partial y^2} + \frac{\partial^2 u}{\partial y\partial x} + \frac{\partial^2 u}{\partial z^2} + \frac{\partial^2 u}{\partial z\partial x}\right) - \mu\!\left(\frac{\partial^2 u}{\partial x^2} + \frac{\partial^2 u}{\partial y^2} + \frac{\partial^2 u}{\partial z^2}\right) + \mu\frac{\partial}{\partial x}\!\left(\frac{\partial u}{\partial x} + \frac{\partial v}{\partial y} + \frac{\partial w}{\partial z}\right)
\]
代入式 (3.37) 简化后得:
\[
\frac{\partial}{\partial t}(\rho\mathbf{v}) + \nabla\cdot[\rho\mathbf{v}\mathbf{v}] = -\nabla p + \mu\nabla^2\mathbf{v} + \mathbf{f}_b
\]
对无黏流动(黏度为零),不可压缩与可压缩无黏流动的动量方程变为:
\[
\frac{\partial}{\partial t}(\rho\mathbf{v}) + \nabla\cdot[\rho\mathbf{v}\mathbf{v}] = -\nabla p + \mathbf{f}_b
\]
3.6 能量守恒(Conservation of Energy)
能量守恒受热力学第一定律支配,该定律指出:能量在过程中既不能被创造也不能被消灭,只能从一种形式(机械、动能、化学等)转化为另一种形式。因此孤立系统中所有形式能量之总和保持不变。
对图 3.8 所示的物质体积,其质量为 \(m\)、密度为 \(\rho\)、以速度 \(\mathbf{v}\) 运动。定义物质体积在时刻 \(t\) 的总能量 \(E\) 为其内能与动能之和:
\[
E = m\!\left(\hat{u} + \tfrac{1}{2}\mathbf{v}\cdot\mathbf{v}\right)
\]
其中 \(\hat{u}\) 为流体比内能(单位质量的内能)。将热力学第一定律应用于物质体积,给出物质体积总能量的变化率等于通过其边界的加热率与做功率的差。数学上写为:
\[
\frac{dE}{dt}_{MV} = \dot{Q} - \dot{W}
\]
本书采用如下符号约定:加入到物质体积的热量与物质体积对外所做的功均为正。为对物质体积应用 Reynolds 输运定理,取 \(B = E\)、\(b = e\)(单位质量的总能量),因此:
\[
B = E \;\Longrightarrow\; b = \frac{dE}{dm} = \hat{u} + \tfrac{1}{2}\mathbf{v}\cdot\mathbf{v} = e
\]
物质微元所获得的总加热率 \(\dot{Q}\) 由两部分组成:通过表面传入的加热率 \(\dot{Q}_S\) 与在物质体积内产生/消耗(如化学反应)的加热率 \(\dot{Q}_V\)。物质体积所做的总功率 \(\dot{W}\) 由表面力所做功率 \(\dot{W}_S\) 与体积力所做功率 \(\dot{W}_b\) 共同贡献。因此第一定律可写为:
\[
\frac{dE}{dt}_{MV} = \dot{Q}_V + \dot{Q}_S - \dot{W}_b - \dot{W}_S
\]
按定义,功是力作用一段距离的结果,功率是单位时间所做的功。因此体积力与表面力所做功率可表示为:
\[
\dot{W}_b = -\int_{V}(\mathbf{f}_b\cdot\mathbf{v})\, dV, \quad \dot{W}_S = -\int_{S}(\mathbf{f}_S\cdot\mathbf{v})\, dS
\]
利用式 (3.26)–(3.28) 中 \(\mathbf{f}_S\) 的等价表达式对表面力所做功率进行展开:
\[
\dot{W}_S = -\int_{S}[(\boldsymbol{\mathcal{R}}\cdot\mathbf{v})\cdot\mathbf{n}]\, dS = -\int_{V}\nabla\cdot[(\boldsymbol{\mathcal{R}}\cdot\mathbf{v})]\, dV = -\int_{V}\nabla\cdot[(-p\mathbf{I}+\boldsymbol{\sigma})\cdot\mathbf{v}]\, dV
\]
整理后 \(\dot{W}_S\) 可重写为:
\[
\dot{W}_S = -\int_{V}\!\left[-\nabla\cdot(p\mathbf{v}) + \nabla\cdot(\boldsymbol{\sigma}\cdot\mathbf{v})\right] dV
\]
设 \(\dot{q}_V\) 为物质体积内单位体积的热源或热汇的加热率,\(\dot{q}_S\) 为通过物质微元表面单位面积的热流率,则 \(\dot{Q}_V\) 与 \(\dot{Q}_S\) 可写为:
\[
\dot{Q}_V = \int_{V}\dot{q}_V\, dV, \quad \dot{Q}_S = -\int_{S}\dot{\mathbf{q}}_s\cdot\mathbf{n}\, dS = -\int_{V}\nabla\cdot\dot{\mathbf{q}}_s\, dV
\]
应用 Reynolds 输运定理并将做功与加热项代入其等价表达式,式 (3.46) 变为:
\[
\frac{dE}{dt}_{MV} = \int_{V}\left[\frac{\partial}{\partial t}(\rho e) + \nabla\cdot(\rho\mathbf{v}e)\right] dV = -\int_{V}\nabla\cdot\dot{\mathbf{q}}_s\, dV - \int_{V}\!\left[-\nabla\cdot(p\mathbf{v}) + \nabla\cdot(\boldsymbol{\sigma}\cdot\mathbf{v})\right] dV + \int_{V}(\mathbf{f}_b\cdot\mathbf{v})\, dV + \int_{V}\dot{q}_V\, dV
\]
合并各项后,上式化为:
\[
\int_{V}\left[\frac{\partial}{\partial t}(\rho e) + \nabla\cdot(\rho\mathbf{v}e) + \nabla\cdot\dot{\mathbf{q}}_s + \nabla\cdot(p\mathbf{v}) - \nabla\cdot(\boldsymbol{\sigma}\cdot\mathbf{v}) - \mathbf{f}_b\cdot\mathbf{v} - \dot{q}_V\right] dV = 0
\]
为使式 (3.52) 中的体积分对任意控制体成立,被积函数必须为零:
\[
\frac{\partial}{\partial t}(\rho e) + \nabla\cdot(\rho\mathbf{v}e) = -\nabla\cdot\dot{\mathbf{q}}_s - \nabla\cdot(p\mathbf{v}) + \nabla\cdot(\boldsymbol{\sigma}\cdot\mathbf{v}) + \mathbf{f}_b\cdot\mathbf{v} + \dot{q}_V
\]
这就是用比总能表示的能量守恒方程(能量方程)。能量方程还可以用比内能、比静焓(或简称比焓)、比总焓表示,并在特定条件下用温度表示。
3.6.1 用比内能表示的能量守恒(Conservation of Energy in Terms of Specific Internal Energy)
将式 (3.53) 用比内能改写,需对动量方程 (3.23) 取与速度矢量的点积:
\[
\left[\frac{\partial}{\partial t}(\rho\mathbf{v}) + \nabla\cdot[\rho\mathbf{v}\mathbf{v}]\right]\cdot\mathbf{v} = \mathbf{f}\cdot\mathbf{v}
\]
经若干整理后,式 (3.54) 变为:
\[
\frac{\partial}{\partial t}(\rho\mathbf{v}\cdot\mathbf{v}) - \mathbf{v}\cdot\frac{\partial}{\partial t}(\rho\mathbf{v}) + \nabla\cdot[\rho(\mathbf{v}\cdot\mathbf{v})\mathbf{v}] - \rho\mathbf{v}\cdot[(\mathbf{v}\cdot\nabla)\mathbf{v}] = \mathbf{f}\cdot\mathbf{v}
\]
重新整理与合并各项得到:
\[
\frac{\partial}{\partial t}(\rho\mathbf{v}\cdot\mathbf{v}) + \nabla\cdot[\rho(\mathbf{v}\cdot\mathbf{v})\mathbf{v}] - \mathbf{v}\cdot\underbrace{\left[\frac{\partial}{\partial t}(\rho\mathbf{v}) + (\mathbf{v}\cdot\nabla)\mathbf{v}\right]}_{=\mathbf{f}\,\text{by Eq. (3.21)}} = \mathbf{f}\cdot\mathbf{v}
\]
注意到左端第三项为 \(\mathbf{v}\cdot\mathbf{f}\),将 \(\mathbf{f}\) 用其等价表达式代入,可得到关于流动动能的方程:
\[
\frac{\partial}{\partial t}\!\left(\tfrac{1}{2}\rho\mathbf{v}\cdot\mathbf{v}\right) + \nabla\cdot\left(\tfrac{1}{2}\rho(\mathbf{v}\cdot\mathbf{v})\mathbf{v}\right) = -\mathbf{v}\cdot\nabla p + \mathbf{v}\cdot[\nabla\cdot\boldsymbol{\sigma}] + \mathbf{f}_b\cdot\mathbf{v}
\]
该方程可改写为:
\[
\frac{\partial}{\partial t}\!\left(\tfrac{1}{2}\rho\mathbf{v}\cdot\mathbf{v}\right) + \nabla\cdot\left(\tfrac{1}{2}\rho(\mathbf{v}\cdot\mathbf{v})\mathbf{v}\right) = -\nabla\cdot(p\mathbf{v}) + p\nabla\cdot\mathbf{v} + \nabla\cdot(\boldsymbol{\sigma}\cdot\mathbf{v}) - (\boldsymbol{\sigma}:\nabla\mathbf{v}) + \mathbf{f}_b\cdot\mathbf{v}
\]
将式 (3.58) 从式 (3.53) 中减去,得到以比内能为主变量的能量方程:
\[
\frac{\partial}{\partial t}(\rho\hat{u}) + \nabla\cdot(\rho\mathbf{v}\hat{u}) = -\nabla\cdot\dot{\mathbf{q}}_s - p\nabla\cdot\mathbf{v} + (\boldsymbol{\sigma}:\nabla\mathbf{v}) + \dot{q}_V
\]
3.6.2 用比焓表示的能量守恒(Conservation of Energy in Terms of Specific Enthalpy)
将能量方程改写为用比焓表示的形式十分直接,依据比内能与比焓的关系:
\[
\hat{u} = \hat{h} - \frac{p}{\rho}
\]
将 \(\hat{u}\) 替换为 \(\hat{h} - p/\rho\) 代入式 (3.59),经代数整理,得到以比焓表示的能量方程:
\[
\frac{\partial}{\partial t}(\rho\hat{h}) + \nabla\cdot(\rho\mathbf{v}\hat{h}) = -\nabla\cdot\dot{\mathbf{q}}_s + \frac{Dp}{Dt} + (\boldsymbol{\sigma}:\nabla\mathbf{v}) + \dot{q}_V
\]
3.6.3 用比总焓表示的能量守恒(Conservation of Energy in Terms of Specific Total Enthalpy)
将比总能 \(e\) 用比总焓 \(\hat{h}_0\) 表达,可推得以比总焓表示的能量方程。首先:
\[
e = \hat{u} + \tfrac{1}{2}\mathbf{v}\cdot\mathbf{v} = \hat{h} - \frac{p}{\rho} + \tfrac{1}{2}\mathbf{v}\cdot\mathbf{v} = \hat{h}_0 - \frac{p}{\rho}
\]
将 \(\hat{h}_0 - p/\rho\) 替换 \(e\) 代入式 (3.53),整理后得到以比总焓表示的能量方程:
\[
\frac{\partial}{\partial t}(\rho\hat{h}_0) + \nabla\cdot(\rho\mathbf{v}\hat{h}_0) = -\nabla\cdot\dot{\mathbf{q}}_s + \frac{\partial p}{\partial t} + \nabla\cdot(\boldsymbol{\sigma}\cdot\mathbf{v}) + \mathbf{f}_b\cdot\mathbf{v} + \dot{q}_V
\]
至此所给出的所有形式的能量方程都是一般性的,适用于 Newtonian 与非 Newtonian 流体的流动。唯一的限制是它们适用于固定控制体。
3.6.4 用温度表示的能量守恒(Conservation of Energy in Terms of Temperature)
要将能量方程写成以温度为主变量的形式,需要附加若干约束。假设 \(\hat{h}\) 为 \(p\) 与 \(T\) 的函数,且流体为 Newtonian;因此后续推导仅适用于 Newtonian 流体。若 \(\hat{h} = \hat{h}(p,T)\),则 \(d\hat{h}\) 可写为:
\[
d\hat{h} = \left.\frac{\partial \hat{h}}{\partial T}\right|_{p} dT + \left.\frac{\partial \hat{h}}{\partial p}\right|_{T} dp
\]
利用下列平衡态热力学关系:
\[
\left.\frac{\partial \hat{h}}{\partial p}\right|_{T} = \hat{V} - T\left.\frac{\partial \hat{V}}{\partial T}\right|_{p}
\]
其中 \(\hat{V}\) 为比容,\(d\hat{h}\) 的表达式可修正为:
\[
d\hat{h} = c_p\, dT + \left[\hat{V} - T\left.\frac{\partial \hat{V}}{\partial T}\right|_{p}\right] dp
\]
将 \(d\hat{h}\) 的表达式 (3.66) 代入式 (3.61) 左端,可将左端用 \(T\) 表达为:
\[
\frac{\partial}{\partial t}(\rho\hat{h}) + \nabla\cdot(\rho\mathbf{v}\hat{h}) = \rho\frac{D\hat{h}}{Dt} = \rho c_p\frac{DT}{Dt} + \rho\!\left[\hat{V} - T\left.\frac{\partial \hat{V}}{\partial T}\right|_{p}\right]\frac{DP}{Dt} = \rho c_p\frac{DT}{Dt} + \rho\!\left[\frac{1}{\rho} - T\frac{\partial(1/\rho)}{\partial T}\bigg|_{p}\right]\frac{DP}{Dt} = \rho c_p\frac{DT}{Dt} + \left[1 + \frac{\partial(\ln\rho)}{\partial(\ln T)}\bigg|_{p}\right]\frac{DP}{Dt}
\]
将式 (3.67) 代入式 (3.61),得到以 \(T\) 为主变量的能量方程:
\[
\rho c_p\frac{DT}{Dt} = -\nabla\cdot\dot{\mathbf{q}}_s - \left[\frac{\partial(\ln\rho)}{\partial(\ln T)}\bigg|_{p}\right]\frac{Dp}{Dt} + (\boldsymbol{\sigma}:\nabla\mathbf{v}) + \dot{q}_V
\]
等价形式为:
\[
c_p\frac{\partial}{\partial t}(\rho T) + c_p\nabla\cdot(\rho\mathbf{v}T) = -\nabla\cdot\dot{\mathbf{q}}_s - \left[\frac{\partial(\ln\rho)}{\partial(\ln T)}\bigg|_{p}\right]\frac{Dp}{Dt} + (\boldsymbol{\sigma}:\nabla\mathbf{v}) + \dot{q}_V
\]
出现在所有能量方程形式中的热通量 \(\dot{\mathbf{q}}_S\) 代表通过扩散进行的热量传递——这是发生在分子层面的现象,由 Fourier 定律支配:
\[
\dot{\mathbf{q}}_s = -[k\nabla T]
\]
其中 \(k\) 为物质的热导率。该方程表明热量沿温度梯度方向流动,并假设材料对热传导没有优先方向,即介质是各向同性的。然而某些固体是各向异性的,式 (3.70) 需替换为:
\[
\dot{\mathbf{q}}_s = -[\boldsymbol{\kappa}\cdot\nabla T]
\]
其中 \(\boldsymbol{\kappa}\) 为二阶对称张量,称为热导率张量。因此在各向异性介质中,热通量方向并不与温度梯度方向一致。后续推导假设介质为各向同性,可应用式 (3.70)。将 Fourier 定律代入后,能量方程 (3.69) 变为:
\[
c_p\frac{\partial}{\partial t}(\rho T) + c_p\nabla\cdot(\rho\mathbf{v}T) = \nabla\cdot[k\nabla T] - \left[\frac{\partial(\ln\rho)}{\partial(\ln T)}\bigg|_{p}\right]\frac{Dp}{Dt} + (\boldsymbol{\sigma}:\nabla\mathbf{v}) + \dot{q}_V
\]
\((\boldsymbol{\sigma}:\nabla\mathbf{v})\) 在三维笛卡尔坐标下用流动变量表示为:
\[
(\boldsymbol{\sigma}:\nabla\mathbf{v}) = \lambda\!\left(\frac{\partial u}{\partial x} + \frac{\partial v}{\partial y} + \frac{\partial w}{\partial z}\right)^2 + 2\mu\!\left[\left(\frac{\partial u}{\partial x}\right)^2 + \left(\frac{\partial v}{\partial y}\right)^2 + \left(\frac{\partial w}{\partial z}\right)^2\right] + \mu\!\left[\left(\frac{\partial u}{\partial y}+\frac{\partial v}{\partial x}\right)^2 + \left(\frac{\partial u}{\partial z}+\frac{\partial w}{\partial x}\right)^2 + \left(\frac{\partial v}{\partial z}+\frac{\partial w}{\partial y}\right)^2\right]
\]
定义 \(W\) 与 \(U\) 如下:
\[
W = \frac{\partial u}{\partial x} + \frac{\partial v}{\partial y} + \frac{\partial w}{\partial z}
\]
\[
U = 2\!\left[\left(\frac{\partial u}{\partial x}\right)^2 + \left(\frac{\partial v}{\partial y}\right)^2 + \left(\frac{\partial w}{\partial z}\right)^2\right] + \left(\frac{\partial u}{\partial y}+\frac{\partial v}{\partial x}\right)^2 + \left(\frac{\partial u}{\partial z}+\frac{\partial w}{\partial x}\right)^2 + \left(\frac{\partial v}{\partial z}+\frac{\partial w}{\partial y}\right)^2
\]
以温度表示的能量方程化为:
\[
c_p\frac{\partial}{\partial t}(\rho T) + c_p\nabla\cdot(\rho\mathbf{v}T) = \nabla\cdot[k\nabla T] - \left[\frac{\partial(\ln\rho)}{\partial(\ln T)}\bigg|_{p}\right]\frac{Dp}{Dt} + \lambda W + \mu U + \dot{q}_V
\]
为后续引用方便,将能量方程展开为:
\[
\frac{\partial}{\partial t}(\rho c_p T) + \nabla\cdot(\rho c_p \mathbf{v}T) = \nabla\cdot[k\nabla T] + \rho T\frac{Dc_p}{Dt} - \left[\frac{\partial(\ln\rho)}{\partial(\ln T)}\bigg|_{p}\right]\frac{Dp}{Dt} + \lambda W + \mu U + \dot{q}_V = \nabla\cdot[k\nabla T] + \mathbf{Q}_T
\]
并改写为:
\[
\frac{\partial}{\partial t}(\rho c_p T) + \nabla\cdot(\rho c_p \mathbf{v}T) = \nabla\cdot[k\nabla T] + \mathbf{Q}_T
\]
能量方程一般不以完整形式求解;针对具体物理情形可建立若干简化版本。耗散项 \(\Phi\) 在速度梯度不大或亚声速情形下通常可忽略。对不可压缩流体,连续方程意味着 \(W = 0\),且密度为常数,从而
\[
\partial(\ln\rho)/\partial(\ln T) = 0
\]
。因此不可压缩流体流动的能量方程 (3.77) 简化为:
\[
\frac{\partial}{\partial t}(\rho c_p T) + \nabla\cdot(\rho c_p \mathbf{v}T) = \nabla\cdot[k\nabla T] + \dot{q}_V + \rho T\frac{Dc_p}{Dt} = \nabla\cdot[k\nabla T] + \mathbf{Q}_T
\]
式 (3.79) 也适用于定压系统下的流体流动。对固体情形,密度为常数、速度为零,若温度变化不大则热导率可视为常数,能量方程变为:
\[
\rho c_p\frac{\partial T}{\partial t} = k\nabla^2 T + \dot{q}_V
\]
对理想气体,
\[
\partial(\ln\rho)/\partial(\ln T) = -1
\]
,理想气体可压缩流动的能量方程简化为:
\[
c_p\frac{\partial}{\partial t}(\rho T) + c_p\nabla\cdot(\rho\mathbf{v}T) = \nabla\cdot[k\nabla T] + \frac{Dp}{Dt} + \lambda W + \mu U + \dot{q}_V
\]
若黏性可忽略(流动为无黏),则式 (3.81) 进一步简化为:
\[
c_p\frac{\partial}{\partial t}(\rho T) + c_p\nabla\cdot(\rho\mathbf{v}T) = \nabla\cdot[k\nabla T] + \frac{Dp}{Dt} + \dot{q}_V
\]
3.7 通用守恒方程(General Conservation Equation)
质量、动量与能量的守恒方程可由"比量"或称强度量(intensive property)——即以单位质量表示的量——表达。例如动量方程用"单位质量的动量"——即速度——来表达线性动量守恒。同样的守恒方程形式可应用于任何强度量 \(\phi\),例如溶液中盐的浓度或某化学组分的质量分数。控制体内 \(\phi\) 随时间的变化可以写成形如下式的平衡方程:
\[
\text{Term I} = \text{Term II} + \text{Term III}
\]
对图 3.9 所示的固定控制体,利用 Reynolds 输运定理,物质体积内 \(\phi\) 随时间的变化可写为:
\[
\text{Term I} = \frac{d}{dt}\int_{MV}(\rho\phi)\, dV = \int_{V}\left[\frac{\partial}{\partial t}(\rho\phi) + \nabla\cdot(\rho\mathbf{v}\phi)\right] dV
\]
其中 \(\rho\) 为流体密度,\(V\) 为控制体体积,\(S\) 为其表面积。项 \(\rho\mathbf{v}\phi\) 表示 \(\phi\) 由流场的输运,称为对流通量(convective flux):
\[
\mathbf{J}_{\text{convection}} = \rho\mathbf{v}\phi
\]
第二项表示由跨越控制体表面的物理过程引起的 \(\phi\) 变化。本书所关注的物理过程中,导致 \(\phi\) 流入/流出的机制是扩散——由分子碰撞产生,记为 \(\mathbf{J}_{\text{diffusion}}\)。设 \(\phi\) 的扩散系数为 \(\Gamma_\phi\),则扩散通量可写为:
\[
\mathbf{J}_{\text{diffusion}} = -\Gamma_\phi\nabla\phi
\]
从而 Term II 变为:
\[
\text{Term II} = -\int_{S}\mathbf{J}_{\text{diffusion}}\cdot\mathbf{n}\, dS = -\int_{V}\nabla\cdot\mathbf{J}_{\text{diffusion}}\, dV = \int_{V}\nabla\cdot(\Gamma_\phi\nabla\phi)\, dV
\]
其中 \(\mathbf{n}\) 为控制面外法向单位矢量,负号源于本书所采用的符号约定(即向内通量为正)。Term III 可写为:
\[
\text{Term III} = \int_{V}Q_\phi\, dV
\]
其中 \(Q_\phi\) 为单位体积内 \(\phi\) 的生成/消耗率,也称为源项。因此守恒方程可表达为:
\[
\int_{V}\left[\frac{\partial}{\partial t}(\rho\phi) + \nabla\cdot(\rho\phi\mathbf{v})\right] dV = \int_{V}\nabla\cdot(\Gamma_\phi\nabla\phi)\, dV + \int_{V}Q_\phi\, dV
\]
整理为:
\[
\int_{V}\left[\frac{\partial}{\partial t}(\rho\phi) + \nabla\cdot(\rho\mathbf{v}\phi) - \nabla\cdot(\Gamma_\phi\nabla\phi) - Q_\phi\right] dV = 0
\]
为使积分对任意控制体为零,被积函数必须为零,得到守恒方程的微分形式:
\[
\frac{\partial}{\partial t}(\rho\phi) + \nabla\cdot(\rho\mathbf{v}\phi) - \nabla\cdot(\Gamma_\phi\nabla\phi) - Q_\phi = 0
\]
为后续引用方便,方程可改写为:
\[
\frac{\partial}{\partial t}(\rho\phi) + \nabla\cdot\mathbf{J}_\phi - Q_\phi = 0
\]
其中总通量 \(\mathbf{J}_\phi\) 为对流通量与扩散通量之和:
\[
\mathbf{J}_\phi = \mathbf{J}_{\phi,C} + \mathbf{J}_{\phi,D} = \rho\mathbf{v}\phi - \Gamma_\phi\nabla\phi
\]
\(\phi\) 输运的通用守恒方程 (3.90) 的最终形式为:
\[
\underbrace{\frac{\partial}{\partial t}(\rho\phi)}_{\text{unsteady term}} + \underbrace{\nabla\cdot(\rho\mathbf{v}\phi)}_{\text{convection term}} = \underbrace{\nabla\cdot(\Gamma_\phi\nabla\phi)}_{\text{diffusion term}} + \underbrace{Q_\phi}_{\text{source term}}
\]
将式 (3.93) 与前面推导出的各守恒方程对照可知:通过为 \(\phi\)、\(\Gamma_\phi\) 与 \(Q_\phi\) 赋予适当的值,式 (3.93) 便可代表任意守恒方程。这是一个非常重要的观察——它使得后续章节中数值方法的发展得以集中在通用方程 (3.93) 上,而非各独立的守恒方程。
3.8 无量纲化过程(Non-dimensionalization Procedure)
守恒定律对应的微分方程很少直接用有因次变量求解。常见做法是选取合适的特征尺度,将这些方程写成无量纲形式。引入无量纲变量有若干优点:减少问题所涉及的参数数量;揭示守恒方程各项的相对量级,并据此判断哪些项可忽略;简化待解方程并保留量级相近的项,从而提高数值精度;此外所得解可适用于所有动力学相似的问题。
将一个有因次变量转换为无量纲量,是将其除以一个具有相同量纲的量(由一个或多个物理性质组成)。例如空间坐标可除以特征长度,速度可除以特征速度或组合量
\(\Gamma_{\text{ref}}/(\rho_{\text{ref}} L_{\text{ref}})\)
(整体量纲与速度 m/s 相同),压力通常除以参考动压 \(\rho_{\text{ref}} v_{\text{ref}}^2\),时间可除以 \(L_{\text{ref}}/v_{\text{ref}}\),等等。理解无量纲化过程的最好方式是举例。考虑一黏度与热导率为常数的不可压缩黏性流动,且体积力仅沿 \(y\) 方向作用(即重力加速度 \(\mathbf{g} = (0, -g, 0)\))。在无内热源的三维笛卡尔坐标系下,质量、动量与能量守恒方程为:
\[
\frac{\partial u}{\partial x} + \frac{\partial v}{\partial y} + \frac{\partial w}{\partial z} = 0
\]
\[
\frac{\partial}{\partial t}(\rho u) + \frac{\partial}{\partial x}(\rho u u) + \frac{\partial}{\partial y}(\rho v u) + \frac{\partial}{\partial z}(\rho w u) = -\frac{\partial p}{\partial x} + \mu\!\left(\frac{\partial^2 u}{\partial x^2} + \frac{\partial^2 u}{\partial y^2} + \frac{\partial^2 u}{\partial z^2}\right)
\]
\[
\frac{\partial}{\partial t}(\rho v) + \frac{\partial}{\partial x}(\rho u v) + \frac{\partial}{\partial y}(\rho v v) + \frac{\partial}{\partial z}(\rho w v) = -\frac{\partial p}{\partial y} + \mu\!\left(\frac{\partial^2 v}{\partial x^2} + \frac{\partial^2 v}{\partial y^2} + \frac{\partial^2 v}{\partial z^2}\right) - \rho g
\]
\[
\frac{\partial}{\partial t}(\rho w) + \frac{\partial}{\partial x}(\rho u w) + \frac{\partial}{\partial y}(\rho v w) + \frac{\partial}{\partial z}(\rho w w) = -\frac{\partial p}{\partial z} + \mu\!\left(\frac{\partial^2 w}{\partial x^2} + \frac{\partial^2 w}{\partial y^2} + \frac{\partial^2 w}{\partial z^2}\right)
\]
\[
c_p\frac{\partial}{\partial t}(\rho T) + \frac{\partial}{\partial x}(\rho u T) + \frac{\partial}{\partial y}(\rho v T) + \frac{\partial}{\partial z}(\rho w T) = k\!\left(\frac{\partial^2 T}{\partial x^2} + \frac{\partial^2 T}{\partial y^2} + \frac{\partial^2 T}{\partial z^2}\right)
\]
当体积力相对于其他力可忽略时,\(\rho g\) 项可置零并从方程中移除。此时流场独立于温度场,可先独立求解速度场,再求解温度场。然而,要使流动存在,必须有外力将流体推过计算域——即流体应具有入口速度。该速度是写成无量纲形式时的重要参数(特征速度),其传热机制称为受迫对流(forced convection)。
另一方面,当流场由域内温差自然引发时,体积力不可忽略。此时温度变化引起密度变化(如前所述),产生浮力驱动流动。此情形下的传热称为自然对流(natural convection)。由于流动由自然引发,不存在明显的特征速度,也不能作为无量纲数的一部分——其尺度未知。为将速度写成无量纲形式,必须使用与速度具有相同量纲的物性组合。后续讨论针对自然对流问题。
设温差 \(\Delta T = T - T_1\)(其中 \(T_1\) 为域内最低与最高温度之间的参考温度,通常取平均值)足够小,使得 \(\Delta T^2\) 及更高阶项可忽略,则任意温度 \(T\) 处的密度可由其在参考温度 \(T_1\) 处取值通过 Taylor 展开的一阶截断表达:
\[
\rho = \rho\big|_{T=T_1} + \left.\frac{d\rho}{dT}\right|_{T=T_1}(T - T_1)
\]
其中略去了 \(\Delta T^2\) 及更高阶项。引入体膨胀系数 \(\beta\):
\[
\beta = -\frac{1}{\rho}\frac{\partial \rho}{\partial T}\bigg|_{p}
\]
密度(或状态方程)变为:
\[
\rho = \rho_1[1 - \beta(T - T_1)]
\]
该式在文献中称为 Boussinesq 近似。仅在体积力项中使用该密度表达式,并以常数 \(\rho\) 表示密度以简化记号,\(y\) 方向动量方程变换为:
\[
\frac{\partial}{\partial t}(\rho v) + \frac{\partial}{\partial x}(\rho u v) + \frac{\partial}{\partial y}(\rho v v) + \frac{\partial}{\partial z}(\rho w v) = -\frac{\partial}{\partial y}(p + \rho g y) + \mu\!\left(\frac{\partial^2 v}{\partial x^2} + \frac{\partial^2 v}{\partial y^2} + \frac{\partial^2 v}{\partial z^2}\right) + \rho g\beta(T - T_1)
\]
这清楚表明:在求解自然对流问题时,动量方程与能量方程相互耦合,需要联立求解两者。
通过定义以下无量纲量,可得到守恒方程的无量纲形式:
\[
\hat{x} = \frac{x}{L}, \quad \hat{y} = \frac{y}{L}, \quad \hat{z} = \frac{z}{L}, \quad \hat{u} = \frac{u}{\Gamma/(\rho L)}, \quad \hat{v} = \frac{v}{\Gamma/(\rho L)}, \quad \hat{w} = \frac{w}{\Gamma/(\rho L)}
\]
\[
\hat{t} = \frac{t}{\rho L^2/\Gamma}, \quad \hat{p} = \frac{p + \rho g y}{\Gamma^2/(\rho L^2)}, \quad \hat{T} = \frac{T - T_1}{T_{\max} - T_1}
\]
其中 \(L\) 为特征长度,\(\Gamma\) 为流体的动力黏度,\(T_{\max}\) 为域内最高温度,上标 ^ 标记无量纲量。守恒方程中的各项可按上述新变量改写。下面通过从每一类中取一个典型项来说明过程。
连续方程中的典型项:
\[
\frac{\partial u}{\partial x} = \frac{\partial[\Gamma\hat{u}/(\rho L)]}{\partial(L\hat{x})} = \frac{\Gamma/(\rho L)}{L}\frac{\partial \hat{u}}{\partial \hat{x}} = \frac{\Gamma}{\rho L^2}\frac{\partial \hat{u}}{\partial \hat{x}}
\]
动量方程中的典型项:
\[
\frac{\partial}{\partial t}(\rho u) = \frac{\partial(\Gamma\hat{u}/L)}{\partial(\rho L^2\hat{t}/\Gamma)} = \frac{\Gamma/L}{\rho L^2/\Gamma}\frac{\partial \hat{u}}{\partial \hat{t}} = \frac{\Gamma^2}{\rho L^3}\frac{\partial \hat{u}}{\partial \hat{t}}
\]
\[
\frac{\partial}{\partial t}(\rho u u) = \frac{\partial[\Gamma^2/(\rho L^2)\hat{u}\hat{u}]}{\partial(L\hat{x})} = \frac{\Gamma^2/(\rho L^2)}{L}\frac{\partial}{\partial \hat{x}}(\hat{u}\hat{u}) = \frac{\Gamma^2}{\rho L^3}\frac{\partial}{\partial \hat{x}}(\hat{u}\hat{u})
\]
\[
\hat{p} = \frac{p + \rho g y}{\Gamma^2/(\rho L^2)} \;\Longrightarrow\; \frac{\partial}{\partial \hat{x}} = \frac{\rho L^3}{\Gamma^2}\frac{\partial p}{\partial x} \;\Longrightarrow\; \frac{\partial p}{\partial x} = \frac{\Gamma^2}{\rho L^3}\frac{\partial \hat{p}}{\partial \hat{x}}
\]
\[
\mu\frac{\partial^2 u}{\partial x^2} = \mu\frac{\partial^2[\Gamma\hat{u}/(\rho L)]}{\partial(L\hat{x})^2} = \mu\frac{\Gamma/(\rho L)}{L^2}\frac{\partial^2 \hat{u}}{\partial \hat{x}^2} = \frac{\Gamma^2}{\rho L^3}\frac{\partial^2 \hat{u}}{\partial \hat{x}^2}
\]
\[
\rho g\beta(T - T_1) = \rho g\beta(T_{\max} - T_1)\hat{T} = \rho g\beta(\Delta T)\hat{T}
\]
能量方程中的典型项:
\[
\frac{\partial}{\partial t}(\rho T) = \frac{\partial[\rho(T_1 + \Delta T\hat{T})]}{\partial(\rho L^2\hat{t}/\Gamma)} = \frac{\Gamma\Delta T}{L^2}\frac{\partial \hat{T}}{\partial \hat{t}}
\]
\[
\frac{\partial}{\partial x}(\rho u T) = \frac{\partial[\Gamma\hat{u}(T_1 + \Delta T\hat{T})/L]}{\partial(L\hat{x})} = \frac{\Gamma T_1}{L^2}\frac{\partial \hat{u}}{\partial \hat{x}} + \frac{\Gamma\Delta T}{L^2}\frac{\partial}{\partial \hat{x}}(\hat{u}\hat{T})
\]
\[
k\frac{\partial^2 T}{\partial x^2} = k\frac{\partial^2(T_1 + \Delta T\hat{T})}{\partial(L\hat{x})^2} = \frac{k\Delta T}{L^2}\frac{\partial^2 \hat{T}}{\partial \hat{x}^2}
\]
将各项替换为其等价表达式,得到连续方程、动量方程与能量方程的无量纲形式:
\[
\frac{\partial \hat{u}}{\partial \hat{x}} + \frac{\partial \hat{v}}{\partial \hat{y}} + \frac{\partial \hat{w}}{\partial \hat{z}} = 0
\]
\[
\frac{\partial \hat{u}}{\partial \hat{t}} + \frac{\partial}{\partial \hat{x}}(\hat{u}\hat{u}) + \frac{\partial}{\partial \hat{y}}(\hat{v}\hat{u}) + \frac{\partial}{\partial \hat{z}}(\hat{w}\hat{u}) = -\frac{\partial \hat{p}}{\partial \hat{x}} + \frac{\partial^2 \hat{u}}{\partial \hat{x}^2} + \frac{\partial^2 \hat{u}}{\partial \hat{y}^2} + \frac{\partial^2 \hat{u}}{\partial \hat{z}^2}
\]
\[
\frac{\partial \hat{v}}{\partial \hat{t}} + \frac{\partial}{\partial \hat{x}}(\hat{u}\hat{v}) + \frac{\partial}{\partial \hat{y}}(\hat{v}\hat{v}) + \frac{\partial}{\partial \hat{z}}(\hat{w}\hat{v}) = -\frac{\partial \hat{p}}{\partial \hat{y}} + \frac{\partial^2 \hat{v}}{\partial \hat{x}^2} + \frac{\partial^2 \hat{v}}{\partial \hat{y}^2} + \frac{\partial^2 \hat{v}}{\partial \hat{z}^2} + \text{Gr}\,\hat{T}
\]
\[
\frac{\partial \hat{w}}{\partial \hat{t}} + \frac{\partial}{\partial \hat{x}}(\hat{u}\hat{w}) + \frac{\partial}{\partial \hat{y}}(\hat{v}\hat{w}) + \frac{\partial}{\partial \hat{z}}(\hat{w}\hat{w}) = -\frac{\partial \hat{p}}{\partial \hat{z}} + \frac{\partial^2 \hat{w}}{\partial \hat{x}^2} + \frac{\partial^2 \hat{w}}{\partial \hat{y}^2} + \frac{\partial^2 \hat{w}}{\partial \hat{z}^2}
\]
\[
\frac{\partial \hat{T}}{\partial \hat{t}} + \frac{\partial}{\partial \hat{x}}(\hat{u}\hat{T}) + \frac{\partial}{\partial \hat{y}}(\hat{v}\hat{T}) + \frac{\partial}{\partial \hat{z}}(\hat{w}\hat{T}) = \frac{1}{\text{Pr}}\!\left(\frac{\partial^2 \hat{T}}{\partial \hat{x}^2} + \frac{\partial^2 \hat{T}}{\partial \hat{y}^2} + \frac{\partial^2 \hat{T}}{\partial \hat{z}^2}\right)
\]
其中 Gr 为 Grashof 数,Pr 为 Prandtl 数,\(\nu\) 为运动黏度,定义如下:
\[
\text{Gr} = \frac{g\beta\Delta T L^3}{\nu^2}, \quad \text{Pr} = \frac{\mu c_p}{k}, \quad \nu = \frac{\mu}{\rho}
\]
Grashof 数与 Prandtl 数是由所涉及的物理性质组合而成的无量纲群。因此,影响解的参数数量被缩减为两个,并可针对这两个参数的不同取值生成解。此外,只要这些物性组合所给出的 Gr 与 Pr 值与解所对应的值一致,单一解对多种物性组合都有效。这两个无量纲数与其他无量纲数的物理意义,将在下面结合其他无量纲化条件予以讨论。
3.9 无量纲数(Dimensionless Numbers)
将守恒方程写成无量纲形式会得到一些无量纲数,这些无量纲数对工程问题的参数研究极为有用。对不可压缩黏性流动,支配自然对流传热的无量纲参数被缩减为 Grashof 数 Gr 与 Prandtl 数 Pr 两个无量纲数。在不同条件下(如可压缩流动、多孔介质流动等),其他类型的流体作用力与耗散项会出现在控制方程中,从而产生不同的无量纲群。例如对多孔介质中的流动,Darcy 数 Da 作为重要参数出现;对自由表面流动,出现 Weber 数 We;对明渠流动,出现 Froude 数 Fr;对可压缩流动,出现 Mach 数 M,等等。下面讨论若干最重要的无量纲群。
3.9.1 Reynolds 数(Reynolds Number)
Reynolds 数定义为:
\[
\text{Re} = \frac{\rho U L}{\mu}
\]
可解释为对流(惯性)动量通量与扩散(黏性)动量通量之比的度量。若动量通量同方向,Reynolds 数揭示流动的边界层特性;若扩散通量定义为横向方向上的,则 Reynolds 数反映流动状态(层流、过渡或湍流),如图 3.10 所示。
图 3.11 展示了 Reynolds 数不同取值下流动场的一个示例。该示例为边长 \(L\) 的方形腔内、由顶壁以速度 \(U\) 驱动的流动。图中流线表明:随着 \(\text{Re} = \rho U L/\mu\) 增大,流动强度增强。
3.9.2 Grashof 数(Grashof Number)
如前推导,Grashof 数由下式给出:
\[
\text{Gr} = \frac{g\beta\Delta T L^3}{\nu^2}
\]
Grashof 数代表浮力与黏性力之比。它在自然对流中所扮演的角色与 Reynolds 数在受迫对流中所扮演的角色相同。图 3.12 展示了 Grashof 数影响的一个示例。该物理情形为水平偏心热圆管与冷方外壳之间环形区域内的自然对流传热。从图中等温线可看出:随着 Gr 值增大,等温线畸变更显著,反映了由更强流场引起的自然对流效应增强。
3.9.3 Prandtl 数(Prandtl Number)
Prandtl 数定义为动量扩散率(运动黏度 \(\nu\))与热扩散率(\(\alpha\))之比:
\[
\text{Pr} = \frac{\mu c_p}{k} = \frac{\mu/\rho}{k/\rho c_p} = \frac{\nu}{\alpha}
\]
Prandtl 数代表水动力学边界层与热边界层厚度之比。如图 3.13 所示:对 \(\text{Pr} < 1\)(图 3.13a),热边界层厚于水动力学边界层;对 \(\text{Pr} > 1\)(图 3.13b),则相反;二者厚度在 \(\text{Pr} = 1\) 时重合。
对前述方形腔内驱动流动问题,在 \(\text{Re} = 100\) 不变的条件下,不同 Pr 值下的等温线示于图 3.14。从图中可看出:随着 Pr 增大,对流相对导热的比例提高。
3.9.4 Péclet 数(Péclet Number)
Péclet 数定义为某物理量的对流输运率与扩散输运率之比。对传热情形,Péclet 数由下式给出:
\[
\text{Pe} = \frac{\rho U L c_p}{k} = \frac{UL}{\alpha} = \text{Re}\,\text{Pr}
\]
在这种情况下,Pe 等于 Reynolds 数与 Prandtl 数的乘积。图 3.15 展示了 Pe 的影响一例——平板上方不同 Péclet 数下的等温线分布。低 Pe 值下,热量传递以导热为主;随着 Pe 增大,对流的作用逐渐增强,并在 \(\text{Pe} = 1000\) 时成为主导传热机制。
对质量输运,Péclet 数由下式给出:
\[
\text{Pe} = \frac{UL}{D} = \text{Re}\,\text{Sc}
\]
其中 \(D\) 为质量扩散率,Sc 为 Schmidt 数。此时 Pe 等于 Reynolds 数与 Schmidt 数的乘积。
大的 Péclet 数意味着流动对下游位置的依赖较低,而对上游位置的依赖较高。因此,对高 Péclet 数的情形可采用相对简化的计算模型。
3.9.5 Schmidt 数(Schmidt Number)
Schmidt 数定义为:
\[
\text{Sc} = \frac{\nu}{D}
\]
Schmidt 数在质量传递中是 Prandtl 数在热量传递中的对偶量,代表动量扩散率(\(\nu\))与质量扩散率(\(D\))之比。物理上,Sc 关联水动力学边界层与质量传递边界层的厚度。图 3.16 展示了 Schmidt 数影响的一个示例。
该图表示两根水平横截面为菱形管之间环形区域内的自然对流质量传递。溶质浓度在外管内壁附近较高。浓度不均匀导致密度变化,进而建立流场。随着 Sc 值增大,流动强度增强——表现为等浓度线畸变更显著,说明对流质量传递相对扩散质量传递增强(后者在低 Sc 值下占主导)。
3.9.6 Nusselt 数(Nusselt Number)
Nusselt 数由下式表达:
\[
\text{Nu} = \frac{hL}{k}
\]
是对流传热系数 \(h\) 的无量纲形式,提供固体表面对流传热强度的度量。Nusselt 数并非在守恒方程无量纲化过程中自然出现的无量纲群;它通常被广泛用于报告对流传热数据。
3.9.7 Mach 数(Mach Number)
Mach 数定义为运动物体相对流体的速度与当地声速之比:
\[
M = \frac{|\mathbf{v}|}{a}
\]
其中 \(|\mathbf{v}|\) 为流体相对其流动介质的当地速度幅值,\(a\) 为声速。声速的一般方程为:
\[
a = \sqrt{\left.\frac{\partial p}{\partial \rho}\right|_{T}}
\]
对理想气体,简化为:
\[
a = \sqrt{\gamma R T}
\]
其中 \(\gamma\) 为定压比热与定容比热之比 \(c_p/c_v\),\(R\) 为气体常数。
Mach 数小于 0.2 的流动可视为不可压缩。\(M < 1\) 时流动为亚声速,\(M = 1\) 时为声速,\(1 < M < 5\) 时为超声速,\(M > 5\) 时为高超声速。从亚声速加速到超声速的流动称为跨声速流动。计算域边界上 Mach 数的大小(小于 1 或大于 1)决定该边界所需的边界条件个数。
图 3.17 通过 Mach 等值线给出亚声速、跨声速与超声速流动场的示例。该物理情形为流体绕过一个最大曲率为通道高度 10%(亚、跨声速情形)或 4%(超声速情形)的圆弧形凸起的流动。从图中可看出:随着 Mach 数从亚声速(\(M < 1\))增加到超声速(\(M > 1\)),流动类型由椭圆型变为双曲型(并以激波形式出现间断)。
3.9.8 Eckert 数(Eckert Number)
Eckert 数是关联流动动能与其焓值的无量纲数,定义为:
\[
\text{Ec} = \frac{\mathbf{v}\cdot\mathbf{v}}{c_p \Delta T}
\]
其中 \(\Delta T\) 为特征温差。该无量纲数在对可压缩能量方程无量纲化时,作为黏性耗散项 \(U\) 的乘子出现。Ec 值大意味着高速流动(高动能)下显著的黏性耗散。当 Eckert 数较小(\(\text{Ec} \ll 1\))时,能量方程中若干项(如黏性耗散项、体积力项等)可忽略,能量方程化为其不可压缩形式(即仅含导热与对流的平衡)。
3.9.9 Froude 数(Froude Number)
Froude 数定义为特征速度 \(U\) 与重力波速度 \(\sqrt{gL}\) 之比:
\[
\text{Fr} = \frac{U}{\sqrt{gL}}
\]
它是流体中部分浸没物体所受流体阻力的一种度量,Fr 值越高表示阻力越大。
对图 3.18 所示的自由表面流动,流动的性质由 Froude 数决定。\(\text{Fr} > 1\) 时流动为超临界(supercritical),\(\text{Fr} < 1\) 时为亚临界(subcritical)。两个区域之间的界面——即所谓的水跃(hydraulic jump)——为临界状态,其 Froude 数值为 1。
3.9.10 Weber 数(Weber Number)
Weber 数定义为:
\[
\text{We} = \frac{\rho U^2 L}{\sigma}
\]
其中 \(U\)(m/s)与 \(L\)(m)分别为特征速度与特征长度,\(\sigma\) 为表面张力(N/m)。Weber 数代表惯性力与表面张力之比,对分析涉及不同流体之间界面、且具有弯曲表面(如液滴、气泡)的多相流动非常有用。
3.10 小结(Closure)
本章说明了多种物理现象可通过守恒方程建模。这些方程由有限体上的平衡出发,基于第一原理推导而出。本章还表明,支配质量、动量、能量以及其他比量输运的守恒方程具有共同的形式——即通用标量输运方程。该方程包含瞬态项、对流项、扩散项与源项,每一项都为方程贡献特征性的作用,需要由离散化过程加以再现。
本章个人批注
本章是 Moukalled FVM 一书中最"硬核"的一章——它将流体动力学所涉及的全部守恒律,从质量、动量到能量,以一致的方式重新整理为通用输运方程的形式。在阅读中我有几点具体的感受与疑问,留待后续章节与文献验证。
第一,3.3.1 节给出的物质导数推导非常标准,但作者强调"将注意力放在固定点/固定体积上"的工程取向值得关注。Lagrange 描述在理论物理(比如将流体作为经典场论对象来分析)时反而更自然;FVM 之所以强制选择 Euler 描述,是因为控制体的边界与网格重合,便于积分。这点在后续第 5 章会反复体现——读者最好从一开始就将"控制体 = 网格单元"作为思维框架。
第二,3.4 节末尾对不可压缩性的讨论写得很克制,但有一个微妙之处:\(D\rho/Dt = 0\) 不等于 \(\rho\) 在空间上均匀——它只意味着每个微元在运动中保持原密度。文中提到水中的盐浓度差与空气中的温差会引起密度差并产生浮力,这正是 3.8 节 Boussinesq 近似的动机。这条线索如果跟下去,会连接到 3.9.2 Grashof 数与自然对流。这里我有一个问题:Boussinesq 近似只对 \(\Delta T^2\) 及更高阶项做截断,但许多工程应用(比如低温流体或大温差)下截断并不成立;Moukalled 在后续章节是否会讨论非 Boussinesq 形式的自然对流?
第三,3.5.3 节对应力张量
\(\boldsymbol{\mathcal{R}} = -p\mathbf{I} + \boldsymbol{\sigma}\)
的分解,是连续介质力学中标准的 Cauchy 分解。但作者在写作上做了一个隐含约定:把"压力"定义为法向应力平均值的负值(式 (3.27))。在静止流体中这与热力学压强一致;在有黏性应力或非平衡情形下,二者可能不同——这在多孔介质流动或非 Newtonian 流体中尤其突出。第 14 章处理源项时会再次涉及这个区别。
第四,3.5.5 节针对 Newtonian 流体将 \(\lambda\) 取为 \(-(2/3)\mu\),这是 Stokes 假设。许多教材(包括 Bird 等的 Transport Phenomena)都明确指出该假设在某些情况下不成立(例如单原子气体的快速压缩)。Moukalled 没有讨论 Stokes 假设的适用条件——读者最好将其作为工程近似,而非严格物理定律。
第五,3.6.4 节推导温度形式的能量方程时,我注意到对密度为常数的情形
\[
\partial(\ln\rho)/\partial(\ln T) = 0
\]
,这一项被略去;但对理想气体
\[
\partial(\ln\rho)/\partial(\ln T) = -1
\]
——式 (3.81) 中 \(Dp/Dt\) 之前出现 \(+\) 号而非 \(-\) 号,这一符号差异是值得在数值实现时仔细核对的。
第六,3.8 节的 Boussinesq 无量纲化示例非常清晰,但作者刻意只对一个特例(自然对流)展开推算。受迫对流与可压缩流动的无量纲化会导出哪些无量纲数,本章没有给完整推导,只是通过 3.9 节以"定义 + 物理意义"的方式介绍。这种取舍对工程读者是合理的——他们关心无量纲数如何读图、如何在参数研究中分组;对理论读者则略显不足。Reynolds 数与 Grashof 数分别对应受迫与自然对流的"惯性 vs. 黏性"或"浮力 vs. 黏性"比值,是非常工整的对偶。
第七,3.9 节关于 Mach 数的讨论有一处隐含的工程规则:\(M < 0.2\) 时可视为不可压缩。这一阈值在后续章节(如第 16 章可压缩流动)应当会被再次提及。Moukalled 在第 1 章说过这本书以不可压缩为主、可压缩为辅,所以 Mach 数在本章出现,更多是为了给后续章节做铺垫——读者不需要在这里做完整推导。
最后,3.7 节的"通用守恒方程"是整本书最具战略价值的一段。作者明确指出:通过为 \(\phi\)、\(\Gamma_\phi\) 与 \(Q_\phi\) 赋值,式 (3.93) 可以代表任意守恒方程。这意味着后续第 5–13 章的离散化方法只需针对一个通用方程展开,不必对每个守恒方程单独推导。这是有限体积法相比有限差分法与有限单元法在工程实现上的一大优势——也是 Moukalled 全书组织的核心论点。
与上下章的衔接(一段话)
第 2 章复习了矢量计算与张量微积分的基础工具,第 3 章则使用这些工具系统地建立了流体动力学全部守恒方程的数学形式:从 Reynolds 输运定理出发推导出质量、动量、能量方程的 Euler 形式,再统一为一个含瞬态、对流、扩散与源项的通用标量输运方程。本章占据了全书约 80 页的篇幅,是后续所有数值方法(离散化、网格、边界条件、对流格式、时间步进)的物理基础——后续章节(第 4 章离散化过程、第 5 章 FVM 的基本方程、第 8 章扩散项的空间离散等)都将回到本章给出的通用方程 (3.93),对其四项分别构造离散近似。第 4 章即将从通用方程的"物理形式"过渡到"离散代数形式",给出从连续 PDE 到代数方程组的整体映射。