1、The Elements of Computational Fluid Dynamics第三章 发展型模型方程的有限差分 和有限体积方法3.1 一阶线性对流方程的差分格式3.2 抛物型模型方程对流扩散方程的 差分格式3.3 有限体积方法3.4 差分格式数值解的性质3.1 一阶线性对流方程的差分格式讨论双曲型模型方程:一阶线性对流方程0(3.1.1)uuatx线性对流方程的差分格式和流体力学中Euler方程的差分格式以及Navier-Stokes方程中对流项的差分格式有密切的关系,因此,掌握其差分格式的构造方法具有非常重要的意义。本节中,介绍的差分格式构造方法包括:(1)基于导数逼近(2)基于特
2、征理论(3)基于时间展开(4)基于算子分裂3.1.1 基于导数逼近的差分格式构造差分格式的最简单的方法。采用前差、后差和中心差等离散方法,直接近似微分方程中的导数项。1.Euler显式格式时间方向:前差。空间方向:中心差。1110(3.1.2)2nnnnjjjjuuuuatx11+1+12+1=021sin()1sin()1jjjjjikxnnjikxikxikxikxnnnnnnnnuA eAeA eA eA eatxa tAAik xxAa tGk xAx稳定性分析:取则格式是无条件不稳定的,没有应用价值。2.Euler隐式格式时间方向:后差。空间方向:中心差。1+1+1110(3.1.3
3、)2nnnnjjjjuuuuatx11+1+1+1+1+12=021+sin()111sin()jjjjjikxnnjikxikxikxikxnnnnnnnnuA eAeA eAeAeatxa tAik xAxAGAa tk xxt稳定性分析:取则格式是无条件稳定的。在数值计算中可以采用较大的,有助于提高收敛速度和计算效率。2.()()LT EOxt 截断误差:3.蛙跳(Leap-Frog)格式时间方向:中心差分。空间方向:中心差分。11110(3.1.4)22nnnnjjjjuuuuatx11+11+1+1+11+1122=0222 sin()02 sin()10sin()1sin()jjj
4、jjikxnnjikxikxikxikxnnnnnnnnnnnnuA eAeAeAeAeatxAAa tik xAAxAAGAAa tGik x Gxa ta tGik xk xxx 稳定性分析:取则可以认为,所以,22211sin()0sin()sin()111sin()0a tk xxGa ta ta tk xk xk xxxx 221CFLcourant.()()ata tccxxLT EOxt 稳定性的条件:。无量纲数,称为数或数。截断误差:在满足稳定性的条件时,放大因子等于1,格式具有零耗散,称为中性稳定的。4.一阶迎风(upwind)和顺风(downwind)格式时间方向:前差。空
5、间方向:前差或后差。110(3.1.5)nnnnjjjjuuuuatx110(3.1.6)nnnnjjjjuuuuatx.()()LT EOxt 截断误差:Courant Friedrichs Lewy101(1)0(3.1.5)0(3.1.5)10(,)()(3.1.5)(2)0(3.1.6)10(3.1.6)nnjjataacxau x tuxatuuatacax稳定性分析:时,无条件不稳定;时,的稳定条件是。时,解析解,波动从左向右传播,不受的影响,违背了波的传播特性,从而造成格式的不稳定。时,的稳定条件是;时,无条件不稳定。1111(upwind scheme)0(0)0(0)1nnn
6、njjjjnnnnjjjjuuuuaatxuuuuaatxatcx迎风格式:稳定的条件:1111(downwind scheme)0(0)0(0)nnnnjjjjnnnnjjjjuuuuaatxuuuuaatx顺风格式:无条件不稳定11100022nnnnnnjjjjjjaauuuuuuaaaatxx为了综合考虑和两种情况,迎风格式可以改写为:3.1.2 基于特征线理论的差分格式,CFL条件特征性质是双曲型方程的重要特点。在构造差分格式时,考虑微分方程的数学物理性质,有助于得到性态较好的差分格式。0(0)uuaatx图3.1:线性对流方程的特征性质。xatcuconst特征线:,沿特征线。1(
7、,(1)(,)njjOPjuu xntuuu xa t n t PPABCD由于 不在网格点上,点的值必须通过、等点插值获得。111(1)()()nnjjnnCBPBjjBCuuuuuuxa tuuxa txx 采用、两点线性插值 11110()(1)nnnnjjjjnnnjjjuuuuatxtcaxuc ucu或记,则 迎风格式.()()LT EOxt 截断误差:01(0)ca稳定条件:1111111111(2)12()022122nnnnnjjjjjDBPBnnnnnjjjjjBDuuuuuuuuuxa taxtxcuuuuu 采用、两点线性插值 或LaxLax-Friedrich)s(或
8、格式格式.()()LT EOxt 截断误差:1c 稳定条件:2012012012201222(3)A()()()()()()2()()22AABBCCCAACBPBBCu xaa xa xaaau xuu xuu xuaaaPu xaa xa xuuuuuuuxa txa txx 采用、三点二次插值设二次曲线,为待定系数首先,由,可以确定,。然后,把 点坐标代入,得 1212212Warming-Bea342()22mnnnnnjjjjjnnnjjjuuuuuatauuutxx即格式22.()()LT EOxt 截断误差:02(0)ca稳定条件:2012012012201222(4)()()(
9、)(),()()2()()22BBCCDDBDCDBPBBCDu xaa xa xaaau xuu xuu xuaaaPu xaa xa xuuuuuuua ta txx 采用、三点二次插值设二次曲线,为待定系数首先,由,可以确定,。然后,把 点坐标代入,得 1211112Lax-Wendroff2()22nnnnjjjjnnnjjjuuuuatauuutxx格式即22.()()LT EOxt 截断误差:1c 稳定条件:1.0warming-BeamLaxLax-Wendroff3.njnjauauaP小结:当时,用及其左侧的点进行插值,得到的是迎风格式。(格式可以看做二阶迎风格式),稳定条件
10、与 的正负有关。2.当用及其与之对称的点进行插值时,得到的是中心型格式(格式,格式),稳定性条件与 的正负无关。由插值点内插出 点的值,差分格式才是稳定的。Lax定理:双曲型方程的差分格式收敛(在等价性定理满足时亦即稳定)的必要条件是差分格式的依赖域包含微分方程的依赖域。An numerical method can be convergentonly if its numerical domain of dependence contains the true domain of dependence of the PDE,at least CFL(Courant-Friin the lim
11、it ase drichs-Lew and g y)o tdxdt:条件o zero.3.1.3 基于时间展开的差分格式112321223Taylor1)()2!1()()()(3.1.13)2!nnjjnnnnjjtjttjtxttxtxxtxnnnnjjxjxxjuuuuututOtuauuaua ua uuua uta utOt 把在作展开,得TaLaylx-WeorndroffCauchy-Kowalewski把展开式中的时间导数项用空间导数代替。是推导差分格式的方常法或方法:用技巧。1211112(3.1.13)222Lax-Wendroffxxxnnnnjjjjnnnjjjuuuu
12、uuatauuutxx和都采用中心差分离散,则改写为格式:0(3.1.1)uuatx12122121212122Warmin(1)034222(2)034g-22B a2e mxxxnnnnnjjjjjnnnjjjnnnnnjjjjjnnnjjjuuauuuuuatauuutxxauuuuuatauuutxx用二阶迎风差分离散,也用迎风离散,则 格式:1123123213Taylor1)()2!1)()2!1Crank-Nicolso11)()22n!2nnjjnnnnjjtjttjnnnntjtjttjtttjnnnttjtjtjuuuuututOtuuututOtutututOt 隐式格式
13、把在作展开,得由得格式:1131311)()22)()2nnnnjjtjtjnnnjxjxjuuututOta tuuuOt 所以,11111111),)22nnxjxjnnnnjjjjnnxjxjuuuuuuuuxx对和用中心差分离散111111111111114444044nnnnnnjjjjjjnnjjnnnnjjjja ta ta ta tuuuuuuxxxxuuaauuuutxx则即Crank-Nicolson格式是无条件稳定的。3.1.4 基于算子分裂方法的格式211111212122212212211Lax-W112endr221112o11f2fnnnnnnnjjjjjjjnn
14、jxxxxjnxxxxjnxxxjnnjxxjnnjxjccuuuuuuuccuuccucccuuccuucu 格式改写成算子:形式令1.MacCormack+格式 两步格式(预测步:校正步)1111112211122nnnjjxxjnnjxjuuccuucu 则1111111112MacC12ormacknnnjjxjnnnnjjjxjnnnjjxjnnnnjjjxjuucuuuucuuucuuuucu 所以,预测步:格一种形式校正步:预测步:另一种形式校正步:式:MacCormackMacCormack对于线性对流方程,格式和Lax-Wendroff格式是等价的。但是,对于非线性问题,如流
15、体力学方程组,二者不等价,格式更为简单有效。Warming-Beam2.两步的格式12122121212122Warmin(1)034222(2)034g-22B a2e mxxxnnnnnjjjjjnnnjjjnnnnnjjjjjnnnjjjuuauuuuuatauuutxxauuuuuatauuutxx用二阶迎风差分离散,也用迎风离散,则 格式:212121222220342221221111(1)(1)2222111(1)(1)222nnnjjjnnnnnjjjjjnxxxjnxxxxjnxxxjauuucuucuuuccuccccucccu 考虑11112Warming1-Beam2n
16、nnjjxjnnnnnjjjxjxjuucuuuucucu 所以,预测步:一种形式校正步:两步的格式:212121222234222122111(1)22nnnjjjnnnnnjjjjjnxxxjnxxxxjuuucuucuuuccucccu 1121Warming-112Beam2nnnjjxjnnnnjjxxjxjuucuuucucu 两步的格所以,预测步:另一种形式校:正步:式3.Runge-Kutta方法,半离散格式112345Taylor111)()2!3!4!1111111234njnnnnnnjjtjttjtttjttttjnjuuuututututOtttttutttt 把作展
17、开,得(1)(2)(1)(3)(2)1(3)4321nnjjjnjjjnjjjnnjjjtuuuttuuuttuuuttuuut 写成分步形式:(1)(2)(1)(3)(2)1(3)4321nnjjjnjjjnjjjnnjjjtuuauxtuuauxtuuauxtuuaut()()0(1,2,3)kkjjuuatxuuaktx Runge-Kutta上式中,空间导数可以采用前面介绍的方法进行离散,基于上述过程的格式(四为)称级方法。Runge-Kutta对于线性偏微分方程,无论空间方向如何离散,四级的方法时间方向具有四阶精度。对于非线性偏微分方程,格式在时间方向一般只有注:二阶精度。()()(
18、1)(2)(1)(3)(2)1(3)Runge-Kutta(,)0(,)0()4)3Runge-Kutta)2)1xxkkjjnnjjjnjjjnjjjnnjjjuP u utuP u utuPttuuPtuuPtuuPtuuP 方法可以用来求解以下形式的方程:对空间导数进行离散后,记为根据利用半离散格式方法,得(3.1.28)()全离散格式1121120Jameson0(3.1.29)()(464)(3.1.30)2jjjjjjjjjjuuatxuPtctPuucuuuuu线性对流方程的半离散格式。提出一种半离散格:式:其中,例434434(3.1.29)()()0()uuuaxtxxuuu
19、axtxx 式相当于小正数 的中心差分离散,与是相容的。添加的目的在于保证格式的稳定性。(3.1.29)(3.1.30)(3.1.2Jameson8)把和代入,就得到了全离散格式,称为格式。Jameson:空间方向有二阶精度;对于线性问题,时间方向有格式四阶精度。2(1)(2)(1)(3)(2)(3.1.29),4(1 cos)sin,()111,432jikxjnnnnAuAetZAtZcik xAAZAAAZAAAZA 把代入,有其中,所:,稳析以定性分1(3)(2)(1)11111211112311111234nnnnnnnnnnnnAAZAAZAZAAZ AZAZAAZAZ AZAZA
20、123411111+12624nnAZZZZA 放大因子:11Jameson2 2nnAcA由得格式的稳定性条件:3.1.5 边界条件的数值处理0(0)uuaatx方程:axx给定物理边界条件12:(2)nnnbMMMxxuuu线性外推边界处格式的精度比内点低一阶时,不会影响数值计算的整体精度。,(0baabnnnnjMjMjjuuxLxxLMxxuatuujxu求解的初值问题。周期性边界求解域为无限大,初值是空间变量 的周期性函数。取计算域的长度等于周期的长度,即,被均匀划分为等份。在计算域的边界上,满足。在处理周期性边界条件时,可以把计算区域内的网格点向左、右两侧做周期性延拓,延拓出的网格
21、点称为。虚拟网格点上的函数值可以利用解的周期性确定虚拟网:条件格点0110)(3.1.32),(3.1.32)Mx xx在点求解差分方程,如果差分方程的模板点中包含虚拟网格点,则用周期性条件确定虚拟网格点上的因变量值。3.2 抛物型模型方程 对流扩散方程的差分格式2222220uuabxyuuatxuuxyuutx对流项的离散:借鉴线性对流方程中空间导数的离散方法;扩散项的离散:借鉴热传导方程中空间导数的离散方法。2222(3.2.1)uuuuuabtxyxyxtyt考虑二维的对流扩散方程在,平面上,方程是抛物型的,可以用时间推进的方法求解。Navier-Stokes方程的物理意义:包含了对流
22、和扩散的过程,可以看作方程的一个模型方程。3.2.1 求解域的离散和边界条件的处理0(,)Mxxy00(,)xyxy(,)MMxyxy(,)ijxy,abcdxybadcxyxxyyxMyMxxyyxyMM 求解域:方向等分为份,方向等分为份。,(,)(,)(,)ijnni jxyi x j ytn tu i x j y n tu0,1,0,0,1,(,)(,)(,(Dirichlet)(Neu)(,mann),)(,)anjjnannjjnnnnjjjju xy tf y tuf y tuxy tb y txuubuuxbx边界条件第一类边界条件条件第二类边界条:求解域的上下左右边界上都需要
23、边界条件。考虑左侧边界,:将边界条件给出的值赋给边界网格点上的数值解:件条件对边界条件的一阶近似:构造边界条件的:二阶近似:201211,22,0122,1,0,(,)(,),(,),(,),423ajnjnnnnjjnjjnjx xnnnjjjnju x y taa xa xu x y tbu x y tuu xy tuxa a auuxbu 在边界附近,令利用求出系数。从而得到边界条件的二阶近似为:3.2.2 差分格式2222(3.2.1)uuuuuabtxyxy1.FTCS格式1,1,1,1,11,1,1,1222222nnnnnnnnnnnni ji jijiji ji jiji ji
24、ji ji ji juuuuuuuuuuuuabtxyxy 122,22212222(2).(,)nnnnnni ji jxi jyi jxi jyi juuuuuuabtxyxyEELT EOxyt用差分算子表示:截断误差:2.一阶迎风格式2222122,22,2222.(,)nnnnnnnni ji jxi jxi jyi jyi jxi jyi juuuuuabtxyxyuuuuuuuuaabbtxxyyxyaaaabbbbaabbLT EOxyt二维对流扩散方程:其中,截断误差:11110022022nnnnnnjjjjjjnnnnjjxjxjuuatxuuuuuuaaaatxxuuu
25、uaaaatxx一阶线性对流方程:2323123123xxxxxxxxDxDx 利用微分算子与差分算子之间的关系:3.二阶迎风格式12222,22,2222+2222.(),(),nni ji jyynnnnxxxi jxi jyi jyi jnnxi jyi juuaabbuuuutxxyyuuxyLT EOxyt 截断误差:4.Crank-Nicolson格式:隐式格式11111111122221,1,00440444nnjjnnnnjjjjnnjjnnxjxjnni ji jnnnxi jxi jyi jyuuatxuuaauuuutxxuuauutxuuuuuabtxyxyuuabuu
26、uutxy一阶线性对流方程:二维对流扩散方程:1,222121,2222(3.2.7)2ni jnnnnxi jyi jxi jyi juuuuxyxy时间方向具有二阶精度!3.2.3 近似因式分解方法1,11,222121,22222244Crank(3.2.7)2,-Nicols4on1nni ji jnnnnxi jxi jyi jyi jnnnnxi jyi jxi jyi jxyxyxxuuabuuuutxyuuuuxyxya tb tttccxxxyc二维对流扩散方程的令则格式:221,22,14221(3.2.8)4422yynxyxyi jyynxxxyxyi jnncuccu
27、AubA利用边界条件,可以写成矩阵形式:五对角矩阵(Approximate Factorization,AF)(3.2.8)近似因式分解方法:近似求解。221,221,22221,221,144221142421688411()4242yynxxxyxyi jyynxxxxyyi jxyxyyxxynxyyxxyi jyynxxxxyyi jccuccuc cccuccuOt 左边:222122,(3.2.8)11142424422(Crank-Nicolson)yyyynnxxxxxxyyi jxyxyi jccccuu在时间方向,具有一阶所以,简化为下列格式:可以分两步求解上式,每一步只需
28、求精度。格式具有二阶精度解一个三对角。方程组。222,21,11424422142yynxxxxxxi jxyxyi jyynyyi ji jcccuucuu第一步:第二步:22,22,2222,1442211424216884yynxxxyxyi jyynxxxxyyi jxyxyyxxynxxyyxxyi jccuccuc cccu 保持时间方向二阶精度的做法:右边221,22,144221(3.2.8)4422yynxxxyxyi jyynxxxyxyi jccuccu由22122,22122,211442244221111424242421688yyyynnxxxxxyxyi jxyx
29、yi jyyyynnxxxxxxyyi jxxyyi jxyxyyxxyccccuuccccuuc ccc 得2221,22122,34111142424242()xynnyxxyi ji jyyyynnxxxxxxyyi jxxyyi juuccccuuOt 22122,111142424242()yyyynnxxxxxxyyi jxxyyi jccccuu在时间方向,具有二阶所以,上式的求解,可以分精度。两步进行。222,21,111424242142yynxxxxxxi jxxyyi jyynyyi ji jcccuucuu第一步:第二步:,21,210,0,10,0,0,0,0,1,1
30、4201(3.2.13)42(3.2.13)i jyyni jyyi jyynjyyjnjjnnjjjnnnjjjucuuicuuuuubuxuuxbu 边界条件:中间变量在边界上的处理由考虑左边界第一类边界条件:已知,利用直接计算。第二类边界条件:已知,未知。首先,由计算0,+11+10,0,0,0,2njnnnnjjjjuuuu,然后,利用计算,210,0,0,1242yynnjyyjjcuuu所以,0,0,0,jjjuuu如果对精度的要求不高,对最简单的边界处理方法是3.2.4 多维问题差分格式的稳定性分析22221222121,2222(),1(1)(01)()()()()1 4(1)
31、sx iyjnnnnnni ji jxi jyi jxi jyi ji k xk ynni jxnnuuutxyuuuuuutxyxyuA eAGA考虑二维热综合格式:把代入上式得传程,导方222222insin2214sinsin221 4(1)sinsin2214yxyyxxyyxxykykxkykxkykxBBGB 令,则1111022214(1)142142114142142BBBBBGGBBBBB 当时,令,差分格式无条件稳定。221110022214(1)142142141421421401140111142(12)2(12)111()()2(12)xyBBBBBGGBBBBBBG
32、BGBtxy 当时,令,当时,格式是不稳定的。当时,格式是稳定的,所以,格式稳定的条件为即220FTCS111()()2txy当时,综合格式化为格式,格式稳定的条件为22()41()FTCS22xxytxt 如果取跟一维问题相比,则,多。(一维问题格式维问题的稳定性条的件稳定性条件:一般更0)为严格。3.3 有限体积方法(Finite Volume):对积分型的守恒方程进行离散,把积分型方程近似为代数方程进行求有限解体积方法的方法。3.3.1 积分型守恒方程对微分型的守恒方程在有限体积上积分,利用Gauss公式,可以得到积分型守恒方程。00()uuatxuffautx考虑一维线性对流方程守恒形
33、式2112211122,0,(),()xxx xudxffff u xff u xt在控制体上积分,可得积分型的守恒方程:3.3.2 空间控制体,baabxxxxMxM 求解域:等分为份,1/21/21/2,(0,1,)aMbjaxxxxxxj xjM 1/21/21/21/2,2jjjjjjxxxxuxjj:数值解存储在控制体的中心,第 个即第 个控制体网格点。3.3.3 有限体积方法的全离散形式1/21/21/2111/21/21/21/21/21/21/21/21/2111/21/20,(,),()0(3.3.4)(jjjnnjnnjjxjjjjjjxnnxttnnjjxttxnxnju
34、dxfffaufautt tuudxfdtfdtu dxuxj对沿时间方向在 之间积分,则定义数值解在单第元内的平个控制体均值代入111/21/21/211/21/23.3.4)0(3.3.6)1,nnnnnnjjjjtnnnjnnjjtuufftxft tffdtt一维问题有限体积方法的标准形式。有限体积方,得全离散有限体积方法的表达式为其中,在之间的平均值,称法的解是因变量在单元上的平均值。数值。通量为111/21/21/21/211/21/20(3.3.6)1,()()(1)(1,2nnnnnnjjjjtnnjjnnjtnnjjnjuufftxffdtt tutttuutuj在中若能得到
35、数值通量,便可计算出下一时刻因变量在单元上的平均值。为了计算,必须求出之间的。已知量只有 时刻单元平均值。为了求出,必须有限体积方法的标解决两个问准形式题:由数值通量111/21/21/21/2,),()(reconstruction)1(2)(),()()(evolution)1.nnnntnnnnnjjjjtMtxuxuxt tutffudtt求出 时刻数值解沿 方向的分布,这一过程称为重构。由求出之间的以及,这一过程称为演化过程。有限体积方法:在两个过程中引入各种近似,从而把积分型方程化为代数方程。重构零阶重构:在每个控制体内1/21/2,()jjnnjxxxjuxu物理量的分布为常数1
36、/21/21/21/21/21/23/2,1,()()()()()1()jjjjjjjnnjjjxxxxnxnjnnjjjnjjxxxxnnjnjjjxjjjnjuxuD xxu dxuxuDuuDuxuD xxjuD xxdxu线性重构:在每个控制体内物理量满足线性分布利用平均值的计算公式:得线性和未重构中,:将延拓到知的确第单元定方法,1/21/21/21,1()()1jjjnjjnnjjjnnjjjxxxnnjjjuDxxuuDxuxuD xxjuuDx将延拓到第单元,1/21/21/211/21/2,1/21/22.(),()0()()(1)()()=jjjnnnnjjnxxxnnjj
37、njutt tnuxautuxa ttutuxa ttu演化过程确定单元界面处的在之间的演化过程,从而计算数值通量。对于线性对流方程,在已知 时刻数值解的分布后,可以利用特征线理论计算数值通量。考虑的情况,根据特征关系采用零阶,重构时所以111/21/21/21/2111/21/2111()1()00()nnnntnnnjjjttnnnjjjtnnnnjjjjnnnnjjjjnnjjff ut dtautff ut dtautuufftxuuuuauutx,同理,时间方向一阶精度代入,得,一阶迎风,空间方向格式一阶精度。1/21/21,1/21/21/221/21/21/2(2)0,()()(
38、)()()()21()+221jjnnnnjjjxxxnjjnnjjjnjnjjntnnjjjjjtnjajuxuD xxutuxa ttuD xa ttxxuDa ttxaff ut dta uDDttf时,如果采用线性重,构则所以121/211111/21/221111()+220+222(3.3.15)nntnjjjjtnnnnjjjjnnnnjjjjjjjjxaf ut dta uDDttuufftxxxauDuDt DDuuatxx代入,得1121212211211112(3.3.15)34222Warming-Beam222Lax-Wendnnjjjnnnnnjjjjjnnnjjj
39、nnjjnnjjjnnnnjjjjnnnjjjuuDxuuuuuatauuutxxuuuuDxuuuuatauuutxx在中,当取时,线性重构中,上式为格式。当取时,上式为roff时间方向二阶精度,空间方向二格式。阶精度。3.3.4 有限体积方法的半离散形式1/21/21/21/20,()jjxjjxudxfffutja第对个控制体1/21/21/21/21/21/21/21/21/21/21/21/21010jjxnxnjjjjjjjjjjju dxuxufftxufftxffffff将数值解的单元平均值代半离散形式的有限体积格式半离散的有限体积格式的入上式,得:,是,的近似值,称作数值通量
40、。,的构造需要经过重构和演化过程,但是演化过程比全离标准形式散方法简单。1/21/21/21/21/20()()jjjnjnjuuatxffaunutux考虑线性对流方程,假定重构过程已经完成。在半离散格式中,时间导数并未离散,所以,只需计算某一时刻的考虑 时刻,则1/21/21/21/21/2,1/21/21/21/211/2()0()()()()()()2jjnjnjnjxxxnnjnjnjnnnnjjjjnjjjnuxautuxutuxa ttufaufauxutuDa tt由重构得到的物理量分布在控制体界面两侧一般存在间断,即是多值的。根据特征理论,当时,采用零阶重构则;采用线性重构1
41、/21/21/21111/2111/21()2221.50.5()()2njnjjnnnnjjjjjjnnjjnnnjjjjnnjjnnnjjjjxutuDxxfa uDfa uDuuDfauuxuuaDfuux则;当时,具有迎风特性当时,中心型格式实际应用中,对半离散方法的时间导数项进行离散后,可以得到用于数值计算的全离散格式。Runge-Kutta时间方向的离散:方法1/21/2(1)(2)(1)(3)(2)1(3)14321jjjnnjjjnjjjnjjjnnjjjPffxtuuPtuuPtuuPtuuP定义求解过程可以写为1.2.3.4.小结:对于一维问题,有限体积方法和有限差分方法是
42、等价的。对于多维问题,当控制体为矩形或长方体时,在某些意义上也是等价的。有限差分方法的相容性、收敛性和稳定性概念以及稳定性分析方法,适用于有限体积方法。有限体积方法对复杂的求解域形状有较好的适应性,更适于计算包含间断的流场。3.4 差分格式数值解的性质0,0 xtxt 差分格式是对微分方程的离散意义上的近似。如果差分格式是收敛的,则当时,差分格式的数值解趋近于微分方程的差精确解。在分格式的修实际计正方程算中,和都是有限的,这时,可以通过分析解的性质,得到差分格式数值解的性质。3.4.1 修正方程111(,)0(,)(,)(,)(,)(,)nnnxjjnnnjjxxxxxnnnjjuQ u uu
43、utQtuutQuQ u uR u utR u uuuQ u uR u utuutQ :设微分方程的差分格式为,如果和差分格式等价的微分方程可以写为,其中只包含 及其空间导数,则称为差分格式的修修正方程正方程。1122322323230(0)0Taylor+0(3.4.3)2!3!2!3!nnnnjjjjnjuuaatxuuuuatxuututuuxuxuatttxxx 例:考虑一维线性波动方程的迎风格式:首先,把格式在处做展开,得到与差分格式等价的微分方程:其次,消去上式中的二阶及二阶以上的时间导数项,化222(3.4.3)0(*)2266(*)(*)mmmmtxttxxtttxxxuuua
44、txxta xta xuauuuuu为下列形式,自循环消元法:把重新写成通过及其各阶导数的线性组合,消去中与时间相关的二阶及以上导数项。220(*)2266txttxxtttxxxta xta xuauuuuutuxuttutxuxxutttuttxutxxu(*)1a2t2a x26t00 xxxu26a x0(*)2tt 2t2a t024t04a t x 0(*)2a tx2a t22at024a t024at x 222(*)12tt212t212a t0022(*)3a tt x 23a t 223at02222(*)34ata t xx 2234ata t x 32234atat
45、x 按列求和1a00(1)2a xc00022(231)6a xcc.(1)()22xxa xT Ec uOx修正方程的作用:1.分析差分格式的截断误差和精度由修正方程,得一阶迎风格式的截断误差:,所以,一阶迎风格式具有一阶精度。.分析差分格式解的特性,差分格式的耗散和频散特性。3.4.2 差分格式的耗散和频散22(1)(231)26txxxxxxa tcxa xa xuauc uccu其中,所以,一阶迎风格式的修正方程为20(3.4.10)(3.4.10)lllluuatxuuuatxx考虑线性波动方程的某一差分格式,修正方程为研究差分格式解的性质,相当于研究方程解的性质。00/20/2 1
46、(,0)()(,0)()Fourier()(0)mMik xmmMu xuxLLMu xuxuxAe设初始条件是周期性函数,周期为。当长度为 的区间被划分为等份时,可以展开为离散的级数02Fourier()(0)(3.4.10)(,)()(3.4.13)(3.4.13)(3.4.10)()()()()()()()()()mmmmmmik xmmik xmik xik xik xlmmmlmmllmmlmlmkuxAeku x tA t edA tea ikA t eikA t edtdA ta ikikA t 利用级数的正交性,可以只考虑波数为的分量式解的波数为的分量为把带入,得2dt22()
47、()()()()(0)(,)(0)lmlmllmlmlma ikiktmma ikiktik xmA tCeCAu x tAee由初始条件,2221211221(1)(1)()(1)2211(,)(0)(0)(1)0llllmlml mllmlll mlikxaktktmikx atmmktllmlmlu x tAeeAg egeaakuuatx数值解的振幅随时间变化。即其中,在相同的初始条件下,的解为()(,)(0)mikx atmu x tAe解析解的振幅不随时间变化。221(1)111lll mlktmmmmgeggg:当时,数值解的振幅随时间衰减,格式有正的耗散,是稳定的。当时,数值解
48、的振幅不变,格式有零耗散,是中性稳定的。当时,数值解的振幅随时间增长,格式数值解与解析解振幅的比值:差分格式的耗散特性(造有负耗散,是不稳定的。利用差分格式的耗散特性,可成数值解与解析解之间幅值的差别)以通过修正方程进行稳定性分析。221(1)42(1,2,)4(1,2,)lll mlktmgennn n由,在修正方程中,只有偶数阶导数项影响差分格式的耗散特性,因此,偶数阶导数项被称为修正方程中的耗散项。一般地,如果修正方程中所有阶项的系数均为非负,所有阶项的系数均为非正,则数值解的振幅不会随时在实际应用中间增长,差分,通常只考虑格式主项是耗散的。的影响。2222222212(1)(231)2
49、6(1)(1)2(1)021Von Neumanntxxxxxxllxxlmla xa xuauc uccuua xc ukxa xcc对于一阶精度的差分格式,如迎风格式,修正方程为耗散项的主项,所以,的符号主要由 决定。如果差分格式是耗散的,必须有,即,与方法得到的稳定性条件是一致的。2322432242414324Lax-Wendroff(1)(1)68(1)(1)8(1)081Von Neumanntxxxxxxxxllxxxxlmla xa xuauc ucc uua xcc ukxa xccc 对于二阶精度的差分格式,如格式,修正方程为耗散项的主项,所以,的符号主要由 决定。如果差分
50、格式是耗散的,必须有,即,与方法的结论相同。()(dispersion)0(,)(0)()(,)(0)mmmmikx atmikx ammkkuuau x tAeatxau x tAg e相速度:波数为波动的传播速度。如果相速度与有关,则称差分格式有色散(或者频散)。对于微分方程,精确解为相速度为不同波数的波传播速度都是,波形在传播过程中保持不变数值解与解析解之间相位的差别,取决于差分格。对于修正方程,解析解式的频散(色散)特性。)t2211(1)lllmlaakaaaa相速度:如果,称差分格式有正的色散;如果,称差分格式有负的色散。因此,对于差分方程的数值解,由于色散作用,不同波数的波的传播