化工流体力学第六章

1、第六章第六章 计算流体力学计算流体力学 本章要点：本章要点： 1. 流体力学主要方程与模型流体力学主要方程与模型 2.数值方法（有限元法）数值方法（有限元法） 3.节典型流场（绕流，管流）的数值计算节典型流场（绕流，管流）的数值计算 4.节鼓泡反应器中流场的数值模拟节鼓泡反应器中流场的数值模拟概概 述述 计算流体力学采用数值模拟来研究流体运动规律，具有耗费少，时间短，重复性好等优点。解决问题思路：遵循以下步骤（1）进行简化以确定合理的物理模型。 （2）根据所确定物理模型建立相应的数学模型。 （3）对数学模型中各偏微分方程组选用合适的方法进行离散成代数方程组，用恰当的格式进行求解。 （4）对所得

2、的数值结果进行验证。6.1 流体力学主要方程与模型不可压N-S方程 本节主要介绍不可压N-S方程的统一形式，并推导出涡流函数方程。并简要介绍湍流模型。在现代计算流体力学里，N-S方程写成如下统一形式：式中各项代表意义与第三章相同。涡流函数方程对于二维不可压黏性流动，忽略体积力，其直角坐标下无量纲方程形式：JzHyGxFtU0yuxuyx动量方程中有压力项，为了消去压力项，引入流函数：这时，连续方程自动满足。将上式带入下二维涡量的定义式中：可得泊松方程：涡量传递方程为：上两式即为涡流函数方程的N-S2222Re11yuxuxpyuuxuutuxxxyxzx2222Re11yuxuypyuuxuu

3、tuyyyyyxyxuyuyx,yuxuxy2222yx2Re1yuxutyx湍流模型：第三章所述的雷诺方程不封闭，因此必须补充方程，所补充的方程有零方程、一方程和二方程模型，还有代数应力模型（ASM）和雷诺应力方程模型（DSM）等。单方程模型：湍流动能方程形式如下：式中 为湍流产生项，k为湍流脉动动能的普朗特数，常取为1，cD为常数，湍流的长度标尺l可由经验公式得到。单方程模型考虑到了湍流动能的对流与扩散，湍流粘度联系到整个流场。双方程模型：lcGxxxutDkiktii2/3jiijijtkxuxuxuG湍流耗散率常用形式如下：湍流动能方程可改为如下形式：模型中各经验常数取值如下：c=0.

4、9，c1=1.44，c2=1.92，k=1.0，=1.3，T=1.0。和是表征湍流发生、发展与耗散的重要特征量，是应用最广泛而有成效的湍流模型。应力输运方程模型（DSM）：雷诺应力输运方程中需要模拟的有湍流扩散项和压力应变项。湍流扩散项其梯度模型为：21cGcxxxutkitiiikiktiiiGxxxut1,1xuuuuxCpupuuuuxDjilkksjikijkkjikijt湍流动能耗散率的输运方程类似于标准k-模型中的方程：式中=1.0，C1=1.44，C2=1.92C3.。湍流的高级模拟：对于高雷诺数湍流流动，涡旋近似于均匀各向同性，当只需描述湍流的统计平均量时，从雷诺平均方程出发，

5、引入雷诺应力的封闭模型进行模拟，称为雷诺平均数值模拟（RANS）。直接模拟为了对湍流的物理性质进行深入的探讨，需要对描述湍流瞬时运动控制方程直接进行最为精细的数值模拟计算，不需要补充模型，这种称为直接数值模拟（DNS）。大涡模拟（LES）综合了直接模拟与雷诺平均法两种模拟思想，是一种介于DNS和RANS之间的方法。对大涡进行直接模拟，对小涡用模型进行模拟。21CPCxxxutiiitiii6.2 数值计算方法目前，计算流体力学常用的数值计算方法有：有限差分法、有限体积法、有限元法、谱方法、摄动法、边界元法及有限分析法等。6.2.1 微分方程的离散和求解数值方法的求解过程包括将定解问题的计算区域

6、用网格划分为有限个网格节点，将微分方程转换为节点所对应的离散代数方程组，利用初始边界条件对离散方程组进行计算求解。要求对定解问题进行离散求解时，对构造离散方程的格式要求具有相容性和稳定性，其数值解具有收敛性。在保证相容性的前提条件下，通过选取适当的数值方法对定解问题控制方程组进行离散，需满足稳定性条件。常用的对离散方程的求解方法通常分两类：直接求解法和迭代法。6.2.2 有限差分法 构造差分方程的方法有很多，最常用的为泰勒级数展开法。用ui，j表示控制区域内点（i，j）在x方向速度，则点（i+1，j）的速度可用泰勒级数展开成关于点（i，j）速度的表达式：上式可写成：式中右端第一项为 在点（i，

7、j）的差分表达式，其他项则为截断误差。在截断误差表达式中，关于x的最低幂次数为1，因而上式又称为一阶精度表达式。写成以下形式：.62)(3,332,22, 1xxuxxuxxuuujijijijiji.62)(2,33,22, 1, xxuxxuxuuxujijijijiji.xu)()(, 1,xOxuuxujijiji通常情况下一阶精度不能满足需要，因而需构造二阶甚至更高阶的差分表达式，有：同样的，两阶偏导的两阶中心差分公式为：用有限差分法构造描述物理问题微分方程的差分离散格式有许多，但通常必须满足相容性与稳定性的要求。6.2.3 有限元法比较适用于求解复杂计算区域的流体流动问题。有限元法

8、将计算区域连续划分成许多离散的网格单元，在网格单元内，将微分方程中的变量改写成由各变量或其导数的节点值与所选用的插值函数组成的线性表达式，借助于加权余量法或变分原理，将控制微分方程离散成代数方程组进行求解。2, 1, 1,)(2)(xOxuuxujijiji22, 1, 1,22)()(2)(xOxuuuxujijijiji 在有限元分析中，插值函数j一般用不同幂次的多项式函数表示，有拉格朗日插值和埃尔米特插值。以一维空间为例采用拉格朗日插值函数如下：式中，n为单元节点数， j为（n-1）次逼近的插值函数。为了简便可在一个单元只取2个节点，则由上拉格朗日插值函数式，可给出下图所示线性分布，在图

9、中单元内，有：因而，对于一维线性插值，在有限元中的近似函数u表达为： njjjjjjnjjjxxxxxxxxxxxxxxxxx.).(.)(1111111111,jjjjjjjjxxxxxxxxjjjjuuu11在选取单元插值函数后，还需选择权函数或变分原理来构建离散方程组。下面介绍权函数的选取方式：配置法：配置法选取Dirac函数作为权函数，其定义为：对于区域中所选的配置点xi，性质如下：设r为任意小的正数，记i为小区域|x-xi|r，则有：可选取n个线性无关的权函数设算子方程为：L（u）=p，u。方程的近似解可以表示为n个基函数的线性组合: iiiixxxxxxxx,0,0,1,iidxx

10、xxdxxxxiiiinixxwii,.,2 , 1, njjjnau1上式中aj是待定系数，令方程的余量分别与wi正交，即它们的内积为零，因而有：方程中的余量一般不为零，将基函数的线性组合代入上式中，求解代数方程组可求得系数aj，可得近似解。6.2.4 有限体积法有限体积法是以守恒型控制方程为出发点，把计算域分成许多控制容积，并对每个控制容积作为积分来构造离散方程。为了建立离散方程，需对求解区域网格及节点相关的几何要素的命名方法规范化。如图所示的一维网格。nidwpuLwini,.2 , 1, 0, 控制体积法中一般将节点用P表示，其网格的相邻节点为E与W。虚线表示控制面，e，w表示相应的界

11、面。相邻两节点间的距离以x表示，x表示相邻两界面间的距离，图中控制体积为x*1*1。以热扩散系数为常数的一维非定常热传导方程为例说明用有限体积法进行方程离散的方法。一维非定常热传导方程如下：上式中，假设热扩散系数和相邻距离x为常数。在控制容积内积分上述方程得：为了导出离散方程，作出两个假设即：假设网格点P上的T值代表整个控制体积上的值，即阶梯式分布；假设T随x呈分段线性分布。则可得到：22xTtT ewttttttewdxdtxTxdtdxtT其中：上式即为一维非定常热传导方程用有限体积法离散成代数方程的通用形式。现针对不可压黏性流N-S方程的二维守恒形式利用有限体积法进行离散，则二维不可压粘

12、性流标准方程为式中S为方程源项；为通用变量。x与y方向的对流-扩散总通量表达式分别如下：0000111pWEpWwwEEEppTaafaTffTaTffTaTa00,pWEppwWeEaaafatxaxaxa SyuyxuxtyxyuJxuJyyxx,用有限体积法来离散偏微分方程，将通式即二维粘性流标准方程在以P点位代表的主控制容积内作时间和空间的积分，并作出假设。将通式在主控制容积面内积分，有：在关于P点的主控制容积面内积分连续性方程，有：式中，Fe，Fw，Fn及Fs分别是通过主控制容积面e，w，n，s的体积流率，既有：上三式合并得到：yxSSJJJJtyxppcsnwepp00snweFF

13、FFxuFxuFyuFyuFssnnwwee, yxSSFJFJFJFJtyxppcpsspnnPwwpeepp0psspssnpNpnnpwwpwwEpEpeeAFJAFJAFJAFJ,0 ,0 ,0 ,0 ,ssssssnnnnnNewwwwweeeeeEFPADPADAFPADPADAFPADpADAFPADPADAxxDxxDyxDyxDssnnwwee ,将以上方程带入合并式中，可得到二维不可压黏性流标准方程在节点P的离散方程。6.3 流场数值计算的应用实例圆柱绕流：低雷诺数圆柱绕流电磁控制的数值计算（1）数学模型：其涡量和流函数的表达式为：考虑壁面无滑移，找出初始条件和边界条件，数

14、值计算时，涡量函数采用ADI格式，流函数采用FFT格式。下图为t=400时，流场的涡量与流线等值分布图，展示了充分发展后的圆柱卡门涡街尾迹情况。HFFNHUUtHr22222122222Re2非牛顿流体管流：非牛顿流体收缩管流与管道输运和血液在动脉中流动等现象息息相关。控制方程组：设流体均匀不可压，血流是周期性的脉动流动，则非牛顿流体无量纲涡流函数方程：式中，和分别为流函数与涡量，St为斯托克斯数。边界条件：脉动流，其壁面上的流函数用下式表示：壁面无滑移，对于定常流动，壁面上=0.5，进口处，流体已充分发展。2222222221112Re1rzrrrzrzrrrzrzrrzzrrtStrrrr

15、zrzrzrrzzrzrzrrtwcos15 . 00zz 在出口处：设管中流动关于轴线对称，只计算上半部分的管流，计算域为-10z10，z是以管入口半径a为特征长度的无量纲数。计算结果：其壁面压力与剪切应力分布如下两图所示，由图可知，狭窄管流在周期脉动过程中，壁面剪切应力在狭窄处的变化很大，壁面压力则是在狭窄端下游变化较大。02222zz 当Re=200时，非牛顿流体在脉动狭窄管中流动的瞬时流线分布如下图所示。下图显示了一个完整周期中非牛顿流体的流线在狭窄管中的变化，t=0时，狭窄区下游有一个大的分离区，t=时，流动被减速，净通量为零。在t=7/5时，流体流动变的平滑，直到t=8/5时，流体

16、被加速，新的涡旋出现并被拉长，新周期开始 伴有化学反应的边界层流：研究对翅片表面燃料气流燃烧反应。控制方程：翅片坐标系统如下图所示，气流通过平板形成边界层，燃料在平板表面反应产生热量。边界层内的运动、能量和组分方程分别为：设定边界条件，以上模型可以用来预测燃料燃烧情况下翅片温度。iiiiyixiiifyxpxxyxxWyYDyyYuxYuhWyTyyTuxTuCyuyyuuxuu 计算结果如下图，下左图与中间图分别为翅片无量纲温度（平均温度）p与主流温度的变化情况。下右图为翅片表面边界层温度的分布图。鼓泡塔中汽液两相流的数值计算：鼓泡塔是汽液两相进行质量、动量和能量传递的重要设备，广泛应用于化

17、工、石油及生物化工等领域中。双流体模型:鼓泡塔一般为圆柱形结构，气体从底部进入鼓泡床后，产生垂直向上运动的气泡，一般情况下气泡在塔内分布不均匀，可用双流体模型（欧拉-欧拉模型）来描述汽液两相流动。作出假设：气泡与液体可视为互相穿透的连续介质，两种流体以一定几率出现，用气含率与液含率表征，汽液两相流的基本方程由连续性方程与动量方程组成。对于含有Np相的流体其连续方程如下：式中右端项描述为由相到的质量传递。r 代表相的体积分数。忽略相间质量传递与源相：pNiiSrmmurxrt1,0,iiurxrt11pNiir动量方程为：方程右端描述了所有作用在控制体上相流体上的力，总压力项，黏性应力项以及重力

18、与相间相互作用项F,i。只考虑阻力的影响，可根据Clift等所得的公式：式中，阻力系数CD=0.44,；气泡直径db=4mm为常数； 为气泡相对于流体的滑移速度； 为流体流动速度。对于连续相可用单相标准湍流模型k-e模型来模拟汽液两相流中连续相中的湍流现象。其湍流动能k与湍流耗散率e方程如下：iiijjiejijiiiFgrxuxurxxpruurxurt,uuuudrCFmmbDi143, uummuk-e模型中的常数可查表得。汽液流中气泡诱导的湍流对流动有一定影响，可通过改进湍流动能与湍流耗散率产生项G体现，该项描述相间的相互作用。如下式：式中， = =1.44与 = = 1.92。Grxkrxkurxkrtieiii,21,CGCkrxrxurxrt