第三章有限差分法

1第3章有限差分法1.1波动方程式的差分法（线性双曲线方程）即前进波波动方程（又称为移动方程或传递方程：convectionequation）00ccuuxt（3-1）)()0,()()0,(00xvxuxuxux（3-2）从此方程的差分求解方式分析常用的差分形式和稳定性条件。理论解：物理意义：波形保守不变，位置随时间以速度c前进。（前进波）f(x)cf(x-ct)ct2)(),(ctxftxu（3-3）ctxdvcctxutxu0001,（3-4）其中：)0,()()0,()(00xuxvxuxux（3-5）例：00000)(;1)(xxxuxxxu（3-6）即)()(00xxxv（3-7）其解为：)(),(0ctxutxu（3-8）3.1.1显式法对于发展型问题而言，当某一程度上的数值解已知，求解下一程度的数值解时，如未知的相只要一个，称为显式法（explicittimeintegrationmethod）。i.FTCS（ForwardinTimeandCentralDifferenceinSpace）方法tuuunjnjt1xuuunjnjx211（3-9）则能产生：xu1x030211xuuctuunjnjnjnj（3-10）变形后:njnjnjnjuuuu11121（3-11）这儿，为Courant数。xtc（3-12）Courant数表示物理的传播速度c和数值传播速度(x/t)的比值。该解的特性如图的三角形所示，nju的值由nju1和nju1所确定。当比值x/t保持一致时，不管x和t取多小，其影响的范围是一样的。当物理传播速度c比数值传播速度大的话，用此方法无法得到下一步长的物理特性。(如图绿线所示)，也即，当数值的影响领域无法包括物理特性领域，数值方法将不安定。1928年Courant、Friedrich、Lewy因此而提出了所谓Courant条件。FTCS法的解的发展xtTan-1(t/x)物理速度走的路程（Cu1）数值速度走的路程物理速度走的路程（Cu1）4即Courant条件为（CFL条件）1（3-13）但是波动方程不能由此方法判别的例子有：)(0),0(),(1),0(00xxxuxxxu（3-14）此问题有理论解，如图。例如，＝0.5时，时间步长为1/2x。解析解FTCS的解表1FTCS的解（＝0.5）xj-2xj-1xjxj+1xj+2t=011000t=115/41/400t=215/1623/169/1600t=353/6449/3275/6413/640其值是振荡不稳定的。随着时间的延续，振幅增加，甚至在正负值间振荡。此问题可从分析其差分方程的Fourior展开来分析。因此类展开可与频率振幅相角等相关。设FTCS格式的解的展开的某一分项为：ijjjegu（3-15）gj表示幅度，为与波数有关的相位角。定义一个时间步长前后的解的幅度比为振幅，用表示。则令它可用振幅和相位差来表示：ie（3-16）但偏离真实解时，振幅产生误差。偏离真实解时表示该波数代表的波不能以正确的速度移动。对于线性双曲波动问题，其理论解的振幅和某一波段上的相位差为1（3-17）且u可空间时间分离变量：即g为时间的变量。代入FTCS格式：xu1x0xu1x05sin11iggnn（3-18）其振幅为：)sin(tan),(sin1),(122（3-19）可见，无论是振幅还是相位都偏离真实解。尤其是振幅，它始终保持1，这说明随时间的发展，差分方程的解的幅度会无限制变大。这种不稳定称为VonNeumann不稳定。FTCS方法的不稳定也可从物理上理解。此差分的Tayor展开为...21)(612tuxucxuctuttxxx（3-20）代入波动方程...21)(6122tucxucxuctuxxxxxx（3-21）可见其截断误差的第一项为负的扩散值。负的扩散项从物理上是不稳定的。ii.Lax差分格式（Lax-Friedrich方法）FTCS式中的unj项用两边的平均值来代替:njnjnjnjnjuuuuu111112121（3-22）其：tantan,sincos,12222（3-23）在的一定范围内小于1。－11的范围内收敛。该方法使数值安定，但代价是解的扩散性增强。即耗散误差大。Lax格式Leap-Frog格式iii.Leap-Frog格式(蛙跳法)空间时间都用中心差分：022111xuuctuunjnjnjnj（3-24）xu1x06njnjnjnjuuuu1:11（3-25）当21时：2/1212sin1sintan,1,（3-26）因此，无论对什么信号，振幅都不变。因此无振幅误差。称为中立稳定。但当很小时，存在延迟相位误差（Laggingphaseerror）。即波的位置在真的波后发生。当21时：为纯虚数。2221sinsin（3-27）即当sin1/时不稳定。此方法分散误差大。iv.Lax-Wendroff格式原始Lax-Wendroff格式二阶Tayor展开，空间中心差分：为FTCS格式的修正。此方法由于有附加的扩散项，使得格式稳定。附加项包括空间和时间的2阶精度。在15.0的条件下，除了很大波数，主要含有延迟相位误差。Lax-Wendroff两步格式第一步（前1/2t）：第二步（后1/2t）：v.MacCormack的方法为航天领域应用最广的格式，为Lax-Wendroff两步格式的一种，但不计算j+1/2点的值。适用于非线性问题。对于线性问题。采用7预测修正因子法(predictor-corrector)。预测阶段：修正阶段结果同Lax-Wendroff两步格式。vi.1阶精度上风法时间向前，空间向后：Taylor展开右为截断误差。utt、uttt用uxx、uxxx表示：时间空间都为1阶精度。当＝0时，截断误差为零。0.51.0时：向前相位误差0.5时，延迟相位误差83.1.2显式法的小结1.粘性的附加i.FTCS格式njnjnjnjuuxtcuu1112(1)ii.Lax格式njnjnjnjnjnjnjuuuuuxtcuu111112212(2)iii.Lax-Wendroff格式njnjnjnjnjnjnjuuuxtcuuxtcuu112111222(3)iv.1次精度上风法njnjnjnjnjnjnjuuuxtcuuxtcuu11111222(4)上述所有的格式都可认为对FTCS的格式的修正，而且都是2阶微分方程utt的差分格式，为扩散项。其中Lax格式扩散最严重，Lax-Wedroff扩散最小，同2次精度中心差分接近。2.有限差分法的一般格式9线性波动微分格式可写成：0xtfu(5)xxcuf一般的单纯的Eular显示法为：njnjnjnjffxtuu111~~(6)f~称为数值流束(numericalflux)j-1/2j+1/2j+3/2j-1jj+1图x.通过界面的流束i.FTCS格式njnjnjfff12/121~(7)ii.Lax格式njnjnjfff1121~12/1(8)iii.Lax-Wendroff格式njnjnjfff111121~12/1(9)iv.1次精度上风法10011~2/1njfnjffnj(10)上式可表示为：jjjjnjuucfff112/121~(11)10此式在后面的高精度上风法的讨论中非常有用。3.1.3隐式法（implicittimeintegrationmethod）i.Grank-Nicolson格式022211111111xuuxuuctuunjnjnjnjnjnj(12)它可化成：njnjnjnjnjnjuuuuuu111111141414141(13)它符合对三角型法则（trapezoidalrule）11它的解同Leap-Frog格式一样振荡。隐式法的特点是，n时刻的任意坐标上的解的变化会影响n+1点上所有的点的解。情报无限制的快速传递。这是由于方法稳定，无CFL条件限制的缘故。事实上，隐式法也有不稳定的时候。ii.Beam-Warming格式：对上方法更一般化：02)1(2211111111xuuxuuctuunjnjnjnjnjnj(14)3.2扩散方程的差分法（抛物型方程）3.2.1一般形式1维波动方程？？（抛物性方程）的一般形式：22xuatu(15)流动和传热中常用形式稳态形式：1维稳定对流/扩散问题xxxu(16)Dirichlet边界条件：LLxxatxat00时，其解为xu1x01201/01LPeLxPeee(17Pe为Peclet数，定义为：uLPe(18由于此问题很简单，常被用于检验离散和求解方法。物理上，它代表了在流线方向对流与扩散间的平衡。事实上，很少有这种平衡起重要作用的流动。通常，对流与压力梯度或垂直流动方向的流动平衡。假定：u0andoL,u0或大值，Pe0,对流项可以忽略。解是线性的。Pe是大的，解在缓慢变化了一段后，在x=L附近很快变到L。此的突然变化往往成为对离散方法的考验。3.2.2差分方法i.Euler显式法显式差分格式：（空间中心差分）21112xatninininini(19)简化成：Pe0Pe0Pe=0LOLx013ninininini1112(20其中2xta(21)形式同波动方程一样，可利用VonNeumann稳定性条件，得到：21(22)ii.通用Crank-Nicolson法2112111111212xxatnininininininini(2314iii.边值问题3点计算分子方法，最后的代数方程形式：iiiWiiEiiPQAAA11(24)扩散项，采用中心差分〔CDS〕：11212121iiiiixxxxxx(25iiiiixxx1121；1121iiiiixxx(26dWdEdPiiiidWiiiidEAAAxxxxAxxxxA11111122(27对流项采用上风法（UDS）：;0if,;0if,1111uxxuuxxuxuiiiiiiii（28cWcEcPiicWiicEAAAxxuAxxuA