CFD基础与Fluent工程应用分析

978-7-115-59286-6
作者: 江帆
译者:
编辑: 李瑾
分类: Ansys

图书目录:

详情

本书是 CFD 基础及 Fluent 工程应用分析的指导性教材或参考书。本书根据ANSYS Fluent 2021编写,并结合Mesh模块的网格划分方法,详细介绍计算流体力学、多相流、动网格、滑移网格、流固耦合等基础知识,涵盖Mesh模块内复杂几何体的网格划分实例,以及Fluent软件在各个工程领域的实际应用分析。全书共 12 章:第 1章为计算流体力学基础,第 2 章为 Fluent 简介,第 3 章为网格划分,第 4 章为稳态与瞬态流动分析,第 5 章为离散相流动分析,第 6 章为传热流动分析,第 7 章为多孔介质流动分析,第 8 章为多相流动分析,第 9 章为动网格流动分析,第 10 章为滑移网格流动分析,第11 章为流固耦合分析,第 12 章为化学反应、燃烧与微流动分析。书中以详解实例的方式来说明 Fluent 软件在各个工程应用领域的详细操作及一些值得注意的问题,设有大量练习题,具有较强的实用性。 本书可作为机械、材料、暖通、水利、动力、能源、航空、冶金、环境、建筑等相关专业的本科生和研究生的教材,也可供上述领域的技术人员,特别是进行 CFD 应用研究的人员参考。

图书摘要

版权信息

书名:CFD基础与Fluent工程应用分析

ISBN:978-7-115-59286-6

本书由人民邮电出版社发行数字版。版权所有,侵权必究。

您购买的人民邮电出版社电子书仅供您个人使用,未经授权,不得以任何方式复制和传播本书内容。

我们愿意相信读者具有这样的良知和觉悟,与我们共同保护知识产权。

如果购买者有侵权行为,我们可能对该用户实施包括但不限于关闭该帐号等维权措施,并可能追究法律责任。


版  权

编  著 江 帆 谢宝山 张冥聪 黄春曼

责任编辑 李 瑾

人民邮电出版社出版发行  北京市丰台区成寿寺路11号

邮编 100164  电子邮件 315@ptpress.com.cn

网址 http://www.ptpress.com.cn

读者服务热线:(010)81055410

反盗版热线:(010)81055315

内容提要

本书是CFD(计算流体力学)基础,Fluent工程应用分析及Mesh模块网格划分实践的指导性教材或参考书。全书共12章:第1章为计算流体力学基础,第2章为Fluent简介,第3章为网格划分,第4章为稳态与瞬态流动分析,第5章为离散相流动分析,第6章为传热流动分析,第7章为多孔介质流动分析,第8章为多相流动分析,第9章为动网格流动分析,第10章为滑移网格流动分析,第11章为流固耦合分析,第12章为化学反应、燃烧与微流动分析。书中以详解实例的方式来说明Fluent软件在各个工程应用领域的详细操作及一些值得注意的问题,设有大量练习题,具有较强的实用性。

本书可作为机械、材料、暖通、水利、动力、能源、航空、冶金、环境、建筑等相关专业的本科生和研究生的教材,也可供上述领域的技术人员,特别是进行CFD应用研究的人员参考。

前  言

CFD(计算流体力学)是融合流体力学、数值数学和计算机科学的一门新兴交叉学科,它从计算方法出发,利用计算机快速的计算能力得到流体控制方程的近似解。CFD兴起于20世纪60年代,随着20世纪90年代后期计算机的迅猛发展,CFD得到了飞速发展,成为很多领域的一个重要分析手段。CFD软件有很多,其中Fluent是比较流行的商用CFD软件包,适用于流动、传热和化学反应等相关领域。它包含基于压力的分离求解器、基于密度的隐式求解器与显式求解器、多求解器等技术,可以用来模拟从不可压缩到高度可压缩范围内的各种复杂流场。它具有强大的前后处理功能,在汽车外形设计、机翼气流、炉内燃烧、石油运输、血液流动、注塑、机加工冷却、暖通空调、建筑、环境保护等方面有着广泛的应用。

为了让读者掌握CFD基础知识与Fluent各方面的应用操作,编者根据ANSYS Fluent 2021编写了本书。本书的内容包含计算流体力学、多相流、动网格、滑移网格、流固耦合等基础知识,以及Fluent软件在各个工程领域的实例分析、Mesh模块的网格划分方法及实例。这些实例涵盖稳态与瞬态流动分析,离散相流动分析,传热流动分析,多孔介质流动分析,多相流动分析,动网格流动分析,滑移网格流动分析,流固耦合分析,化学反应、燃烧与微流动分析等方面。部分实例还给出了实验结果对比,便于读者进一步研究。

全书分12章,第1章为计算流体力学基础,第2章为Fluent简介,第3章为网格划分,第4章为稳态与瞬态流动分析,第5章为离散相流动分析,第6章为传热流动分析,第7章为多孔介质流动分析,第8章为多相流动分析,第9章为动网格流动分析,第10章为滑移网格流动分析,第11章为流固耦合分析,第12章为化学反应、燃烧与微流动分析。书中讲解大量工程实例及Fluent使用过程中值得注意的一些问题,具有较强的实用性。书中还设有大量练习题(配套资源有相关的模型、网格与CAS文件)供读者强化练习,熟悉Fluent相关操作。

本书由江帆、谢宝山、张冥聪、黄春曼编著,江帆构思全书框架并筛选实例,江帆、黄春曼编写第1章,谢宝山编写第2~8章及第12章,张冥聪编写第9~11章,全书由江帆、谢宝山校对与统稿。

本书多数实例来源于编者的科研项目和与企业合作的项目,部分材料来源于网络资源,感谢相关人员所做的工作。

本书的编写得到了广州市科技计划项目“稠油运输柔性管内油水环状流稳定性研究”(No.202102010386)和广州大学机械与电气工程学院老师们的支持,还得到了编者家人的理解与大力支持,在此一并致以深深的谢意!

本书既可作为使用Fluent软件进行工程应用分析的技术人员的参考书,又可作为高等院校相关专业本科生和研究生的CFD分析教学、复杂流动问题研究的参考书,书中实例的SAT、STEP、MSH、CAS、DAT等文件及相关的源程序、操作视频已上传人民邮电出版社异步社区,扫描书中二维码即可获得相关资料。

限于编者的水平,书中难免有不当之处,还请广大读者给予指正,来信请致jiangfan2008@126.com或xiebaoshan2022@163.com,不胜感激。

江帆

2022年1月于广州

资源与支持

本书由异步社区出品,社区(https://www.epubit.com)为您提供相关资源和后续服务。

配套资源

本书提供如下资源:书中实例的SAT、STEP、MSH、CAS、DAT等文件,以及相关的源程序和操作视频。

要获得以上配套资源,请在异步社区本书页面中单击,跳转到下载界面,按提示进行操作即可。注意:为保证购书读者的权益,该操作会给出相关提示,要求输入提取码进行验证。

提交错误信息

作者和编辑尽最大努力来确保书中内容的准确性,但难免会存在疏漏。欢迎您将发现的问题反馈给我们,帮助我们提升图书的质量。

当您发现错误时,请登录异步社区,按书名搜索,进入本书页面,单击“提交勘误”,输入错误信息,单击“提交”按钮即可。本书的作者和编辑会对您提交的错误信息进行审核,确认并接受后,您将获赠异步社区的100积分。积分可用于在异步社区兑换优惠券、样书或奖品。

扫码关注本书

扫描下方二维码,您将会在异步社区微信服务号中看到本书信息及相关的服务提示。

与我们联系

我们的联系邮箱是contact@epubit.com.cn。

如果您对本书有任何疑问或建议,请您发邮件给我们,并请在邮件标题中注明本书书名,以便我们更高效地做出反馈。

如果您有兴趣出版图书、录制教学视频,或者参与图书翻译、技术审校等工作,可以发邮件给我们;有意出版图书的作者也可以到异步社区在线投稿(直接访问www.epubit.com/ contribute即可)。

如果您所在的学校、培训机构或企业,想批量购买本书或异步社区出版的其他图书,也可以发邮件给我们。

如果您在网上发现有针对异步社区出品图书的各种形式的盗版行为,包括对图书全部或部分内容的非授权传播,请您将怀疑有侵权行为的链接发邮件给我们。您的这一举动是对作者权益的保护,也是我们持续为您提供有价值的内容的动力之源。

关于异步社区和异步图书

“异步社区”是人民邮电出版社旗下IT专业图书社区,致力于出版精品IT专业图书和相关学习产品,为作译者提供优质出版服务。异步社区创办于2015年8月,提供大量精品IT专业图书和电子书,以及高品质技术文章和视频课程。更多详情请访问异步社区官网https://www.epubit.com。

“异步图书”是由异步社区编辑团队策划出版的精品IT专业图书的品牌,依托于人民邮电出版社数十年的计算机图书出版积累和专业编辑团队,相关图书在封面上印有异步图书的LOGO。异步图书的出版领域包括软件开发、大数据、人工智能、测试、前端、网络技术等。

异步社区

微信服务号

第1章 计算流体力学基础

计算流体力学(Computational Fluid Dynamics,CFD)将流体力学控制方程中的积分、微分项近似地表示为离散的代数形式,使其成为代数方程组,然后通过计算机求解这些离散的代数方程组,获得离散的时间点与空间点上的数值解。CFD是融合流体力学、数值数学和计算机科学的交叉学科,已广泛应用于热能动力、土木水利、汽车、铁道、船舶工业、航空航天、石油化工、流体机械、环境工程等领域。工欲善其事,必先利其器。下面先简单介绍CFD基础知识。

1.1 流体力学的基本概念

1.1.1 流体的连续介质模型

流体质点:几何尺寸同流动空间相比是极小量,又含有大量分子的微元体。

连续介质:质点连续地充满所占空间的流体或固体。

连续介质模型:把流体视为没有间隙地充满它所占据的整个空间的一种连续介质,且其所有的物理量都是空间坐标和时间的连续函数的一种假设模型,即u=u(t,x,y,z)。

1.1.2 流体的性质

惯性:流体不受外力作用时,保持其原有运动状态的属性。惯性与质量有关,质量越大,惯性就越大。单位体积流体的质量称为密度(Density),用表示,单位为kg/m3。对于均质流体,设其体积为V,质量m,则密度为

(1-1)

对于非均质流体,密度因点而异。若取包含某点在内的体积与质量,则该点密度用极限方式表示,即

(1-2)

压缩性:作用在流体上的压力变化可引起流体的体积变化或密度变化。压缩性可用体积压缩率k来量度,即

(1-3)

式中,p为外部压强。

在研究流体流动的过程中,若考虑流体的压缩性,则为可压缩流动,相应的流体称为可压缩流体,例如高速流动的气体;若不考虑流体的压缩性,则为不可压缩流动,相应的流体称为不可压缩流体,如水、油等。

黏性:在运动的状态下,流体所产生的抵抗剪切变形的性质。黏性大小由黏度来量度。流体的黏度是由流动流体的内聚力和分子的动量交换所引起的。黏度有动力黏度和运动黏度之分。动力黏度由牛顿内摩擦定律导出

(1-4)

式中,为切应力,单位为Pa;为动力黏度,单位为Pa•s;为流体的剪切变形速率。

运动黏度与动力黏度的关系为

(1-5)

式中,为运动黏度,单位为m2/s。

在研究流体流动的过程中,若考虑流体的黏性,则为黏性流动,相应的流体称为黏性流体;若不考虑流体的黏性,则为理想流体的流动,相应的流体称为理想流体。

根据流体是否满足牛顿内摩擦定律,可以将流体分为牛顿流体和非牛顿流体。牛顿流体严格满足牛顿内摩擦定律且保持为常数。非牛顿流体的切应力与速度梯度不成正比,非牛顿流体一般又分为塑性流体、假塑性流体、胀塑性流体3种。

塑性流体,如牙膏等,有一个保持不产生剪切变形的初始应力,只有克服了这个初始应力,其切应力才与速度梯度成正比,即

(1-6)

假塑性流体,如泥浆等,其切应力与速度梯度的关系是

(1-7)

胀塑性流体,如乳化液等,其切应力与速度梯度的关系是

(1-8)

1.1.3 流体力学中的力与压强

质量力:与流体微团质量大小有关并且集中在微团质量中心的力。在重力场中有重力;直线运动中,有惯性力。单位质量力是一个矢量,一般用单位质量所具有的质量力来表示,其形式如下

(1-9)

式中,ijk为单位质量力在坐标轴上的投影。

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

对于理想流体的流动,流体质点所受到的作用力只有正压力,没有切向力。对于黏性流体的流动,流体质点所受到的作用力既有正压力,又有切向力。

作用在静止流体上的表面力只有沿表面内法线方向的正压力。单位面积内所受到的表面力称为这点的静压强。静压强具有两个特征:①静压强的方向垂直指向作用面;②流场内某点静压强的大小与方向无关。

液体的表面张力:作用于液体表面的液体间的相互作用力。液体表面有自动收缩的趋势,收缩的液面存在与该处液面相切的拉力。正是这种力的存在,才会出现弯曲液面内外出现压强差及常见的毛细现象等。

实验表明,液体表面张力的大小T与液面的截线长度L成正比,即

(1-10)

式中,称为表面张力系数,它表示液面上单位长度截线上的表面张力,其大小由液体性质与接触相温度、压力等决定,其单位为N/m。

标准大气压的压强是101325 Pa(760 mmHg),记作atm。若压强大于标准大气压,则以此压强为计算基准得到的压强称为相对压强,也称表压强,通常用表示。若压强小于标准大气压,则压强低于标准大气压的值就称为真空度,通常用表示。以压强0 Pa为计算的基准,得到的压强称为绝对压强,通常用表示。这几者的关系如下

(1-11)

(1-12)

在流体力学中,压强都用符号p表示,但对于液体,用相对压强;对于气体,特别是马赫数大于0.1的流动,应视为可压缩流体,用绝对压强。

压强的单位较多,一般用Pa,也可用单位bar,还可以用汞柱、水柱,换算关系如下

1 Pa=1 N/m2

1 bar=100 kPa

1 atm=760 mmHg=10.33 mH2O=101325 Pa

对于静止状态的流体,只有静压强。对于流动状态的流体,有静压强、动压强、测压管压强和总压强之分,可从伯努利方程中分析它们的意义。

伯努利方程阐述一条流线上流体质点的机械能守恒。对于理想流体的不可压缩流动,其表达式如下

(1-13)

式中,为压强水头,也是压能项,为静压强;为速度水头,也是动能项;为位置水头,也是重力势能项。这3项之和就是流体质点的总机械能。H为总的水头高。

将式(1-13)的等号两边同时乘以,则有

(1-14)

式中,为静压强,简称静压;为动压强,简称动压;为总压强,简称总压。对于不考虑重力的流动,总压就是静压和动压之和。

1.1.4 流体运动的描述

1.描述流体运动的方法

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

拉格朗日描述也称随体描述,它着眼于流体质点,并认为流体质点的物理量是随流体质点及时间变化的,即把流体质点的物理量表示为拉格朗日坐标及时间的函数。设拉格朗日坐标为(a, b, c),以此坐标表示的流体质点的物理量,如矢径、速度、压强等,在任一时刻t的值,都可以用abct表示出来。

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

(1-15)

例如,设时刻t流体质点的矢径(即t时刻流体质点的位置)为r,则其拉格朗日描述为

(1-16)

同样,流体质点的速度的拉格朗日描述是

(1-17)

欧拉描述也称空间描述,它着眼于空间点,认为流体的物理量随空间点及时间变化,即把流体物理量表示为欧拉坐标及时间的函数。设欧拉坐标为(q1, q2, q3)[欧拉坐标可以用直角坐标(x, y, z)、柱坐标(r, θ, z)或球坐标(r, θ, ϕ)来表示],用欧拉坐标表示的各空间点上的流体物理量,如速度、压强等,在任一时刻t的值,都可以用q1q2q3t表示出来。由数学分析可知,当某时刻一个物理量在空间的分布一旦确定,该物理量便会在此空间形成一个场。因此,欧拉描述实际上描述了一个个物理量的场。

若以f表示流体质点的某一物理量,其欧拉描述的数学表达是(设空间坐标用直角坐标)

(1-18)

同样,流体质点的速度的欧拉描述是

(1-19)

2. 拉格朗日描述与欧拉描述之间的关系

拉格朗日描述着眼于流体质点,将物理量视为随体坐标与时间变化的函数;欧拉描述着眼于空间点,将物理量视为随空间坐标与时间变化的函数。它们可以描述同一物理量,必定互相相关。设表达式表示流体质点(a, b, c)在t时刻的物理量;表达式表示空间点(x, y, z)于t时刻的同一物理量。如果流体质点(a, b, c)在t时刻恰好运动到空间点(x, y, z)上,则应有

(1-20)

(1-21)

将式(1-20)代入式(1-21)左边,即有

(1-22)

或者反解式(1-20),得到

(1-23)

将式(1-23)代入式(1-21)的右边,也应有

(1-24)

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

3.随体导数

流体质点物理量随时间的变化率称为随体导数,又称物质导数、质点导数。

按拉格朗日描述,物理量f表示为f的随体导数就是跟随质点(a, b, c)的物理量f对时间t的导数。例如速度是矢径对时间的偏导数

(1-25)

即随体导数就是偏导数。

按欧拉描述,物理量f表示为,但并不表示随体导数,它只表示物理量在空间点上的时间变化率。随体导数必须跟随t时刻位于空间点上的流体质点。由于该流体质点是运动的,即xyz是变化的。若用abc表示该流体质点的拉格朗日坐标,则xyz将根据式(1-16)变化,从而f=F(x, y, z, t)的变化满足链式法则。因此,物理量f=F(x, y, z, t)的随体导数是

(1-26)

其中,表示随体导数。

从中可以看出,对于流体质点物理量的随体导数,其欧拉描述与拉格朗日描述大不相同。前者是两者之和,而后者是偏导数。

4. 定常流动与非定常流动

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

定常流动:流体流动过程中各物理量均与时间无关。

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

5. 迹线与流线

常用迹线和流线来描述流体的流动。

迹线:随着时间的变化,空间某一处的流体质点在流动过程中所留下的痕迹。t=0时,位于空间坐标(a, b, c)处的流体质点的迹线方程为

(1-27)

式中,uvw分别为流体质点速度的3个分量,xyzt时刻此流体质点的空间位置。

流线:在同一个时刻,由无数个流体质点组成的一条曲线,曲线上每一点的切线与该质点处流体质点的运动方向平行。流场在t时刻的流线方程为

(1-28)

对于定常流动,流线的形状不随时间变化,而且流体质点的迹线与流线重合。在实际流场中,除驻点或奇点外,流线不能相交,不能突然转折。

6. 流量与净通量

流量:单位时间内流过某一控制面的流体体积,用Q表示,单位为m3/s。若单位时间内流过的流体是以质量衡量的,则称为质量流量,不加说明时“流量”一词概指体积流量。在曲面控制面上,有

(1-29)

在流场中取整个封闭曲面作为控制面A,封闭曲面内的空间称为控制体。流体经一部分控制面流入控制体,同时也有流体经另一部分控制面从控制体中流出。此时用流出的流体减去流入的流体,所得出的流量称为流过全部封闭控制面A的净通量(或净流量),用符号q表示,可通过下式计算

(1-30)

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

7. 无旋流动与有旋流动

由速度分解定理可知,流体质点的运动可以分解为:

1)随同其他质点的平动;

2)自身的旋转运动;

3)自身的变形运动(拉伸变形和剪切变形)。

在流动过程中,若流体质点自身做无旋运动,则称流动是无旋流动,也就是有势流动,否则就称流动是有旋流动。流体质点的旋度是一个矢量,通常用ω表示,其大小为:

(1-31)

ω=0,则称流动为无旋流动,否则就是有旋流动。

ω与流体的流线或迹线形状无关;黏性流动一般为有旋流动;对于无旋流动,伯努利方程适用于流场中的任意两点之间。对于无旋流动,存在一个势函数,满足

(1-32)

(1-33)

8. 层流与湍流

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

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

(1-34)

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

对于圆形管内的流体流动,特征长度L取圆形管的直径d。一般认为临界雷诺数为2320,即

(1-35)

Re<2320时,管中是层流;当Re=2320时,管中流动处于临界状态;当Re>2320时,管中是湍流。

对于异型管道内的流体流动,特征长度L取水力直径dH,则雷诺数的表达式为

(1-36)

异型管道水力直径的定义如下

(1-37)

式中,A为过流断面的面积;S为过流断面上流体与固体接触的周长。

临界雷诺数会因管道形状的不同而有所差别。通过实验得到几种异型管道的临界雷诺数,见表1-1。

对于平板的外部绕流,特征长度取沿流动方向的长度,其临界雷诺数为5×105~ 3×106

表1-1 几种异型管道的临界雷诺数

管道截 面形状

正方形

正三角形

正六边形

矩形

偏心缝隙

a

2070

1930

2190

2260

1000

1.2 CFD基本模型

流体流动所遵循的物理定律,是建立流体运动基本方程组的依据。这些定律主要包括质量守恒定律、动量守恒定律、动量矩守恒定律、能量守恒定律、热力学第二定律。在实际计算时,还要考虑流体的不同流态,如层流与湍流。

1.2.1 基本控制方程

1. 系统、控制体与常用运算符

在流体力学中,系统是指某一确定流体质点集合的总体。系统以外的环境称为外界,分隔系统与外界的界面称为系统的边界。系统通常是研究的对象,外界则用来区别于系统。系统将随系统内质点一起运动,系统内的质点始终包含在系统内,系统边界的形状和所围空间的大小,则可随运动而变化。系统与外界无质量交换,但可以有力的相互作用,及能量(热和功)交换。

控制体是指在流体所在的空间中,以假想或真实流体边界包围,固定不动、形状任意的空间体积。包围这个空间体积的边界面,称为控制面。控制体的形状与大小不变,并相对于某坐标系固定不动。控制体内的流体质点的组成并非不变。控制体既可通过控制面与外界进行质量和能量交换,也可以直接与控制体外的环境有力的相互作用。

本书用到的一些数学符号如下所示

梯度       (1-38)

散度       (1-39)

旋度       (1-40)

式中,

(1-41)

(1-42)

(1-43)

式中,称为拉普拉斯算子。

2. 连续性方程

在流场中,流体通过控制面流入控制体,同时也会通过另一部分控制面流出控制体,在这期间控制体内部的流体质量也会发生变化。按照质量守恒定律,流体流入的质量与流出的质量之差应该等于控制体内部流体质量的增量,由此可导出流体流动连续性方程的积分形式为

(1-44)

式中,V表示控制体,A表示控制面。等式左边第一项表示控制体V内部质量的增量,第二项表示通过控制表面流入控制体的净通量。

根据数学中的奥-高公式,在直角坐标系下可将其化为微分形式

(1-45)

对于不可压缩均质流体,密度为常数,则有

(1-46)

对于圆柱坐标系,其微分形式为

(1-47)

对于不可压缩均质流体,密度为常数,则有

(1-48)

3. 动量方程(运动方程)

动量守恒是流体在运动时应遵循的另一个普遍定律,其描述为:在一给定的流体系统中,其动量的时间变化率等于作用于其上的外力总和。其数学表达式即为动量方程,也称运动方程或N-S方程,其微分形式如下

(1-49)

式中,分别是单位质量流体上的质量力在3个方向上的分量,是流体内应力张量的分量。

运动方程在实际应用中有许多表达形式,其中比较常见的有如下几种。

1)可压缩黏性流体的运动方程

(1-50)

2)常黏性流体的运动方程

(1-51)

3)常密度常黏性流体的运动方程

(1-52)

4)无黏性流体的运动方程(欧拉方程)

(1-53)

5)静力学方程

(1-54)

6)相对运动方程。在非惯性参考系中的相对运动方程是研究大气、海洋及旋转系统中流体的运动所必须考虑的。由理论力学得知,绝对速度为相对速度与牵连速度之和,即

(1-55)

其中为运动系中的平动速度,是转动角速度,为质点矢径。

而绝对加速度为相对加速度、牵连加速度及科氏加速度之和,即

(1-56)

其中

将绝对加速度代入运动方程(1-49),得到流体的相对运动方程

(1-57)

4. 能量方程

将热力学第一定律应用于流体运动,把上式各项用有关的流体物理量表示出来,得到能量方程

(1-58)

式中,是有效热传导系数,,其中是湍流热传导系数,根据所使用的湍流模型来定义;是组分j的扩散流量;包括了化学反应热及其他用户定义的体积热源项。式(1-58)等号右边括号中的3项分别描述了热传导、组分扩散和黏性耗散带来的能量输运。

1.2.2 湍流模型

湍流是自然界广泛存在的流动现象。大气层中气体的流动和海洋环境中水流的流动,飞行器和船舰的绕流,叶轮机械、化学反应器、核反应器中的流体运动都是湍流。湍流流动的核心特征是其在物理上近乎无穷多的尺度和数学上强烈的非线性,这使得人们无论是通过理论分析、实验研究还是计算机模拟都很难彻底认识湍流。20世纪80年代,相关人员提出和发明了一大批高精度、高分辨率的计算格式,相当成功地解决了欧拉方程的数值模拟,可以说欧拉方程数值模拟方法的精度已接近于它有效使用范围的极限;欧拉方程也能适用于各种实践所需。在此基础上,研究人员同时进行了求解可压缩雷诺平均方程及其三维定态黏流流动的模拟。20世纪90年代开启了一个非定常黏流流场模拟的新局面,黏流流场具有高雷诺数、非定常、不稳定、剧烈分离流动的特点,需要继续探求更高精度的计算方法和更实用可靠的网格生成技术,因此研究湍流机理,建立相应的模式,并进行适当的模拟仍是解决湍流问题的重要途径。

1.湍流模型的分类

湍流模型的种类有很多,大致可以归纳为以下3类。

第1类是湍流输运系数模型,即将速度脉动的二阶关联量表示成平均速度梯度与湍流黏性系数的乘积,用笛卡儿张量表示为

(1-59)

模型的任务就是给出计算湍流黏性系数的方法。根据建立模型所需要的微分方程的数目,还可将其细分为零方程模型(代数方程模型)、单方程模型和双方程模型。

第2类是抛弃了湍流输运系数的概念,直接建立湍流应力和其他二阶关联量的输运方程,称为直接数值模拟。

第3类是大涡模拟。前两类模型是以湍流的统计结构为基础,对所有涡旋进行统计平均。大涡模拟把湍流分成大尺度湍流和小尺度湍流,通过求解三维修正的N-S方程,得到大涡旋的运动特性,而对小涡旋运动还采用第1类湍流模型。

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

Fluent 提供的湍流模型包括:单方程模型(Spalart-Allmaras模型)、双方程模型[标准(Standard)模型、重整化群(Renormalization Group,RNG)模型、可实现的(Realizable)模型]、雷诺应力(Reynolds Stress)模型和大涡模拟(Large Fddy Simulation,LES)。

2. Fluent中的湍流模型选择策略

选取湍流模型需要考虑的因素有:流体是否可压、精度的要求、计算机的能力、时间的限制等。Fluent中湍流模型详解如图1-1所示,各湍流模型的应用范围及特点说明如下。

图1-1 湍流模型详解

(1)Spalart-Allmaras模型

Spalart-Allmaras模型由Spalart-Allmaras提出,用来解决因湍流动黏滞率而修改的数量方程,主要是壁面约束的流动,且已经显示出很好的效果。

应用范围:多用于航空领域、透平机械,对于有壁面约束的空气流动问题应优先选用此模型。

该模型有以下特点。①Spalart-Allmaras模型是相对简单的单方程模型,只需求解湍流黏性的输运方程,不需要求解剪切层厚度的长度尺度;由于没有考虑长度尺度的变化,对一些流动尺度变换比较大的流动问题不太适合。例如平板射流问题,从有壁面影响的流动突然变化到自由剪切流,流场尺度变化明显。②该模型的输运变量在近壁处的梯度要比模型中的小,表现出对网格粗糙带来的数值误差不太敏感。③该模型不能预测均匀衰退、各向同性湍流等复杂的工程流动问题。

(2)模型

1)标准模型。

在Fluent中,标准模型的系数由经验公式给出,计算比较稳定,自从被Launder和Spalding提出之后,标准模型就变成了工程分析计算中的主要工具。现有研究表明,对于简单的充分发展湍流问题,该模型完全适用,而且比其他湍流模型更快收敛。

应用范围:适合完全湍流(充分发展的湍流)的流动过程模拟,特别是对高雷诺数的湍流有效,包含黏性热、浮力、压缩性等选项,但模拟旋流和绕流时有缺陷。

2)重整化群模型。

重整化群模型来源于严格的统计技术,与标准模型相比有一些改进:①重整化群模型在方程中加了一个条件,有效地改善了精度;②考虑到了湍流漩涡并提高了在这方面的精度;③为湍流普朗特数提供了解析公式,而标准模型是用户提供的常数;④标准模型是一种高雷诺数的模型,重整化群模型提供了一个考虑低雷诺数流动黏性的解析公式,这些公式的作用取决于如何正确地对待近壁区域。

应用范围:除强旋流过程无法精确预测外,其他流动都可以使用此模型,如模拟射流撞击、分离流、二次流和旋流等,能取得较标准模型更高的精度,但中等复杂的湍流受到各向同性涡旋黏度假设限制。

3)可实现的模型。

可实现的模型与标准模型相比,主要有两个不同点:①可实现的模型为湍流黏性增加了一个公式;②可实现的模型为耗散率增加了新的传输方程,该方程是为层流速度波动而作的精确方程。

应用范围:除强旋流过程无法精确预测外,其他流动都可以使用此模型,包括有旋均匀剪切流、自由流(射流和混合层)、腔道流动和边界层流动。

(3)模型

模型是双方程模型中的一种模型。模型的应用范围为:对自由剪切湍流、附着边界湍流、表现出强曲率的流动以及适度分离湍流都有较高的计算精度。

1)标准模型。

标准模型通过两个输运方程求解。对于有界壁面和低雷诺数的可压缩性和剪切流动,能取得较好的模拟效果,尤其是圆柱绕流问题、放射状喷射、混合流动等,包含转捩、自由剪切和压缩性选项。

2)剪切应力传输模型。

剪切应力传输模型合并了来源于方程的交叉扩散,湍流黏度考虑到了湍流剪应力的传波,而且模型常量不同,这些特点使得该模型在近壁区域的自由流、逆压梯度流动、翼型、跨音速激波等方面有着更高的精度。

(4)雷诺应力模型

雷诺应力模型同时考虑了连续性方程、运动方程、输运方程和各向异性湍流剪切应力方程,但计算收敛较难,计算耗时较长。

应用范围:主要应用于复杂3D流动(如弯曲管道、旋转、旋流燃烧、旋风等流动),尤其是强旋流运动,若要考虑雷诺应力各向异性,则必须应用该模型。

(5)大涡模拟

大涡模拟将湍流过程分为大尺度脉动和小尺度(小于惯性区尺度)脉动,把小尺度脉动认为是湍流耗散[作为额外的应力项(亚网格应力)],通过运动方程进行传递,将计算的重点放在大尺度脉动上,可以表现出高精度的湍流细节。

应用范围:在学术界应用广泛,但需要更为精细的网格及较长的计算耗时。

(6)尺度自适应模拟(Scale-Adaptive Simulation,SAS)

尺度自适应模拟的计算精度和大涡模拟接近,但其计算量明显更少,网格精细程度也更低。尺度自适应模拟的优点是模型的雷诺平均N-S方程(RANS)部分不受网格尺寸的影响,因此不会出现精度下降的情况。

应用范围:通常适用于分离区域、航空航天等问题的分析。

(7)分离涡模拟(Detached Eddy Simulation,DES)

分离涡模拟通过比较湍流的长度与网格间距的大小来实现雷诺平均N-S方程与大涡模拟之间的切换。该模型选择两者的最小值,从而在雷诺平均N-S方程和大涡模拟之间切换。它通过在近壁面采用大涡模拟、在湍流核心区采用雷诺平均N-S方程进行计算,可以阻止模型应力损耗及网格导致的分离。

应用范围:适用于外流空气动力学、气动声学、壁面湍流等问题的计算。

另外还有转捩模型,用于模拟流动从层流转变到湍流的过程、转捩不连续等计算情况。

壁面对湍流流动影响较大,壁面的不光滑及其他因素均会对流动产生影响,在离壁面很近的地方,黏性力将抑制流体切线方向的速度变化,同时壁面也阻碍正常的波动。重整化群模型在近壁面的外部区域,湍流动能受平均流速的影响而增大,湍流运动加剧。而模型、雷诺应力模型、大涡模拟仅适用于湍流核心区域(远离壁面处),故需要考虑如何使这些模型能够用于壁面边界层处的流动。Fluent中的湍流模型提供了3种近壁处理方法:标准壁面函数法、非平衡壁面函数法和加强壁面函数法。具体选择可以参考下面的说明。

标准壁面函数法:应用广泛,计算精度与经济性合理,但属于高雷诺数的经验公式,不适用于低雷诺数区的流动,且没有考虑质量蒸发、和体积力的影响,三维计算也较差。

非平衡壁面函数法:考虑了和非平衡方程及分离、再附着、撞击现象;但在低雷诺数区的流动、伴随质量蒸发的流动、过大的流动、强体积力流动及三维流动模拟方面存在不足。

加强壁面函数法:不依赖于经验壁面法则,应用于复杂流动及低雷诺数区的流动,但需要精细的网格、较大的CPU与内存。

1.2.3 初始条件和边界条件

在计算过程中,初始条件和边界条件的正确设置是关键的一步。Fluent提供了现成的各种类型的边界条件。

1.初始条件

初始条件是计算初始给定的参数,即时,给出各未知量的函数分布后,用户根据实际情况设置的条件。当流体运动定常时,无初始条件。

2.边界条件

边界条件是流体力学方程组在求解域的边界上,流体物理量应满足的条件。例如,流体被固壁所限,流体就不应有穿过固壁的速度分量;在水面边界上,大气压强可认为是常数(一般在不大的范围内可如此认为);在流体与外界无热传导的边界上,流体与边界之间无温差等。由于具体问题不同,边界条件一般要保持恰当:①保持在物理上是正确的;②要在数学上不多不少,刚好能用来确定积分、微分方程中的积分常数,而不是矛盾的或有随意性的。

Fluent中的基本边界条件包括以下几种。

(1)入口边界条件

入口边界条件就是指定的入口处流动变量的值。常见的入口边界条件有速度入口边界条件、压力入口边界条件和质量流量入口边界条件。

1)速度入口边界条件:用于定义流动速度和流动入口的流动属性相关的标量。该边界条件适用于不可压缩流动,用于可压缩流动会得到非物理结果,因为它允许驻点条件浮动。应注意不要让速度入口靠近圆形障碍物,因为这会导致流动入口驻点属性的非一致性太高。

2)压力入口边界条件:用于定义流动入口的压力及其他标量属性。它既可以用于可压缩流动,也可以用于不可压缩流动。压力入口边界条件可用于压力已知但是流动速度或速率未知的情况。这一情况可用于很多实际问题,例如浮力驱动的流动。压力入口边界条件也可用来定义外部或无约束流体的自由边界。

3)质量流量入口边界条件:用于已知入口质量流量的可压缩流动。在不可压缩流动中不必指定入口的质量流量,因为密度为常数时,速度入口边界条件就确定了质量流量入口边界条件。当要求达到的是质量和能量流速而不是流入的总压时,通常就会使用质量流量入口边界条件。

注意: 调节入口总压可能会导致解的收敛速度较慢,当压力入口边界条件和质量流量入口边界条件都可以接受时,应优先选择压力入口边界条件。

(2)出口边界条件

1)压力出口边界条件:需要在出口边界处指定表压。指定的表压只用于亚声速流动。如果流动变为超声速,就不再指定表压,此时压力要从内部流动中求出,其他的流动属性也应如此处理。

在求解过程中,如果压力出口边界处的流动是反向的,回流条件也需要指定。如果指定了比较符合实际的回流值,计算会比较容易收敛。

2)质量出口边界条件:当流动出口的速度和压力在解决流动问题之前未知时,可以使用质量出口边界条件来模拟流动。需要注意的是,如果模拟可压缩流体或者包含压力出口的流体,则不能使用质量出口边界条件。

(3)固体壁面边界条件

对于黏性流动问题,可为壁面设置无滑移边界条件;也可以指定壁面切向速度分量 (壁面平移或者旋转时),给出壁面切应力,从而模拟壁面滑移;还可以根据流动情况计算壁面切应力和与流体换热情况。固体壁面热边界条件包括固定热通量、固定温度、对流换热系数、外部辐射换热、对流换热等。

(4)对称边界条件

对称边界条件应用于计算的物理区域对称的情况。在对称轴或者对称平面上,没有对流通量,所以垂直于对称轴或者对称平面的速度分量为0。因此在对称边界上,垂直边界的速度分量为0,任何物理量的梯度都为0。

(5)周期性边界条件

如果流动的几何边界、流动和换热是周期性重复的,则可以采用周期性边界条件。

边界类型的改变是有一定限制的,不能随意进行修改。边界类型可以分成四大类(面边界、双面边界、周期性边界、单元边界),所有边界类型都可以被划分到其中一个大类中。边界类型的改变只能在大类中进行,而分属不同大类的边界类型不能相互替换。

1.3 有限体积法

在有限体积法中将所计算的区域划分成一系列控制体积,每个控制体积都有一个节点作为代表。将守恒型的控制方程对控制体积做积分来导出离散方程。在导出过程中,需要对界面上的被求函数本身及其一阶导数的构成作出假定,这种构成的方式就是有限体积法中的离散格式。用有限体积法导出的离散方程可以保证具有守恒特性,而且离散方程系数的物理意义明确,是目前流动与传热问题的数值计算中应用最广的一种方法。

Phoenics是最早投入市场的有限体积法软件。Fluent、STAR-CD、OpenFOAM、XFlow和CFX等都是常用的有限体积法软件。它们在流动、传热传质、燃烧和辐射等方面应用广泛。

1.3.1 计算区域离散化与控制方程离散化

有限体积法是一种分块近似的计算方法,其中比较重要的步骤是计算区域的离散化和控制方程的离散化。

1.计算区域的离散化

计算区域的离散化就是用一组有限个离散的点来代替原来的连续空间。一般的实施过程是:把所计算的区域划分成许多个互不重叠的子区域,确定每个子区域中的节点位置及该节点所代表的控制体积。区域离散后,得到以下4种几何要素。

节点:需要求解的未知物理量的几何位置。

控制体积:应用控制方程或守恒定律的最小几何单位。

界面:定义了与各节点相对应的控制体积的界面位置。

网格线:连接相邻两个节点形成的曲线簇。

一般把节点看成是控制体积的代表。在离散过程中,将一个控制体积上的物理量定义并存储在相应节点处。图1-2所示为一维的有限体积法网格,图1-3所示为二维的有限体积法网格。

图1-2 一维的有限体积法网格

图1-3 二维的有限体积法网格

用于计算区域离散化的网格有两类:结构化网格和非结构化网格。

结构化网格的节点排列有序,当给出了一个节点的编号后,可以立即得出其相邻节点的编号,所有内部节点周围的网格数目相同。结构化网格具有实现容易、生成速度快、网格质量好、数据结构简单等优点,但不能实现复杂边界区域的离散化。

非结构化网格的内部节点以一种不规则的方式布置在流场中,各节点周围的网格数目不相同。这种网格虽然生成过程比较复杂,但却有极高的适应性,对复杂边界的流场计算问题特别有效。

2.控制方程的离散化

前面给出的关于流体流动问题的控制方程,无论是连续性方程、运动方程,还是能量方程,都可写成

(1-60)

式中,div表示散度,计算方法如式(1-39)所示,grad表示梯度,计算方法如式(1-38)所示。

对于一维稳态问题,其控制方程为

(1-61)

式中从左到右各项分别为:对流项、扩散项和源项。方程中的是广义变量,可以为速度、温度或浓度等一些待求的物理量。是相应的广义扩散系数,是广义源项。变量在端点AB的边界值为已知。

有限体积法的关键一步是在控制体积上建立积分控制方程,以在控制体积节点上产生离散的方程。通过一维模型方程(1-61),在图1-2所示的控制体积P上积分,可得

(1-62)

式中,是控制体积的体积值。当控制体积很小时,可以表示为A是控制体积界面的面积。从而有

(1-63)

从上式看到,对流项和扩散项均已转化为控制体积界面上的值。有限体积法最显著的特点之一就是离散方程中具有明确的物理插值,即界面的物理量要通过插值的方式由节点的物理量来表示。

为了建立所需形式的离散方程,需要表示出式(1-63)中的界面e和界面w处的。有限体积法中规定,等物理量均是在节点处定义和计算的。因此,为了计算界面上的这些物理参数(包括其导数),需要一个物理参数在节点间的近似分布。可以想象,线性近似是可以用来计算界面物性值的最直接,也是最简单的方式,这种分布叫作中心差分。如果网格是均匀的,则单个物理参数(以扩散系数为例)的线性插值结果为

(1-64)

的线性插值结果为

(1-65)

与梯度项相关的扩散通量的线性插值结果为

(1-66)

对于源项S,它通常是随时间和物理量变化的函数。为了简化处理,可以将S转化为如下线性方式

(1-67)

式中,是常数,是随时间和物理量变化的项。将式(1-64)~式(1-67)代入方程(1-63),可得

(1-68)

整理后得

记为

(1-69)

式中,

(1-70)

对于一维问题,控制体积界面e和界面w处的面积均为1,即单位面积,因此,式(1-70)中各系数可转化为

(1-71)

方程(1-69)即为方程(1-61)的离散形式,每个节点上都可建立此离散方程,通过求解方程组,就可得到各物理量在相应节点处的值。

为了后续讨论更方便,需定义两个新的物理量FD,其中,F表示通过界面上单位面积的对流质量通量,简称对流质量流量,D表示界面的扩散传导性。定义表达式如下

(1-72)

这样,FD在控制界面上的值分别为

(1-73)

在此基础上,定义一维单元的佩克莱数(Peclet number)

(1-74)

表示对流与扩散的强度之比。当为0时,对流-扩散演变为纯扩散问题,即流场中没有流动,只有扩散;当>0时,流体沿x方向流动,当很大时,对流-扩散问题演变为纯对流问题。一般在中心差分格式中,有< 2的要求。

将式(1-72)、(1-73)代入方程(1-71),可得

(1-75)

瞬态问题与稳态问题相似,主要是瞬态项的离散。其一维瞬态问题的通用控制方程为

(1-76)

该方程是一个包含瞬态及源项的对流-扩散方程。从左到右,方程中的各项分别是:瞬态项、对流项、扩散项及源项。方程中,是广义变量,如速度分量、温度、浓度等;为相应的广义扩散系数;为广义源项。

对于瞬态问题的有限体积法求解,在将控制方程对控制体积做空间积分的同时,还必须对时间间隔做时间积分。对控制体积所做的空间积分与稳态问题相同,这里仅叙述对时间的积分。

将式(1-76)在一维计算网格上对时间及控制体积进行积分,可得

(1-77)

改写后,为

(1-78)

式中,A是图1-2中控制体积P在界面处的面积。

在处理瞬态项时,假定物理量在整个控制体积P上均具有节点值,并用线性插值来表示,源项也分解为线性方程,对流项和扩散项的值按中心差分格式通过节点处的值来表示,则

(1-79)

假定变量对时间的积分为

(1-80)

式中,上标0代表t时刻;是时刻的值;f为0与1之间的加权因子,当f=0时,变量取旧值进行时间积分,当f=1时,变量采用新值进行时间积分。

采用类似式(1-80)进行时间积分,式(1-79)可写成

(1-81)

整理后得

(1-82)

同样引入稳态中关于符号FD的定义,并将原来的定义做一定扩展,即乘以面积A,可得

(1-83)

将式(1-83)代入方程(1-82),得

(1-84)

同样,也像稳态问题一样引入apawaE,上式变为

(1-85)

式中,

(1-86)

根据f的取值,瞬态问题对时间的积分有以下几种方案:当f=0时,变量的初值出现在式(1-85)的右边,从而可直接求出在现时刻的未知变量值,这种方案称为显式时间积分方案;当0<f<1时,有现时刻的未知变量出现在式(1-85)的两边,需要解若干个方程所组成的方程组才能求出现时刻的变量值,这种方案称为隐式时间积分方案;当f=1时,为全隐式时间积分方案。当f=1/2时,称为Crank-Nicolson时间积分方案。

进一步将一维问题扩展为二维与三维问题。二维问题的计算网格及控制体积如图1-3所示。只增加第二坐标y,控制体积增加的上、下界面,分别用ns表示,相应的两个邻点记为NS。在全隐式时间积分方案下的二维瞬态对流-扩散问题的离散方程为

(1-87)

(1-88)

三维问题的计算网格及控制体积(两个方向的投影)如图1-4所示。在二维的基础上增加第三坐标z,控制体积增加的前、后界面,分别用tb表示,相应的两个邻点记为TB。在全隐式时间积分方案下的三维瞬态对流-扩散问题的离散方程为

(1-89)

(1-90)

图1-4 三维问题的计算网格及控制体积(两个方向的投影)

1.3.2 常用的离散格式与建立离散方程应遵循的原则

1. 常用的离散格式

有限体积法常用的离散格式有:中心差分格式、一阶迎风格式、混合格式、指数格式、乘方格式、二阶迎风格式、QUICK格式。各种离散格式对一维、稳态、无源项的对流-扩散问题的通用控制方程(1-91)均能得到式(1-92)的形式,其高阶情况如式(1-93)所示。

(1-91)

(1-92)

(1-93)

式中,若为一阶迎风格式,;若为二阶迎风格式,。其中,系数(高阶还有)取决于所使用的离散格式,为了便于比较和编程计算,将不同离散格式下的系数aWaE的计算公式总结于表1-2中。

表1-2 不同离散格式下的系数的计算公式

离散格式

系数

系数

中心差分格式

一阶迎风格式

混合格式

指数格式

乘方格式

二阶迎风格式

QUICK格式

2. 建立离散方程应遵循的原则

在利用有限体积法建立离散方程时,必须遵守以下4个基本原则。

(1)控制体积界面上的连续性原则

当一个面为相邻的两个控制体积所共有时,在这两个控制体积的离散方程中,通过该界面的通量(包括热通量、质量通量、动量通量)的表达式必须相同。例如,通过某特定界面从一个控制体积所流出的热通量必须等于通过该界面进入相邻控制体积的热通量,否则,总体平衡就得不到满足。

(2)正系数原则

中心节点系数和相邻节点系数必须恒为正值。该原则是求得合理解的重要保证。若违背这一原则,得到的结果往往是不真实的解。例如,如果相邻节点的系数为负值,就可能出现边界温度的增加引起的相邻节点温度降低的情况。

(3)源项的负斜率线性化原则

源项的斜率为负可以保证遵循正系数原则。从式(1-71)中看到,即使相邻节点的系数皆为正值,但只要有源项的存在,中心节点系数仍有可能为负。规定可以保证为正值。

(4)系数等于相邻节点系数之和原则

当源项为0时,可以发现中心节点系数等于相邻节点系数之和,而当有源项存在时也应该保证遵循这一原则,如果不能满足这个条件,可以取为0。

1.3.3 离散方程的解法

建立了控制方程的离散方程即建立了可以进行数值计算的代数方程组。常规解法只能应付已知速度场求温度场分布这类简单的问题,所以需要对生成的离散方程进行调整,并且对各未知量(速度、压力、温度等)的求解顺序及方式进行特殊处理。为此,流场数值计算是面对常规解法存在的主要问题进行改善所形成的一系列方法集。流场计算的基本过程是在空间上用有限体积法(或其他类似方法)将计算区域离散成许多小的体积单元,在每个体积单元上对离散后的控制方程组进行求解。其本质是对离散方程进行求解,一般可以分为分离解法和耦合解法两大类,其各自又根据实际情况扩展出具体的计算方法,如图1-5所示。

图1-5 流场数值计算方法的分类

1. 分离解法

分离解法不直接求解联立的方程组,而是顺序地、逐个地求解各变量代数方程组。分离解法中应用广泛的是压力修正法,其求解基本过程如下。

1)假定初始压力场。

2)利用压力场求解运动方程,得到速度场。

3)利用速度场求解连续性方程,使压力场得到修正。

4)根据需要,求解湍流方程及其他标量方程。

5)判断当前时间步上的计算是否收敛。若不收敛,返回第2步,迭代计算;若收敛,重复上述步骤,计算下一个时间步上的物理量。

在使用分离求解器时,通常可以选择3种压强与速度的关联形式,即 SIMPLE、SIMPLEC 和 PISO。SIMPLE 和 SIMPLEC通常用于定常流动计算,PISO用于非定常流动计算,但是在网格畸变很大时也可以使用PISO格式。

2. 耦合解法

耦合解法可同时求解离散方程组,联立求解出各变量,其求解过程如下。

1)假定初始压力和速度等变量,确定离散方程的系数及常数项等。

2)联立求解连续性方程、运动方程、能量方程。

3)求解湍流方程及其他标量方程。

4)判断当前时间步上的计算是否收敛。若不收敛,返回第2步,迭代计算;若收敛,重复上述步骤,计算下一个时间步上的物理量。

Fluent采用的离散格式包括一阶迎风格式、指数律格式、二阶迎风格式、QUICK格式、中心差分格式等。Fluent中各流场变量的迭代都由松弛因子控制,因此计算的稳定性与松弛因子紧密相关。在大多数情况下,可以不对松弛因子的默认设置进行修改,因为这些默认值是根据各种算法的特点得出的。而对某些复杂流动,默认设置可能不能满足稳定性要求,计算过程中可能出现振荡、发散等情况,此时需要适当减小松弛因子的值,以保证计算收敛。

在实际计算中,可以先用默认设置进行计算,如果发现残差曲线向上发展,则中断计算,适当调整松弛因子后再继续计算。在修改计算控制参数前,应该先保存当前计算结果,调整参数后,计算需要经过几步调整才能适应新的参数。

在计算发散时,可以考虑将压强、动量、湍流动能和湍流耗散率的松弛因子的默认值分别减小为0.2、0.5、0.5、0.5。在计算格式为SIMPLEC 时,通常没有必要减小松弛因子的值。

1.4 多相流

多相流是在流体力学、传热传质学、物理化学、燃烧学、微流动、气固运输等学科的基础上发展起来的一门新兴学科,它广泛应用于能源、动力、核能、石油、化工、冶金、制冷、运输、环境保护及航天技术等许多领域,对国民经济的发展有十分重要的作用。

1.4.1 VOF模型的概述和局限

VOF模型通过求解单独的运动方程和处理穿过区域的每一流体的体积分数来模拟两种或多种不能混合的流体。其典型应用包括射流破碎的预测、流体中大泡的运动模拟、决堤后水的流动模拟和气液界面的稳态与瞬态处理。

在Fluent中使用VOF模型存在以下一些限制。

1)必须使用基于压力的求解器,VOF 模型不能用于基于密度的求解器。

2)所有的控制体积必须充满单一流体相或者混合相;VOF模型不允许在空的区域(其中没有任何类型的流体)中存在。

3)只有一相是可压缩的。

4)周期流动(比质量流率或比压降)问题不能和VOF模型同时计算。

5)二阶隐式的时步(Time-Stepping)公式不能用于VOF模型的显示算法。

6)当并行跟踪粒子时,若共享内存启用,VOF模型不能与离散相模型(Discrete Particle Model,DPM)同时使用。

关于稳态和瞬态的 VOF计算,有如下需要注意的问题。

在Fluent中,VOF公式通常用于计算时间依赖解,但是对于只关心稳态解的问题,它也可以执行稳态计算,例如当解是独立于初始时间并且有明显的单相流入边界时,可用稳态VOF计算;又如,渠道内顶部有空气的水的流动和分离的空气入口可以采用稳定状态(Steady-State)公式求解;而在旋转的杯子中自由表面的形状依赖于流体的初始水平,这样的问题必须使用时间依靠公式。

1.4.2 VOF模型的控制方程

VOF模型中的两种或多种流体(或相)没有互相穿插,模型中每增加一个相,就增加了一个变量,该变量为相的体积分数。在每个控制体积内,所有相的体积分数的和为1。所有变量及其属性的区域被各相共享并且代表了体积平均值,且每一相的体积分数在每一位置是可知的。在任何给定单元内的变量及其属性可能纯粹是其中一相的变量及其属性或者多相混合的变量及其属性,具体取值取决于各相的体积分数。即在单元中,如果第相流体的体积分数为,那么就会出现下面3个可能的情况。

1):第相流体在单元中是空的。

2):第相流体在单元中是充满的。

3):单元中包含了第相流体和一相或其他多相流体的界面。

基于的局部值,属性和变量会被适当地分配给每一控制体积。

1. 体积分数方程

在VOF模型中,跟踪相与相之间的界面是通过求解相的体积分数的连续性方程来完成的。第相的体积分数的连续性方程为

(1-94)

式中,是第相到第相的质量传递,是第相到第相的质量传递。默认情况下,式(1-94)右边的源项为零,但可根据情况为每一相设置常数或自定义质量源项。体积分数方程不求解初始相,而是给出各相体积分数的约束

(1-95)

体积分数方程可通过欧拉显示方案或隐式方案求解。

(1)欧拉显式方案

在欧拉显式方案中,Fluent的标准有限差分插值方案被用于对前一时间步的体积分数进行计算

(1-96)

式中,n+1为新时间步的指针;n为前一时间步的指针;为一阶迎风、二阶迎风、QUICK算法中的第相体积分数的数值;V为单元的体积;为以法向速度通过面的体积流量。

这个公式在每一时间步上都不需要输送方程的迭代解,但在隐式方案中是需要的。

使用欧拉显式方案时,时间依赖解必须通过计算。

(2)隐式方案

在隐式方案中,Fluent的标准有限差分插值方案用于获得所有单元的面通量,包括界面附近的

(1-97)

由于这个方程需要当前时间步的体积分数值(而不是上一时间步,如欧拉显式方案),因此,在每一时间步内标准的标量输送方程为每一个第二相的体积分数迭代性地求解。

隐式方案可用于时间依赖和稳态的计算。

说明:Fluent中还有其他界面处理方法,如CICSAM、BGM等,详见Fluent理论书册。

2.材料属性

出现在输运方程中的属性是由存在于每一控制体积中的分相决定的。例如,在两相流系统中,相用下标1和2表示,如果第二相的体积分数被跟踪,那么每一单元中的密度为

(1-98)

通常,求解n相系统的体积分数的平均密度采用如下形式

(1-99)

所有的其他属性(例如黏度)都以这种方式计算。

3. 运动方程

通过求解整个区域内的单一的运动方程得到的速度场为各相共享,运动方程通过材料属性受到所有相的体积分数的影响。

(1-100)

速度场在共享区域是近似的,它的局限是在各相之间存在大的速度差异的情形下,靠近界面的速度的计算精确性会受到不利影响。

4. 能量方程

能量方程也是各相共享的,表示如下

(1-101)

VOF模型处理的能量E和温度T,采用质量平均变量

(1-102)

这里对每一相的均基于对应相的比热容和共享温度。

属性(有效热传导)是被各相共享的。源项包含辐射的贡献,也包含其他体积热源。和速度场一样,在相间存在大的温度差时,靠近界面的温度的精确度会受到限制。在属性有几个数量级的变化时,这样的问题还会增大。例如,如果一个模型包括液态金属和空气,材料的导热性有4个数量级的差异,如此大的差异使得方程需要设置各向异性的系数,以达到收敛和所需的精度范围。

另外还有附加的标量方程、界面附近的插值、表面张力和壁面黏附等因素,详情请查看Fluent理论手册。

1.4.3 混合模型的概述和局限

混合模型是一种简化的多相流模型,它用于模拟各相有不同速度的多相流,但是假定了在短空间尺度上局部的平衡,相之间的耦合是很强的。它也用于模拟有强烈耦合的各向同性多相流和各相以相同速度运动的多相流。

混合模型可以通过求解混合相的连续性方程、动量方程和能量方程,第二相的体积分数方程,以及相对速度的代数式来模拟n相(流体或粒子)。典型的应用包括沉降、旋风分离器、低载荷粒子流,以及气相体积分数很低的泡状流。

混合模型在以下几种情形中可很好地代替欧拉模型:①当存在大范围的颗粒相分布或界面的规律未知或者它们的可靠性有疑问时,完善的多相流模型是不易实现的;②当求解变量的个数小于完善的多相流模型时,像混合模型这样简单的模型能和完善的多相流模型一样取得好的结果。

在Fluent中使用混合模型存在以下一些限制。

1)必须使用基于压力的求解器。

2)只有一相是可压缩的。

3)周期流动(比质量流率)问题不能和混合模型同时计算。

4)凝固与熔化模型不能与混合模型同时使用。

5)湍流大涡模拟不能使用在混合模型及Singhal空穴模型中。

6)多参考系(Multiple Reference Frame,MRF)不能用在混合模型中。

7)混合模型不能用于无黏性流体。

8)并行示踪粒子时,如果启用了共享内存选项离散相模型不能与VOF模型一起使用。

1.4.4 混合模型的控制方程

与VOF模型一样,混合模型使用单流体方法,它有两方面不同于VOF模型。

1)混合模型允许相之间互相贯穿。所以控制体积的体积分数可以是0和1之间的任意值,具体取决于相q和相p所占有的空间。

2)混合模型使用了滑流速度的概念,允许相以不同的速度运动。(注:相也可以假定以相同的速度运动,混合模型就简化为均匀多相流模型。)

混合模型可求解混合相的连续性方程、混合的动量方程、混合的能量方程、第二相的体积分数方程,还有相对速度的代数表达式(如果相以不同的速度运动)。

1.混合模型的连续性方程

混合模型的连续方程为

(1-103)

式中,是质量平均速度,为

(1-104)

是混合密度,为

(1-105)

式中,是第相的体积分数。描述了由于气穴或用户定义的质量源的质量传递。

2. 混合模型的运动方程

混合模型的运动方程可以通过对所有相各自的运动方程求和来获得。它可表示为

(1-106)

式中,是相数,是体积力,是混合黏性,且

(1-107)

是第二相的漂移速度

(1-108)

3.混合模型的能量方程

混合模型的能量方程采用如下形式

(1-109)

是有效热传导率。等式右边的第一项代表了由传导造成的能量传递。包含了所有的体积热源。

若式(1-109)中的Ek

(1-110)

则是对可压缩相的;而是对不可压缩相的,是第相的显焓。

另外还有相对(滑流)速度和漂移速度等方程,请查看Fluent理论手册。

1.4.5 欧拉模型的概述及局限

采用欧拉模型,第二相的数量仅因为内存要求和收敛行为而受到限制,只要有足够的内存,任何数量的第二相都可以模拟。然而,对于复杂的多相流流动,解由于收敛性而受到限制。

Fluent的欧拉多相流模型对液—液和液—固多相流动不区分。颗粒流是其中一种简单的流动,计算中至少有一相被指定为颗粒相。

Fluent 18.0以后的欧拉模型的解基于以下假设。

1)单一的压力是被各相共享的。

2)动量和连续性方程是对每一相求解。

下面的参数对颗粒相是有效的。

1)颗粒温度(固体波动的能量)是对每一固体相进行计算的。这是基于代数关系的。

2)固体相的剪切和可视黏性是把分子运动论用于颗粒流而获得的。摩擦黏性也是有效的。

3)几个相间的曳力系数函数是有效的,它们适合于不同类型的多相流系(也可以通过用户定义函数修改相间的曳力系数,详细描述见Fluent UDF手册)。

4)所有的湍流模型都是有效的,可以用于所有相或者混合相。

在Fluent中,欧拉多相流模型受到以下限制。

1)雷诺应力模型不可用。

2)颗粒跟踪(使用Lagrangian分散相模型)仅与主相相互作用。

3)周期流动(比质量流率)问题不能和欧拉模型同时计算。

4)无黏流是不允许的。

5)熔化与凝固是不允许的。

6)若共享内存选项开启,并行示踪粒子时,不能用离散相模型。

1.4.6 欧拉模型的控制方程

要从单相模型转变为多相模型,需经过一系列步骤实现,首先要设置混合求解及多相求解。由于多相问题是密切联系的,最好先用初始守恒参数集(一阶时间和空间)直接求解多相问题。同时还需知道多相的体积分数及相间动量、热量和质量交换的机制。

1.体积分数

互相贯穿连续的多相流动的描述组成了相位体积分数(Volume Fraction)的概念。体积分数代表了每相所占据的空间,表示为,并且每相独自地满足质量和动量守恒定律。守恒方程可以使用全体平均每一相的局部瞬态平衡或者混合理论方法获得。

相的体积定义为

(1-111)

式中       (1-112)

相的有效密度为

(1-113)

式中,相的物理密度。

2. 守恒方程

欧拉模型的守恒方程(Conservation Equation)也是由连续性方程、运动方程、能量方程组成的。

(1)连续性方程

相的连续性方程为

(1-114)

相的速度;表示从第相到第相的质量传递;表示从第相到第相的质量传递;为源项,默认为0,可根据情况为每相设置常数或自定义质量源项。在运动方程和焓方程中也有类似的变量。

的关系为       (1-115)

q等于p时         (1-116)

(2)运动方程

相产生的动量平衡为

(1-117)

式中,是外部体积力;是升力;是虚拟质量力;是相之间的相互作用力;是所有相共享的压力;是相之间的速度,定义为如果(即p相的质量传递到相), 则,如果(即相的质量传递到相),则,如果,则是第相的压力应变张量,即

(1-118)

其中,是第相的剪切黏度和体积黏度。

式(1-117)中必须有合适的表达使相间作用力封闭。这个力依赖于摩擦、压力、内聚力和其他影响,并服从条件,

Fluent使用下面形式的相互作用力

(1-119)

式中是相间动量交换系数。

(3)能量方程

欧拉多相流模型中的能量方程体现为如下的每相的分离焓方程

(1-120)

式中,为第相的显焓;为热量通量;为焓源项;为第相与第相的热交换强度;为相间焓。相间的热交换必须符合局部平衡条件,即

欧拉模型应用范围比较广泛,详情请查阅Fluent理论手册。

1.5 动网格

在Fluent中,动网格模型可以用来模拟由流域边界运动引起流域形状随时间变化的流动情况。这种流动情况既可以是一种指定的运动(例如,指定围绕物体重心随时间变化,以一定的线速度或角速度运动),又可以是一个未确定的运动(例如,这种运动的围绕物体重心的线速度或角速度是由流域中固体上的受力平衡得出的),下一时间步的运动情况是由当前时间步的计算结果确定的。各个时间步的体网格的更新是基于边界条件新的位置由Fluent自动来完成的。

1.5.1 动网格模型

对于通量,在任一控制体V内,其边界是运动的,守恒方程的通式为

(1-121)

式中,是液体的密度;是液体的速度矢量;是动网格的网格变形速度;是扩散系数;是通量的源项代表控制体V的边界。

在式(1-121)中,第一项可以用一阶向后差分形式表示为

(1-122)

式中的nn+1代表当前和下一时间步的数值。第n+1步的体积可以由下式计算得出

(1-123)

式中的dV/dt是控制体的体积导数,为了满足网格守恒定律,控制体的体积导数由下式算出

(1-124)

式中,是控制体面的数目,表示j面的面积向量。每个控制体面的数量积由下式计算得出

(1-125)

式中,在时间步长的基础上,是由控制卷j面所卷出的体积。

通过二阶迎风向后差分公式,式(1-121)中的时间导数能用下式表述

(1-126)

式中,n+ 1、nn-1分别表示连续时间步的数值,n+1表示当前时间步的数值。

在差分格式下,控制体积的时间导数与第一阶方案的计算方式相同,如式(1-124)所示。对于二阶差分格式,每个控制体面上的点积由下式计算得出

(1-127)

式中,是控制卷在当前时间段和前一时间步长上被扫描出的卷。

Fluent理论手册给出了六自由度(Six Degrees of Freedom,6DOF)求解理论,有需要的读者可查阅Fluent理论手册第3章内容。

1.5.2 动网格更新方法

当运动条件定义在边界条件上时,Fluent提供了3种动网格运动的方法来更新变形区域内的体网格:①基于弹性变形的网格调整;②动态的网格层变;③局部网格重构。

1. 基于弹性变形的网格调整

对于三角形或四边形网格的流体区域,可以通过基于弹性变形的网格调整方法,根据边界节点上的已知位移来光滑地调整流域内节点的位置;可以通过基于弹性变形的网格光滑地更新体网格且不改变网格之间的连通性。

在基于弹性变形的网格调整方法中,网格上任意两节点之间的连线被理想化成互相联结的弹簧。边界节点上给定的位移将产生一个与联结到这个节点所有弹簧位移成比例的力,这样边界上节点的位移通过体网格就在流域中传播过去。从平衡角度来看,对于每一个节点,联结到其上的所有弹簧的合力必为零。这个条件可以用以下的迭代方程表示

(1-128)

式中,是节点i的位移;是与节点i相邻的节点的数量;是节点和与之相邻的节点之间的弹性系数,其中,弹性系数可以定义为

(1-129)

由于边界上的位移是已知的,方程通过雅可比(Jacobian)矩阵对流域内部所有节点进行扫掠体求解。在求解过程中,更新后的节点位置可以表示为

(1-130)

图1-6所示为活塞缸网格弹性变形前后的效果图。

图1-6 活塞缸网格弹性变形前后

2.动态的网格层变

在棱柱(六面体或楔形)形网格区域,动态的网格层变方法可以根据与运动的物面邻近的网格层的高度来决定增加或减少网格的层数。在Fluent中,动网格模型可以指定一个理想的高度。邻近边界的网格层(图1-7中的j层)根据j层的单元层高度h来分裂出新的单元层或与邻近的i层合并成一个新层。

图1-7 动态网格层

如果j层中单元体积是处于膨胀状态的,Fluent允许它们膨胀到

(1-131)

式中,是理想单元高度,是全局单元层的分裂或合并因子。

时,单元将根据预定义的高度条件进行分裂,也就是说,在j层中的单元将分裂成一个具有理想高度的单元层和一个单元高度为的单元层。如果j层中的单元体积是被压缩的,当压缩到时,这个被压缩的单元层将与邻近的单元层合并成一个新层;也就是说,在i层和j层中的单元将被合并。

图1-8所示为活塞的单元层的网格层变。

图1-8 活塞的单元层的网格层变

3.局部网格重构

当边界位移相对局部单元尺寸较大时,单元质量将恶化或单元将退化,这样会导致下一时间步的求解收敛困难。为了避免出现这个问题,Fluent把这些网格质量差的单元(单元尺寸太大或太小的单元,或者高度变形的单元)聚集成团,再对这个网格聚团区域进行重构。

Fluent检查流域内每一个单元的质量,并且对满足下列一个或多个条件的单元给予标注。

1)单元尺寸小于指定的单元尺寸最小值。

2)单元尺寸大于指定的单元尺寸最大值。

3)单元的畸变度大于指定的最大畸变度。

除了对体网格进行重构外,Fluent还对变形边界上网格的三角形面或线性面随同体网格一起进行重构,如图1-9所示。只有那些与运动节点相邻的网格面在考虑范围之内。这种网格重构技术类似于前文介绍的动态的网格层变技术,只不过它是应用于与运动网格节点相邻面网格的高度的网格运动。

这些与运动的物面邻近的网格单元面(图1-9中的j层)根据其高度h来确定是分裂出新的单元层还是与邻近的层(图1-9中的i层)合并成一个新层。

图1-9 变形区域的重划分

如果j层中单元面是处于膨胀状态的,则Fluent允许它们膨胀到

(1-132)

式中,是理想单元面高度,是高度系数。

时,单元将根据预定义的高度条件进行分裂,这时,在i层中的单元面的高度将正好是理想高度。相反,如果j层中的单元体积是被压缩的,当压缩到时,这些被压缩的单元面将与邻近层的单元面合并成一个新的单元层。

图1-10所示为乒乓球局部单元重构前后的对比。

(a)0 s

(b)0.01 s

(c)0.015 s

图1-10 乒乓球局部单元重构前后

1.6 滑移网格

在Fluent中,除了动网格能够描述计算区域的瞬态流动外,滑移网格也能描述计算区域的运动,这是在多参考系和混合平面法基础上发展起来的。本节主要介绍滑移网格的基本概念及在Fluent中进行相关计算的设置。

1.6.1 滑移网格的应用、运动方式及网格分界面形状

1. 滑移网格的应用

1)图1-11所示的是相遇的两动车。

图1-11 相遇的动车

2)图1-12所示为轴承的运动,图1-13所示为风扇叶片的运动。

图1-12 轴承的运动

图1-13 风扇叶片的运动

2.滑移网格的运动方式

在两个或更多的单元区域(如果在每个区域独立划分网格,则必须在开始计算前合并网格)应用滑移网格技术,则每个单元区域至少有一个分界面,且该分界面区域和另一单元区域相邻,相邻的单元区域的分界面互相联系形成“网格分界面”,这两个单元区域相对网格分界面移动。

在计算过程中,单元区域沿着网格分界面滑动(如旋转或平移),而两个区域的网格不会发生变化。图1-14和图1-5所示为两个网格在不同时刻的位置:前者是初始时刻两个网格区域的位置,后者为旋转后两个网格区域的位置。

图1-14 轴承初始时刻的网格

图1-15 轴承相对旋转后的网格

3. 网格分界面形状

当两个分界面是基于同样的几何体时,网格分界面和相联系的分界面可以是任何形状的。图1-16所示为线性网格分界面,网格分界面以虚线表示;图1-17所示为圆形网格分界面。

如果把图1-16伸展为三维时网格分界面将是平面矩形,如果把图1-17伸展为三维时网格分界面将是圆柱体。

图1-16 线性网格分界面

图1-17 圆形网格分界面

图1-18所示为圆锥形网格分界面,虚斜线表示在三维平面上的圆锥形网格分界面的交线。

对于轴向的转子、定子构造,其旋转和静止部分沿轴线成一行,而不是同中心的,如图1-19所示,网格分界面是扇形的。该扇形垂直于旋转轴的截断面区域,旋转轴为转子和定子之间的轴线。

图1-18  圆锥形网格分界面

图1-19 扇形网格分界面

1.6.2 滑移网格的原理

1. 滑移网格数据传递原理

滑移网格模型允许相邻网格之间相对滑动,因而网格面不需要在分界面上排列,滑移网格的关键是计算流进每个网格分界面的两个非一致的分界面区域。

界面流动时,每隔一个新的时间步长,在该交界面就会产生内部区域(Interior,在两边都有流体单元的区域)和一个或多个周期区域。如果不是周期性的问题,那么交界面会产生一个内部区域和两个壁面区域(Interface),如果两个交界面完全贴合则没有壁面区域,如图1-20所示,壁面区域需要采用一些适当的边界类型。重叠的分界面区域对应所产生的内部区域;不重叠的交界面区域对应所产生的周期性区域。在这些交界面区域内的面的数目随着分界面的相对移动而不同。理论上,网格分界面的流量应该根据两分界面的交叉处所产生的面来计算,而不是根据它们各自的分界面的面来计算。

图1-20 非周期性交界面相交区域

在图1-21中,分界面区域由A-B面和B-C面及D-E面和E-F面组成,交叉处产生a-d面、d-b面和b-e面等。在两个单元区域重叠处产生d-b面、b-e面和e-c面而组成内部区域,剩余的a-dc-f面成对形成周期性区域。例如,计算从分界面流入Ⅳ单元区域的流量时,用d-b面和b-e面代替D-E面,并从Ⅰ单元区域和Ⅲ单元区域传递信息到Ⅳ单元区域。

图1-21 二维网格交界面

2. 滑移网格控制方程

滑移网格通用控制方程与动网格类似。由于滑移网格的网格运动是刚性的,网格不变形,则式(1-123)可以简化为

(1-133)

式(1-122)可以简化为

(1-134)

式(1-124)可以简化为

(1-135)

在应用了上述简化后的公式之后,即可采用动网格的控制方程式(1-121)计算滑移网格问题。

1.7 ANSYS流固耦合分析

流固耦合问题是流体力学与固体力学交叉生成的一门力学分支,它是研究变形固体在流场作用下的各种行为及固体位形对流场影响这二者相互作用的一门科学。流固耦合力学的重要特征是两相介质之间的相互作用,变形固体在流体载荷作用下会产生变形或运动。变形或运动又反过来影响流体运动,从而改变流体载荷的分布和大小,正是这种相互作用使得在不同条件下产生形形色色的流固耦合现象。所以,近年来流固耦合分析在工程设计中,特别是虚拟设计和数值模拟中的应用越来越广泛和深入。

ANSYS很早便开始进行流固耦合的研究和应用,目前ANSYS中的流固耦合分析算法和功能已相当成熟,在ANSYS Workbench协同仿真平台上便可以进行相关流固耦合的分析计算,ANSYS Workbench 14.0加入的耦合(Coupling)模块,使得用Fluent和Structure进行双向流固耦合分析成为可能。

从算法上讲,ANSYS(也包括其他大型商业软件)主要采用分离解法也就是载荷传递法求解流固耦合问题。但从数据传递的角度出发,流固耦合分析还可以分为两种:单向流固耦合分析(One-way Fluid-Solid Coupling)和双向流固耦合分析(Two-way Fluid-Solid Coupling)。其中,双向流固耦合分析因为求解顺序的不同又可分为顺序求解法和同时求解法,图1-22简单概括了基于ANSYS的流固耦合分析类型。

图1-22 基于ANSYS的流固耦合分析类型

1.7.1 理论基础

(1)流体控制方程

流体控制方程前文中已经详细介绍,此处不再赘述。

(2)固体控制方程

固体部分的守恒方程可以由牛顿第二定律导出

(1-136)

其中,是固体密度;是柯西应力张量;是体积力矢量;是固体域当地加速度矢量。

(1-137)

式1-137是传热方程。其中,表示导热系数;表示能量源项;htot为比焓;T为温升。

对于固体部分,增加了由温差引起的热变形项, 为

(1-138)

其中,是与温度相关的热膨胀系数。

(3)流固耦合方程

同样,流固耦合遵循最基本的守恒原则,所以在流固耦合交界面处,应满足流体与固体应力()、位移()、热流量()、温度()等变量的相等或守恒,即满足如下方程组

(1-139)

(1-140)

(1-141)

(1-142)

其中,下标f表示流体;下标s表示固体。

上述就是流固耦合分析所采用的基本控制方程,为便于分析,可以建立控制方程的通用形式,然后给定各参数及适当条件的初始条件和边界条件,统一求解。当前,用于解决流固耦合问题的方法主要有两种:直接耦合解法和分离解法。直接耦合解法通过把流固控制方程耦合到同一个方程矩阵中求解,也就是在同一求解器中同时求解流体和固体的控制方程。

(1-143)

其中,表示迭代时间步;分别表示流场的系统矩阵、待求解值和外部作用力;分别对应固体区域的各项;代表流固的耦合矩阵。

由于同时求解流固的控制方程,不存在时间滞后问题,所以直接耦合解法在理论上非常完善。但是,在实际应用中,直接耦合解法很难将CFD和计算固体力学(Computational Solid Mechanics,CMS)技术真正结合到一起,同时考虑到同步求解的收敛难度和耗时问题,直接耦合解法目前主要应用于如压电材料模拟等电磁-结构耦合和热-结构耦合等简单问题中,对流动和结构的耦合只能应用于一些非常简单的研究中,还没有在工业应用中发挥重要的实际作用。

流固耦合的分离解法不需要耦合流固控制方程,只需要按设定顺序在同一求解器或不同的求解器中分别求解流体控制方程和固体控制方程,通过流固交界面把流体和固体的计算结果互相交换传递,待此时刻的收敛达到要求,再进行下一时刻的计算,直到求得最终结果。相比于直接耦合解法,分离解法有时间滞后性和耦合界面上的能量不完全守恒的缺点,但是这种方法的优点也显而易见,它能最大化地利用已有的CFD和CSM的方法和程序,只需对它们做少许修改,从而保持程序的模块化;另外分离解法对内存的需求大幅降低,因此可以用来求解实际的大规模问题。所以,目前几乎所有的商业计算机辅助工程(Computer-Aided Engineering,CAE)软件的流固耦合分析都采用的是分离解法。

1.7.2 单向流固耦合分析

单向流固耦合分析在耦合交界面处的数据传递是单向的,一般是指把CFD分析计算的结果(如力、温度和对流载荷)传递给固体结构分析,没有固体结构分析结果传递给流体分析的过程。即只有流体分析对固体结构分析有重大影响,而固体结构分析的变形等结果非常小,以至于对流体分析的影响可以忽略不计。单向流固耦合的现象和分析非常普遍,例如热交换器的热应力分析、阀门在不同开度下的应力分析、塔吊在强风中的静态结构分析、旋转机械的结构强度分析等都属于单向耦合分析,图1-23所示的风力发电机叶片和支架的受力分析便为典型的单向耦合分析。

图1-23 风力发电叶片和支架的受力分析

另外,已知运动轨迹的刚体对流体的影响分析在某种程度上也可以看作是一种单向耦合分析。如汽车通过隧道时对隧道内部气流的影响分析,快启阀在开启过程中对流体流动的瞬间影响等。由于固体运动已知,且固体变形忽略不计,所以此类问题一般可以单独在CFD求解器中完成,但是运动轨迹需要用户通过自定义函数设置。

1.7.3 双向流固耦合分析

双向流固耦合分析是指数据交换是双向的,也就是既有流体分析结果传递给固体结构分析,又有固体结构分析的结果(如位移、速度和加速度)反向传递给流体分析。此类分析多用于流体和固体介质密度相差不大或者高速、高压下固体变形非常明显及其对流体的流动造成显著影响的情况。常见的分析有挡板在水流中的振动分析、血管壁和血液流动的耦合分析、油箱的晃动和振动分析等。一般来讲,大多数耦合作用现象,如果只考虑静态结构性能,采用单向耦合分析便足够,但是如果要考虑振动等动力学特性,双向耦合分析必不可少。也就是说双向耦合分析是为了解决振动和大变形问题而进行的,最典型的例子莫过于深海管道的激振问题。同理,如前所述,塔吊在强风中的静态结构分析属于单向流固耦合分析,但是若要考虑塔吊在强风中的振动情况,就需要采用双向流固耦合进行分析,图1-24所示的薄片在流场中的分析便为典型的双向流固耦合分析。

图1-24 薄片在流场中的分析

1.7.4 耦合面的数据传递

流固耦合中耦合面的数据传递是指将流体计算结果和固体计算结果通过交界面相互传递。不管是完美对应的流固网格还是相差很大的非对应网格(Dissimilar Mesh),通过严格设置,ANSYS多场求解器MFS和MFX都能很好地完成传递。但是对于非对应网格的数据传递,传递前的插值运算必不可少。

多场求解器MFS提供两种插值方式,分别是Profile Preserving插值法和Globally Conservative插值法。在Profile Preserving插值法中,数据接收端的所有节点映射到数据发射端的相应单元上,要传递的参数数据在发射端单元的映射点完成插值后,传递给接收端,是一种主动式传递。与之相反,Globally Conservative插值法先把发射端的节点一一映射到接收端单元上,然后把要传递的参数数据按比例切分到各个节点上,对接收端而言,这属于被动式传递。

两种插值法的出发点和原理不同,所以效果也相差甚远。例如使用Profile Preserving插值法传递参数数据(如力、热通量等)时,发射端和接收端的数据有可能不守恒;使用Globally Conservative插值法在局部也有类似的不守恒情况,但是可以保证全部交界面数据的总体守恒。从物理角度出发,力、热通量等参数在耦合交界面处保持总体守恒更有意义,但是对于位移和温度等参数,保持整体上的守恒不是很有意义,反而局部的分布轮廓更需要精确传递。所以一般情况下,对力、热通量等参数的传递,可以根据网格情况采用Globally Conservative 插值法或者Profile Preserving插值法;对位移和温度的传递,一般采用Profile Preserving插值法。

与MFS相似,多场求解器MFX同样提供两种插值方式,分别是Profile Preserving插值法和Conservative插值法。MFX中的Profile Preserving插值法与MFS中的完全相同。虽然MFX中的Conservative插值法与MFS中的Globally Conservative插值法只有一个单词的区别,但其原理、方法完全不同。MFX中的Conservative插值法采用单元分割、像素概念、桶算法及新建控制面等多种方式和方法完成数据传递,只要确保流固耦合面能重合,交界面上的参数数据从全局到局部都能得到精确传递。对于流固耦合面不完全对应的情况,Conservative插值法会通过在不对应区域设置0值、特殊边界条件来忽略此区域数据的传递,从而保持严格的守恒传递。

1.8 CFD软件计算结果的验证及修正

CFD软件要想在工程中得到广泛的应用,必须要考虑准确性与可信性。目前常用的CFD 软件几乎都采用有限体积法(除了CFX采用混合有限元法与有限体积法外, Fluent、 STAR-CD、 Phoenics、 Flow-3D、OpenFOAM等都采用有限体积法)。在工程实践中为了提高CFD计算结果的精度与可信度,通常要对计算结果进行验证。实验是最好的验证手段,但是普遍存在实验过程中的参数很难与输入的参数完全吻合的情况。一般来说,在工程中,数值计算结果与实验结果误差在10%以内是被允许的。在数值计算结果与实验结果存在很大差异时,一般进行以下检查与分析。

1)几何模型的修正。分析是否忽略了关键几何特征、检查边界位置是否恰当。边界位置选择不当,可能会导致计算振荡,不收敛等情况发生。同时由于不同的软件对于不同的边界组合方式的处理方法存在差异,因此需要根据模拟软件自身的特点选择合适的边界组合方式(如Fluent中压力边界与自由出流边界最好不要同时出现,同时出现可能导致计算不收敛问题。流量入口边界收敛要比压力入口困难)。

2)物理模型的修正。分析是否选用了不合适的模型。每一种模型都有一定的使用范围,模型的选择要基于对适用范围和使用条件的深刻认识。例如在Fluent中,湍流模型有很多,一般的工程流动问题适合采用标准模型,但是其应用在强旋流中误差较大;重整化群 模型常用于旋转流动湍流计算;Spalart-Allmaras模型常用于航空外流计算;模型则适合边界层计算;雷诺应力模型适合各向异性湍流的计算,但是计算量大不易收敛。因此在进行模型的选择时需要仔细地考虑选择的模型是否适合当前问题,一旦模型选择错误,轻则造成大的偏差,重则不收敛、计算出错。

3)检查是否忽略了某些重要的物理现象。在计算结果出现较大偏差时,应考虑是否有些物理过程被忽略或者过于简化。例如复杂几何模型中出现大的负压区,是否需要考虑空化;计算高压气体时,是否考虑可压缩性,是否考虑黏性热;还有一些情况下,是否考虑蒸发、冷凝等相变情况。有时候这些物理现象会导致计算的不收敛乃至计算错误。

4)优化网格。高质量的网格能够增强收敛、提高计算精度、减少计算时间。因此在时间充裕的情况下,需尽可能地提高网格质量;同时,对流动情况复杂的区域进行网格加密处理。在计算结果达到要求后,还需要进行网格独立性验证。

5)边界条件分析。分析测量精度是否满足要求;若边界信息来自计算,那么分析计算时的条件是否满足。有时由于微流体自身的复杂性及实验设计与测量的现实困难,实验数据间存在着相互矛盾的现象,此时进行深入的研究是非常必要的。

6)模拟过程各主要环节的匹配。建立几何模型要考虑所建模型的复杂程度,要尽可能与实际物体一致,又要有简化提高,做到模拟计算与几何模型生成方法合理匹配,几何模型网格划分疏密程度与计算机内存合理匹配,几何模型与边界条件合理匹配,物理模型与模拟软件合理匹配。其中,合理建立几何模型对达到模拟精度要求起决定作用。

相关图书

ANSYS Workbench中文版超级学习手册
ANSYS Workbench中文版超级学习手册
ANSYS Fluent中文版超级学习手册
ANSYS Fluent中文版超级学习手册
ANSYS有限元分析完全自学手册
ANSYS有限元分析完全自学手册
ANSYS Workbench有限元分析实例详解(动力学)
ANSYS Workbench有限元分析实例详解(动力学)
CAE分析大系 ANSYS nCode DesignLife疲劳分析基础与实例教程
CAE分析大系 ANSYS nCode DesignLife疲劳分析基础与实例教程
ANSYS 19.0有限元分析完全自学手册
ANSYS 19.0有限元分析完全自学手册

相关文章

相关课程