计算流体力学基础

1、计算流体力学基础The Base of Computational Fluid Dynamics绪论10.1流体力学的基本方程10.1.1守恒原理10.1.2 质最守恒方程20.1.3动量守恒方程20.1.4 其他标量的守恒方程30.2无因次化方程303流体力学方程的简化模型40.3.1不可压缩流动40.3.2无粘流动40.3.3蠕变流40.3.4白然对流40.3.5边界层流动50.4流动的数学分类5第一章计算流体力学概论71.1 引言71.2 CFD的定义813CFD的优、缺点81.4 CFD的基本耍素91.4.1数学模型91.4.2离散方法91.4.3 坐标和向量系统91.4.5数值的网格

2、91.4.6有限近似101.4.7求解方法101.4.8收敛判据101.5 数值方法的基本特性101.5.1相容性101.5.2稳定性101.5.3收敛性101.5.4守恒性111.5.5.有界性111.5.6真实性111.5.7轿确性111.6 离散方法111.6.1 有限差分法(Fmhe Diffeience Method, FD)点近似111.6.2有限体积元法(Fuute olunie Method, FV) 控制体内的平均近似.121.6.3 有限元法(Finite Eleinenf Method )函数逼近12思考题:12笫2章有限差分法142.1引言142.2基本概念142.3

3、阶导数的近似152.3.1 Taylor级数展开法152.3.2多项式拟合法162.3.3迎风格式162.4二阶导数的近似182.5混合导数的近似182.6边界条件192.7时间项的处理192.7.1初值问题的的基木解法192.7.2多层格式和预估校正格式202.7.3 Riinger Kutta 法212.8有限差分法的应用实例212.9差分格式的稳定性分析242.9.1 Von Neiunann 方法24笫3章有限体积法263.1有限体积法基本原理263.2而积分的近似263.3体积分的近似273.4函数的插值283.4.1 迎风插值（UDS） 283.4.2 线性插值（CDS） 283.

4、4.3三阶迎风格式（QUICK）293.4.4高阶格式（4阶精度CDS） 293.5边界的处理293.6有限体积法应用举例293.6 SIMPLE 方法31第四章边界元法344.1线性方程的血加原理（SupeiposMon） 344.2奇点的速度势3443叠加原理在势流数值计算中的应用举例一回转体绕流354.4 Gieen 定理364.5无界流中无升力问题的数值解法Hess-Snuch方法394.5.1 Hess-Sniith 法的棊本原理394.5.2物体表面单元的局部坐标系及坐标转换404.5.3平而奇点的积分公式414.5.4 Hess-Snuth法应用举例圆球的附加质最414.6无厚度

5、机翼的升力问题421The base of Computational Fluid Dynamics计算流体力学基础教材:J.H. Ferziger, M Penc, Computational Methods for Fluid Dynamics参考资料：苏铭徳，计算流体力学，清华人学出版社周雪漪，计算水力学，清华大学出版社绪论流体是这样一种物质：它的分子结构不能承受任何的外部剪切应力，即使是微小的剪力也会引 起流体微团的变形。在宏观上，流体通常被看作是一种连续的介质(连续统假设，continuum) 流体所受的外力包括：表面力:如压力,剪切力,表面张力等体枳力：如重力，旋转产生的向心力等在

6、很多方而，流体的流动速度决定了流动的特点低速流动时(Re<30丿，流体的惯性力可以忽略，这样的流动称为蠕变流(Creepingflow) 端变流(Creepingflow)包括微小颗粒的绕流(如粉尘的降落)，多孔介质中的流动(过滤过程)以 及狭窄通道的流动(纺织技术，微型机械，微型生物)速度较人时(3(XRe<10$10°),惯性力不能被忽略，但足流体微团的运动轨迹为平滑曲线，这样的流 动称为层流(Laminae)速度很人时，流动的紈迹带有很人的随机性，这样的流动称为湍流(turbulent)流速和音速的比值称为Mach数，Ma数决定了是否耍考息流体动能和内能的交换过程M

7、a<0.3,流体可以看作是不可压缩的(incompressible),否则必须考虑压缩性 0.3<Ma<l,流动为亚咅速流(subsonic)Ma> 1，诒i动为超h速流(supeisomc)，此时会ill现激波(shock uave)。Ma>5,流体被压缩产生的高热会使流体的化学性质发生变化，这样的流动称为崗超音速流 (hypersonic)0.1流体力学的基本方程0.1.1守恒原理守恒律可以从给定控制物质(control mass)的外在特性(extensive properties)导出。如质量动 和能磺的守恒。研究固体力学时，控制物质是很容易确泄的，但是在

8、研究流体力学时，处理给定空 间区域(controlvolume)的流动更为方便。设0为任意内在的守谊吊”相W的外在特性为1=J pga上式中，。纽为控制物质的体积，根据迁移定理(01)dp(jdQ + jp(/>(y - vb) -ndS descv(0 2)上式中为控制体（CV）体枳，Sc为控制体表仏n为控制体表面的外法线方向.v是控制体衷而的流体运动速度，Vb是控制体衷而的运动速度，多数情况卜，vb-0o0-1.2质量守恒方程在方程（0.2）中，収尸1,可得质最守恒方程:f pdG - pdG + f pv - ndS = 0 (，t 盒 dt O；v 盒(0.3)写成微分形式为&a

9、mp;p+ V (pv) = 0(04)0-1.3动量守恒方程在方程（0.2）中，取尸v,可得质駅守恒方程:ddij pvdG = J pvdG + j* pvN ndS =工 /Ocvscr(0 5)作用在控制体上的外力包括：表面力（压力，应力，表面张力等）体积力（重力,科氏力Coriolis forces，电磁力等）从微观的角度來讲，压力和应力來源通过衷而的微观动就交换。 对J牛顿流体（Newtonianfluids）»剪应力张龟/zV-v ll+ 2/zD其中，卩为动力学粘性系数.p为静用力 I为单位张帚，D为变形张吊D = （Vv+ （Vv）7 ）/2写成直角坐标系卜的标量形

10、式为(0.6)(0.7)(08)32 %=一 P + 三"3 dx. J /(09)52 %(0.10)注意，所冇的公式都采用了 Emstrd求和约定，即所仃的卜标如果在一项中出现两次表示对所仃的卜标进行求和，如乞=些+殂+些。dxi dVj dx2 dx3由粘性引起的粘性应力可表示为：2% = 2小厂寸貯r如果用b农示体积力，则动駅方程可写成如I、形式：(0.12)J pvt/Q + j pw ndS = J T ndS + | pbdGQcvScvS(vCcv0丄4其他标量的守恒方程能呈守恒jphdQ + jphv-ndS = JVT ndS + J(v- V/? + S: Vv

11、)rf+ JpdG (0.13) '5$c%QcvQc v其中h为烁T为温度，k热传导系数,S粘性剪应力张量0.2无因次化方程流动的试验研究通常用到模型试验.并把试验结果用无因次的形式表达出来.最终换算到实际的流 动条件。这种手段也町丿UJ：数值分析。分别对时间空间坐标X,速度q,压力p,温度T进行无因次化，=丄，xtf = "7_ »W,=，升厶，二茫勒。/。Lq%p'oT-Tq如果流体的物性参数是常数，则连续性方程，动彊方程，温度方程的无因次形式为：d-二 0dx i(0.13)S且+如J1 dzu. dp1 _(014)dtRe dx; dxi Fr1

12、s匹+材).=1 vT(0.15)dt西Re Pi#各无因次参数依次为S2亘,=Fr=-fi=, Pr = o叽“曲k0.3流体力学方程的简化模型守恒方程是耦合非线性方程组,求解十分困难。在很多情况下,方程中的某些部分筹丁o,或 者是影响很小，可以忽略不计，通过对方程进行简化，可以人人的降低求解的难度。0.3.1不可压缩流动如果流体是不可压缩的，则质駅守恒方程和动届守恒方程可简化为：(016)(017)V v = 0也 + V (mv) = V _ 丄越 + 5dtp dXf0.3.2无粘流动(0 18)(019)(0.20)如果不考虑流体的粘性，则动昴守恒方程可简化为：dtoxi如果流体又是

13、不可压缩的无旋流动，则可进一步简化为势流：V2<D= 0V = 0.3.3蠕变流当流动的Re数很小时，惯性力和非定常力可以忽略不计，则动昼方程可简化为：V 閃冷)-？ + 皿=0(0.22)0.3.4自然对流由热传递过程中温度差形成的流体密度的微小差异也可导致流动的产生，这样的流动成为自然对 流。在处理自然对流时，流体依然是看作不叮斥缩的，釆用Boussmesq假设，认为密:度随温度的变 化是线性的，这时动量方程可写成：i E = V 心)-阳 PE P°E 人)(。R丿0.3.5边界层流动(0 24)込+些迎+沁么變dtdxdx2dv/ dxt0.4流动的数学分类考還含两个f

14、t变吊的线性偏微分方程G'rr + blixy + CU yy + x + eU y + fU = S（0.25）如果系数ag仅是fl变駅xy的函数，则偏微分方程是线性的，如果g0,则偏微分方程（0.25丿是 齐次的。（0.25丿还可以表达为算子的形式：U"） = g对线性的算子L.符介以卜运算规律（superposeon）：L （坷）+ 厶（”2）=厶（坷 +“2）L(anJ = aL(uj , oc为常数L(ai( + fluz)=血(“J + J3L(u2)如果系数aF是x. y以及函数u的两数，则偏微分方程（0 25）是卄线性的.但是如果方程授高阶导 数的系数中不包含

15、最高导数.则称为拟线性方程。例：+ y”“ = X2 :非齐次二阶线性偏微分方程wwxr + wvv = 0:拟线性二阶偏微分方程叱-竹=S111X:II：线性二阶偏微分方程二阶（拟）线件偏微分方程的数学分类是依据最高阶（二阶）导数的系数来划分的。<0椭圆型（elliptic）b2 - 4«c- = 0抛物型（parabolic）> 0双曲型（hvpeibolic）不同类型的偏微分方程的解有不同的特征，在实际数值求解时也碍要不同的数值方法。比如，双曲型方程存在特征线（解的特征传递的方向），在流体力学中，超音速流的激波农而就是特 征线（而）但是椭圆熨方程就没冇特征线，这样.

16、双曲蚁方程和椭圆也方程需耍不同的边界条件。流动的数学分类是依据具体流动的控制方程的数学分类进行的。如势流的控制方程为：v:=o,为椭圆型方程，势流就划分为椭圆型流动；超音速流动的控制方程(l_Ma')20+20=O为双曲型方程，因此划分为双曲型流动: dr Or边界层方程(0.24丿是抛物熨方程，则将边界层流动归入抛物型流动。9第一章计算流体力学概论本章教学目标：介绍计算流体力学(CFD)的一些基本概念，计算流体力学的基本内容和数值计算的主耍特性。 使学生对计算流体力学有一个总体的认识。本章的重点和难点：CFD的定义；CFD的基本构架：数值方法的巫耍特征；常用的CFD方法。学时安排：2

17、学时本章主要的主页外语词汇：计算流体力学：Computational Fluid Dynamics, (CFD)结构化网格：(Structured /Regular Grid)块结构化网格：(Block- strucmred Grid)lh结构化网格：(Unstructured Gnd)Yf 限差分法：Finite DifFerence Method (FDM)限体枳尤法：(Finite Volume Method (FVM)仃限元法：Finite Element Method (FEM)1.1引言尽管流体力学的皋本方程早在一个多世纪就已经发现了，但是奁到Id前为止，仍只能用得到一 些简单流动

18、的解析解。这些优先的结果虽然仃助我们深入认识流动的木质，但是很难将ItfLhV 用到工程实践当中。因此人们不得不采用苴他的手段來解决问题。在很多情况F.人们采用堆近似分析和无因次分析以及经验公式的圾础上建芒的简化方程。例如，阻力的定义Fd = CDS ,S为迎流面积，p为密度，丫为流速，Cd为阻力系数，是一些 无鼠纲数的函数。这一公式是皐无因次分析的结果。通常情况卜，CD可以通过模型实验得到，在 相同无危纲参数的情况卜，可用实际程问题。这样的公式对仅有一两个无杲纲数的复杂形体的 流动是非常有用的(很多情况卜°, Re数是主要的无暈纲数)。但是对/依赖序个无砒纲参数的流动，实验中要同时

19、满足所仃的无鉄纲参数就常的困难了。 比如在空气动力学中Re和Ma的矛盾以及水动力学中Re数和Fr数的矛盾。在另外一些情况中，实验很难或者是儿乎不可能进行。比如实验测彊设备给流场带来的扰动， 或者是流动是不町解触的(液态硅的结晶过程丿；何些物理届用H前的设备难以测彊或者是测吊是难以 得到很好的精度等。实验在测杲一些总体参数，如阻力，升力，压差等是非常有效的，但是在很多情况卜，流场的 细节也非常匝耍，比如在技术改进和优化设计的时候。这时实验手段就显得非常的昂贵而H费时费 力了，因此，寻找另一种替代方案是非常必要的。随着电子计算机的出现，数值方法变成了-种现实的替代方案。尽管这些偏微分方程数值解法

20、的某木思想在一个世纪以前就被人发现了，但是由其庞人的计算杲，在计算机问世以前，很少有 人问津。1950s.计算机出现，Z后，计算机的性价比以级数的方式増长。计算速度和存储能力得到 了飞速的提高。但人们认识到了计算机的强人功能后，对数值求解流体力学问题产生了浓厚的兴趣, 数值方法在丁程实践中的地位变得越來越重耍了，逐渐形成了流体力学的第三个分支一计算流体 动力学(Computational Fluid Dynamics)简称 CFD。1.2 CFD的定义很多流体力学问题都可以用偏微分方程或(积分一微分方程)來描述，但是只冇在非常特殊的 情况卜才能得到解析解。为了得到流体力学微分方程的数值解，我们

21、用离散的方法，把原来的微分 方程近似成一个代数方程组，使其能在计算机上进行求解。近似公式应用在空间和时间的小域上， 从而数值解在离散的空间上给出数值结果。这门学科称为计算流体力学。实於流体力学的轿度取决 于测駅手段，而CFD的精度则取决于离散的质駅。1-3 CFD的优.缺点在实验流体力学中许多困难的事情，在CFD 可以得到岡满的解决。比如边界条件的模拟和尺 度效应:，在风洞中模拟汽车行驶过程,当空气运动时，要求地板也应以空气的速度移动,这在实验 中是难以做到的，至于复朵边界的模拟，在实於流体力学屮就更加闲难了，但是在CFD中可以得到 很好的解决；尺度效应也是如此。此外，CFD还可以得到流场的详

22、细倍息，便J：人们对流场仃更加 清晰的认识。CFD的缺点：由J：绝人多数流动都难以得到精确解，而数值结果都是近似解，因此人们很难对 数值结果进行评价，在数值结果和真实流动Z间是有差距的。数值计算的谋差在每一个环节都可能产生： 差分方程的近似和理想化 离散误差 方程求解过程中的迭代谋差 机器的舍入误差差分方程的近似和理想化：虽然从理论上说，只要徽分方程足够精确,理论上可以得到任意精度的数值解,但是实际上, 对许多流动很难得到粘确解，或者是数值解法不可行，这就必耍用些简化的模型。即使我们能 够获得微分方程的精确解，数值解也不是实际流动的克实反映。有的时候，模熨中还包含了试验数 据，这样，模型的真实

23、性就存在疑问。离散谋差可以通过构造高粘:度的插值和近似或者是采用更粘:细的网格來提岛粘度，但是这些都 需要更多的计算时间。在计算过程中，直接求解比迭代法精确，但是需要花费更多的资源，采用迭代方法则冋产生迭 代误差。计算结果的可视化是非常巫耍的，但是漂壳的图形可能并不能反映其实的流动，这-点对商 用软件更应引起重视。1-4 CFD的基本要素1.4.1数学模型数学模型是CFD的基础，通常需耍为特泄的工程问题建立专门的偏微分方程和边界条件。这 些特定的模型通常包含一些合理的简化和理想化。14.2离散方法在姓立数学模型Z后，需耍将微分方程或积分方程转化成仃限点上函数值的代数方程组后才能 在计算机上进行

24、数值求解，离散方法解决了如何将将微分方程或积分方程转化成代数方程组的问题。 常用的离散方法有有限差分法.有限体积法和冇限元法以及苴他的特殊方汰如边界尤法.有限解析 法，谱方法等。1.4.3坐标和向量系统流体力学的基本方程是坐标无关的，但是其貝体形式在不同的坐标系卜旳不同的表达形式，因 此在进行数值计算时必须选择一个合适的坐标系，此外，向吊在该坐标系卜的表达形式也必须事先 予以定义。曲线坐标系中矢量的两种表达形式:投影法则和平行四边形法则1.4.5数值的网格数值网格定义了所求物理量在空间的位置。数值网格是求解域的离散化表达，它将求解域划分 成若干个小的子域。数值网格的种类大致可分为： 结构化网格

25、(structured /regular grid)：结构化网格由多族网格线构成，同族的网格线互 不相交，并且和其他族网格线的任意一条有且只有一个交点。同族的网格线可用连续的整 数编号，求解域内的任意一个网格点都可以用一组(2维或3维指标唯一的标定.结构化网格在逻轲匕相当一个Cartesian网格，二维情况卜毎个网格节点令四个相邻节点, 每一个子域都是四边形。相邻巧点的连续件简化了编程并使得代数方程族11仃规则的給构 性，使得代数方程组的求解牢常的高效，苴缺点在丁只能用J几何形状简单的求解城.并 fl很难控制网格点在空间的分布，不仅浪费计算资源，由此产生的狭长单尤还会影响迭代 的收敛性。13

26、块结构化网格(Block- structured grid)：块结构化网格JUi 结构，在粗网格I .,网格(人的区域)町以是非规则的，而且允许相互重叠；在细网格上，每个区域采用结构化网 格。在区域的界面上需要进行特殊的处理。 非结构化网格(unstructured grid): nf以灵活的任意的求解域边界，尽管从原理匕来说，非 结构化网格町以用任何离散方法，但是采用有限体枳法和有限尤法故为恰当。1.4.6有限近似在选定数值网格以后，还必须确定数值离散过程中的近似方法。例如，何限差分法必须选择节点 处导数的近似公式，仃限体积法必须选择面积分和体积分的近似方法，仃限元法必须选择单尤函数 和权函

27、数。近似程度决定了数值求解的精度以及求解的难度和费用。通常高精度格式的方程中包含了更多网 格节点数，因此求解的匸作杲和难度也相应的増加。因此必须在两者Z间找到平衡点。1.4.7求解方法偏微分方程离散化后生成的非线性代数方程组的求解方法是多种多样的，不同的问题会有不同 的解法。由PCFD的代数方程组系数矩阵为人型稀疏矩阵，通常采用迭代方法求解。1.4.8收敛判据通常采用迭代法求解代数方程组，因此必须给出迭代终止条件，即收敛判据。1-5数值方法的基本特性15.1相容性当网格尺寸趋于o时,离散方法应该是梢确的。离散方程和精确解之间的误差称为截断误差。对于 相容的差分格式，当网格尺寸趋丁()时，截断误

28、差也趙丁0。1.5.2稳定性如果说在数值求解过程中出现的误差不被放大，则称该离散格式是稳定的.1.5.3收敛性如果当网格尺寸趋于o时，离散解和精确解的谋差也趋于o,则称该离散格式是收敛的。Lax定理：对J适定的线性偏微分方程的初值问题，与之相容的羞分格式的稳定性和收敛性是等价 L0的。 对j:非线性方程组，分析离散格式的收敛性非估困难，通常离散格式的收敛性需要用实验加以证明。1.5.4守恒性由于流体力学的基本方程是守恒方程，在数值求解过程中也应满足守恒关系。也就是说，在泄 常状态以及无外部源汇的情况f,守恒就流出和流入封闭体的吊保持守恒。对r强守恒型格式和仃 限体积法，每个控制体上都满足守恒律

29、，在求僻域内也自动满足，对其他离散方法，只要仔细选 择近似方法，也能满足守恒律。源汇的处理应满足一致性，即域内源汇的总杲应和守恒杲流出边界 的总量相等。由于数值计算过程中存在谋差,守恒性也是求解方法的匝要特性。如果守恒定律没有得到满足,则 会导致相用守恒駅在求解域内的不正确分布；年守恒也方程还存在人工源汇，由-方程本身是相容 和稳定的，在细网格上，解仍然是守怛的，但是在粗网格上会出现非守恒现彖。15.5 有界性数值解应在适当的界限之内，非负物理杲的计算结果不W出现负值，百分比应在0和100%之 间，对r拉普拉斯方程，最人和最小值应出现在边界上，诸如此类的特性应被数值近似所继承。1.5.6真实性

30、过复朵的模型难以被直接处理(如湍流，爆炸，影相流等)，所设计的模型和离散算法应能反 映物理的真实，不正确的模型可能导致物理上不存在的解或使计算发散。1.5.7精确性数值解永远是一种近似解，误差在整个数值计算过程屮都可能出现。误差人致町分为三类： 模型误差：实际流动和数学模型精确解之间的误差。离散误差：数学模型精确解和离散代数方程组精确解之间的误差.迭代误差：代数方程组迭代解和精确解之间的误差.通常高権度和求解效率是才盾的，册耍在两者Z间保持平衡。1.6离散方法1.6.1 有限差分法(Finite Difference Method, FD)点近似令限差分法是偏微分方程(parted diffe

31、irential equation,简称PDE)数值求解的瑕为占老的方法, 是Euler在18世纪提出的。对丁简单几何形体，这也是最简单的方法。冇限差分法的出发点是守恒型方程的微分形式，求解域用网格覆盖，在每一个网格节点上，通过将 偏导数近似为节点函数值的代数方程方式来近似原PDE。最终在毎一个节点匕都有一个代数方程，方程的未知数是中心节点以及相邻节点上的变最值。虽然FD理论上可用J：任意网格，但是采用结构化网格址为方便，在结构化网格屮，网格线可作为 局部坐标系的作坐标线。1.6.2有限体积元法(Finite Volume Method, FV)控制体内的平均近似FV方法的出发点是守恒型方程的

32、枳分形式，求解域被分成若干连续的控制体(CVs)。在毎一 个控制体上满足守恒方程。在每一个控制体的屮心作为讣算节点，计算该点上的物理量。控制体边 界上的函数值用巧点函数值的插值获得。体积分和而积分丿【适十的求枳公式近似。給果在每个控制 体上都有一个代数方程，未知数是中心节点以及相邻节点上变鼠的值。1.6.3 有限元法(Finite Eltnwnl Method)函数逼近仃限元法在很多地方和FV非常的类似，仃限尤法也将空间分成连续的控制体，但是通常采用 卄结构网格。FE法的最明显特点是方程在整个域内枳分Z前被乘匕了一个权两数，在最简单的FE 法中，每个控制体内的两数被假设成线性的，并IL保证解在

33、边界上的连续性。思考题：(1) 计算流体力学是如何定义的，计算流体力学的技术械础是什么？(2) 计算流体力学有那些基本耍素？(3丿为什么要对计算域划分数值网格，数值网格有什么作用？(4丿什么是数值方法的相容性？什么是数值方法的收敛性和稳定性，收敛性和穏定性Z间仃和关系？(6丿为什么说CFD结来一定是近似解？1113第2章有限差分法2.1引言所有的守恒方程都JUf+H似的结构.而丄都町以看作是输运方程的特殊形式。本章以一般输运方程 在Caneaan坐标系下的衷达式为例來讲述有限差分法。舟(“)+ot兀巧0)(2.1)#在上述方程中，除0外都认为是已知函数。基本概念#FD方法的第一步是离散求解域，

34、也就是定义数值网格。在FDM屮，网格是局部结构化的，每 个网格节点都可以行作是局部朋标系的原点，网格线则是局部坐标系的坐标线。同族的网格线两两 互不相交。每一个网格节点都可用一组指标唯一的标定.差分形式的标最守恒方6(2.1)是FD法的 原始方程。并被近似为以网格卩点上的守恒战为未知数的代数方程系统。代数方程组的解近似为原 微分方程的解。每一个带冇未知数的节点都必须冇一个代数方程，在节点以及相邻节点上的未知数之间建立联 系。这个代数方程用在接点处用何限差分近似代替偏导数的形式获得。対J* Dinchlet边界条件，边 界上不需耍代数方程，对其他边界条件，则必须将边界条件离散以得到所岛的代数方程

35、。有限差分的概念是从导数的定义中得到的：(2.2)闻=伽飒兀+山)一0(七) dx). aitoA.v几种常用的羞分格式：向前差分(forward difference)(23)切=0 - 0t向后差分(backward diffr ence)(2.4)中心差分(central difference)(2.5)#精确解2.3 一阶导数的近似（2.1）式中的对流项亚型需要对一阶导数进行离散dx2.3.1 Taylor级数展开法任意的连续丙数Kx）可在X】的领域内展开成Taylor级数 如咖+（-£）+呼隍+(x-xj 3!2I (x_xjn+ H(26)H为高阶项，用x宀代桥x这些点处

36、的两数值在X附近的展开式。利用这些展开式可得一阶和 高阶导数在£的近似表达式。例如(2.7)当网格间距较小时,高阶项很小，在各级数中截取第一项,得：(28)(29)152如 + 30 - 60t + 如6A.V+ 0(灯)(2.14)U.10)(2.12)上而三个公式依次FDS(forward difference scheme), BDS(backward difference scheme) CDS(central difFerence scheme)被截去的项称为截断谋差。FDS和BDS为一阶精度差分格式当网格均匀时，CDS貝有二阶精度。2.3.2多项式拟合法在构造差分格式时.

37、也可以用务项式曲线(或样条曲线)来拟介函数.然厉用拟介曲线的导数来近 似原旳数的导数。例如采用抛物线拟合Xbgz三点,可得：(2.13)缎 _ 0屮(£ )' + 0T(A兀£)' + 0 k尸采用不同的插值曲线可获得不同的差分格式,如(2.17)L17(2.17)L#(喲=20出3?； 60+ _ 0+ + 0(厶1丁)(2.16)对于发展方程:du du 八 + a=0 dldx< u(.v,0) = f(x)"(0,/) = 0(/)(0<x <l,a > 0),v>0/>0(2.15)(2.14r2.16

38、)依次为三阶精度的BDS. FDS和4阶精度的CDS2.3.3迎风格式(2.17)L#(2.17)L#方程的解为:(2.17)L#(2.27)L#f(x-at) x-at>0(2.17)"(心煦丄)必0Uj =你利用线性插值：“厂心+丸土(町-吩)XJ 一由此可得差分格式：(218)(2.19)(2.20)当W0时,同理可得:(2.21)(2.22)Aa(2.23)结合(2.22)和(223)可得:（2.24丿式屮，。相当J：流体力学动彊方程屮的速度m差分格式（2.24）反映了上游対卜游影响较人这一 事实，这种差分格式称为迎风格式（upwind scheme）2.4二阶导数的近

39、似方程（2.1）中的扩散项为二阶导数项，二阶导数可以用一阶导数的差分得到，也可以用Taylor级数法得到二阶导数的差分。采用Tavlor级数法可得二阶导数的差分(2.25)也可以用一阶偏导数的差分得到:(2.26)去丿j斗H 一兀*出一兀妇 - 00 - 0T为网格不等距时.后两种方法获得的差分格式是不相容的。0划=如（X,二心）+ 0T险 一 XJ- 0（X出 一 J ） 肝 J,（祐-V1X心-兀）（兀-V1）/ 22.5混合导数的近似混介偏导数的差分可以用一阶导数的差分得到dxdy dxdyJ(2.27)L19(234)2.6边界条件在网格的内点，用差分格式來离散PDE,为了得到唯一的解

40、，还岛耍在边界点上给出边界信息。 对J - Dmchlet边界条件，边界上的值直接给出，无需用新的方程。对J; Neunnnn边界条件和混含边 界条件，则需在边界上对边界条件进行离散。对邻近边界的内点，当采用高阶精度格式时，貝差分格式也应作特殊处理。- g + 6Q + 1803 + 100 - 33060A.V+ O仏卄)(2.28)-210 + 960+180-2400+1470 十。(心),)60Aa(229丿在边界上徑十皿=01办丿1X2 - X.(2.30)采用高阶差分格式仃一 0 (心 一 X J' + 0 (“ 一 X J 01K" X J (X? X(2.31

41、)#(234)#(234)如果网格是均匀的,(Z31)简化为:(232丿_ 03 + 4g 一 302A.V2.7时间项的处理由定解条件的不同，时间项和空间项的处理细节也仃差异。时间项的定解条件是初始条件， 这决定了时间项的处理必须采用时间推进解法，另外网格的划分是从空间上进行划分的，网格数员 和代数方程的数昴是对应的，为了保证方程的封闭熨，毎次求解时，未知数只能是一个时间层面上 的物理氐而其他时间层血上的物理磺必须是匕知值。在现实屮，现在时刻的世界是继承历史的， 可以对未來产生影响，但是不会依赖未來事物的发展趋孙。因此在仃限差分法中，仅仃最后个 时间层而上的物理起是未知的，苴余时间层面上的物

42、理员必须为已知最。2.7.1 初值问题的的基本解法考老初值问题：学 5,0(/) dl0(0) = g对（2.33）式积分可得:这里0叫0"分别表示如）和做f”），而不是0的幕。在具体数值计算中,方程（234丿的右端积分可用多种方式近似。例如宀 0W)/(235)(2.36)宀 0* 2,2)$n+-(237)在（235）中，当0"已知时,0弼可以直接求得,因此称为显格式（explicit method）,也称为forwardEuler method. （2.36）为Baclvard Euler method （2.37）为中值定理，都盅要通过求解方程的形式求 出0“T,因

43、此都是隐格式。2.7.2多层格式和预估校正格式为了提高计算精度，可以采用多层格式或预估校正格式，前者是采用多时间层来构造高阶近似 公式，后者则通过预测步提高积分表达式的计算粘度。预估一校止格式：0二=0" + /（/”，0"）/（预佔步，佔计的近似值，0二是0"的近似值）0网=0" +斗（几“，0"）+几"0二）心（校正步，用预测值來提高计算精度） 多层格式：宀 0" +寻（23他0）- 16他十严）+ 5仏2，严）采用了三个时间层的偏心差分来提高计算精度。212.7.3 Runger Kutta 法RungerKutta法

44、是一种高阶楷度的预佔一校止方法.但是它用了两个时间层的屮间节点的你 而不是多个时间层的值二阶 Runger Kutta0= 0"+分”,硏）卄一 n+2 2卩q 阶 Runger Kutta 法：心=0” + ¥他，0“) 叫20：孑0" +什几灵,喰) 0；+i = 0" + A/a 1，0"1)n+ rr+f(tn,l)+2f(tJ+2/(f“ + 兀卄心)n+ 刀+n-h- n+777y12.8有限差分法的应用实例osx do(238)例i一维扩散方程 如仇 0(x，O) = /(.V) 处 0,r)= A(t) 0(/,/)= B(r)

45、釆用均匀网格，将直线段均匀划分为N份，网格节点X. = /,心0,l,2,.,N + l,网格间 N距Ax = N,时间步长为"。空间扩散项采用中心差分格式:(23刃#(240)对J时间项则分别采用三种不同格式：A： FTCS格式(时间向前差分,空间屮心差分)A/(A.V)2上式中物理杲0的上标表示时间步。令O二丄丝(3严=昭+ (1 - 2咖+叭边界条件：(2.41)23(240)#(240)(242)山采用时间推进算迄 在计算n+1时间层时，11时间层已求解完毕，(2.34丿式可以根据上一时 间层的结果rt接获得半前时间层的物理鼠而无需求解代数方程组.这样的差分格式称为显格式。B

46、： CTCS格式(2.43)“姐-20” +始，0嘗_2犷+绘 + 2( 側(时整理后得：刖T + (1+ 層 _ 号娟=+ Q-60： + 号鑑 (2.44) 边界条件为(2.42)式。注意到(2.42丿式屮每个方程仃3个未知肥 0嘗、0"和0嘗.需耍求解联立代数方程组。当差分格式离散代数方程组中具有多个未知元素而襦要通过求解代数方程组的形式来获得数值解时, 这种差分格式成为隐格式.(2.42丿式写成矩阵形式为:(10000a1+(7O2000a21+ aG00a1+ aa000O1+ O< 000000、'血)、0耳久+(1_60； +专 0；2 20妒1 +(1

47、-0=00N-1y.v-2 +(1 - 6札- + 耳 0；ajtn-H(0叶丿2 2y.v-i + (1 - 60： + yv-H1丿<*(心)>#(240)(2.45) 系数矩阵为三对角矩阵，可通过追赶法求解。通过计算机编程的方式对上而的两种差分格式进行求解 取 / = !, /(x) = sin2 N = 10. A= B = 0. r = l ,这样 2=0.1,取不同的条件数#"器进行数值求解,貝体计算结果见程序“CFD271”的运行结来。（此时运“CFD271S通过实际算例演示两种总分格式数值解的星别）结论：当7 0.5时，FTCS格式计算发散，CTCS格式在

48、任意条件数时计算结果都收敛。（通过对实际计算发现问题，引出卜一节：差分格式的稳定性分析） 例2.对流扩散方程:(2.46)j = 0丄2,N + 1,网格间0 + "0 -血=0 0(兀0) = fW 0(OJ)= 4(/) 如)=B(f)采用均匀网格，将直线段均匀划分为N份，网格节点X,=/,N距Ax = N,时间步长为&。空间扩散项采用中心差分格式：U2A/(如 一 0t) - F(2.47)#"器进行数值求解,貝体计算结果见程序“CFD271”的运行结来。#"器进行数值求解,貝体计算结果见程序“CFD271”的运行结来。对丁时间项则分别采用三种不同格

49、式：A：时间向前差分.空间中心差分：即欧拉显格式（Expliat methods）=_(0：i - 姐)+ TA/妬-20： +姐（心尸(2.48)#"器进行数值求解,貝体计算结果见程序“CFD271”的运行结来。#"器进行数值求解,貝体计算结果见程序“CFD271”的运行结来。上式中物理最0的上标农示时间步。令7=丄耸，5 = ,则有（心）Ax(249)铲=叽+ （1 一 2咖+皿；-沖-C）边界条件：(2.50)由J采用时间推进算法.在计算n+】时间层时.n时间层己求解完毕，（2.34丿式可以根据上一时 间口的结果直接获得当前时间口的物理最血无需求解代数方程组.这样的差

50、分格式称为显格式。B: Crunk Nicolson Method(2.51)T（姐-20” +始侧_2妒+ 0嘗 十 2（呦3整理后得:#-(y + 才)0豐 +(1 + 60：" _(y- j)0常=$ + 才M爲 + 0- 60； + 号 _ 了)姐(2.52)2.9差分格式的稳定性分析2.9-1 Von Neumann 方法VonNeumann方法的莫木思想足将方程的误差看作是一系列正弦波的吾加，若迭代收敛则谋差不 应放人。这样，每个子波的振幅都应趋于6或至少保持不变。例；讨论差分格式(2.41)的稳定性；记n时间层的误差£：=供-$：，其中0为准确解，吋为近似解，

51、两者都满足(2.41)式。宀昭+ (1-2咖+咙两式相减后得:旷=昭 + (1 - 2crk： + g；：(253)将求解域向x轴止负方向作周期性延拓，把如£看作周期函数，展开成富氏级数:N/2宀°严A.V/2(2.54)代入(2.39)得,N/2NilN/2=a £<e+(l-2o-)工 °；/虬 + <7 工必”虹田RN/2A-V/24N/2(255)A/2 工甲严kN/2利用富氏级数的正交性可得:(2.56)(257)谋差放人肉子1- 4(7Sill<1(2.56)即 <7 <1/225(257)再讨论隐格式（2.37

52、）的稳定性，谋差方程为：- + （1+ 0n+1 - = y+ （1_ 6稈 + 珞 将（2.40州入（2.45）得:27(257)#(257)(258)N/2-工(ctcos/Ay - 1 - a)a"'He'kXiN/2N/2=工(<TCOsRAx + 1 - (7)(1：严W/2#(257)#(257)放人因子1- 2a sin <1(2.59)2(a- 2 kz1+ zcrsin 2羞分格N无条件稳定。同样，采JU Von Neumann方法町以讨论基分格式（2.47）和（2.52）的稳泄性。 （比较显格式和隐格式的有缺点）#(3.7)#第3章有限体积法3.1有限体积法基本原理上一章讲到的有限差分法将数值网格的节点上定义为计算节点，并在网格节点上对微分形式的 流体基本方程进行离散，用网格节点上的物理彊的代数方程作为原PDE的近似。在木章所要学习的冇限体积法则采用了不同的离散形式。首先，冇限体积法离散的是积分形式 的流体力学基本方程：£ p(jN - nds = i'll V 羽s + LgdG(3.1)计算域用数值网格划分成疗干小控制体&