離散方程的誤差與物理特性分析-熱流問題的數(shù)值計算-課件-03

1、熱流問題的數(shù)值計算熱流問題的數(shù)值計算Numerical Simulations of Thermal & Fluid Problems第3章 離散方程的誤差與物理特性分析 主講主講 李炎鋒李炎鋒2008年7月 北京nin, i22n, i)Sxxut()(L離散方程的截斷誤差TE,指其微分算子與相應(yīng)的差分算子之間的差,例如：對一維非穩(wěn)態(tài)對流擴散問題，其微分算子為:離散方程的精確解,指在離散方程的求解過程中不引入舍入誤差，記為:3.1.1 截斷誤差及相容性3.1 離散方程的相容性、收斂性及穩(wěn)定性ni2n1inin1in1in1ini1ninit,xSx2x2ut)(Ln , init,

2、x)(L)(LTE則其截斷誤差為:其差分算子為:)x(O) t(OSxxutx2x2ut2n , in , i22n , in , i2n1inin1in1in1ini1ni)x, t(O)(L)(LTE2n, init,x離散方程的截差可以通過對差分方程的精確解作Taylor 展開來導(dǎo)出，對一維模型方程的顯式格式,把 ， 在點(i, n)的Taylor展開式代入該離散方程并整理，得1nin1i3.1.1.3.離散方程的相容性（Consistency） )x,t(Onmt/x定義：當(dāng)時間、空間的網(wǎng)格步長趨近于零時，如果離散方程得截差趨近于零，則稱此離散方程與微分方程相容。即離散方程逼近微分方程

3、。當(dāng)截差表達式呈 得形式時（m, n均大于零），該方程具有相容性。當(dāng)截差中含 有項時，相容性在一定條件下滿足。3.1.2.1 數(shù)值解的離散誤差定義：定義：在網(wǎng)格節(jié)點(i, n)上離散方程的精確解 偏離該點上相應(yīng)的微分方程精確解(i, n)的值,用 表示。3.1.2 離散誤差與收斂性nininini)n, i (即即: : 3.1.2.2 影響離散誤差的因素:離散誤差與離散方程的截斷誤差有關(guān)離散誤差與離散方程的截斷誤差有關(guān)。在相同的網(wǎng)格步長下，截斷誤差階數(shù)越高，離散誤差隨之減小；對同一離散格式，網(wǎng)格加密，離散誤差也會減小.注意注意: :網(wǎng)格不能無限加密，得出網(wǎng)格獨立解即可。原因：原因：1. 經(jīng)濟

4、原因； 2. 整個數(shù)值解的精度取決于求解區(qū)域上各個節(jié)點離散方程的截差，而對于鄰接邊界的節(jié)點，往往難以得到高階截差的表達式；3. 對于對流換熱問題，除了上述數(shù)學(xué)上的一些誤差需考慮外，更要顧及離散格式在物理特性上的表現(xiàn)；4. 過分細密的網(wǎng)格要大大增加計算機的運算次數(shù)，舍入誤差因而增加。當(dāng)時間與空間步長均趨近于零時，如果節(jié)點上的離散誤差都趨近于零，則稱該離散方程是收斂的。3.1.2.3 離散方程的收斂性:對于線性初值問題，離散格式的收斂性可由其穩(wěn)定性保證.3.1.3 舍入誤差與初值問題的穩(wěn)定性nini對給定的物理問題，其數(shù)值解的舍入誤差的大小取決于所采用的計算方法及所用的計算機的字長.若由計算機求得

5、的解為 ,則節(jié)點(i, n)上的舍入誤差 為：3.1.3.1數(shù)值解的舍入誤差nnniiinininininini) n , i () n , i (離散方程的數(shù)值解偏離相應(yīng)精確解的總誤差由離散誤差離散誤差與舍入誤差舍入誤差兩部分組成。設(shè)所得到的數(shù)值解 不存在迭代不完全誤差，則:ni3.1.3.2 數(shù)值解誤差的組成誤差的主要來源為離散誤差誤差的主要來源為離散誤差. .3.1.3.3 初值問題離散格式的穩(wěn)定性一個初值問題的離散格式，如果可以確保在任一時層計算中所引入的誤差都不會在以后各時層的計算中被不斷放大，以致變得無界，則稱此離散格式是穩(wěn)定的。問題的提出：對非穩(wěn)態(tài)問題，在初始條件給出過程中引入的

6、誤差或在某一時層計算中引入的誤差會不會在以后的各時層的計算中被逐漸放大？對于適當(dāng)?shù)木€性初值問題所建立起來的相容格式，穩(wěn)定性是收斂性的充分必要條件；對于非線性問題，離散方程的相容性與穩(wěn)定性僅是獲得收斂解的必要條件而非充分條件。 3.1.3.4 聯(lián)系收斂性與穩(wěn)定性的Lax原理 凡是穩(wěn)定的格式，任何一個信息或擾動在計算過程中被放大的程度總是有限的；凡是不穩(wěn)定的格式，無論什么誤差都會在計算過程中被不斷放大，以致當(dāng)計算的時間層足夠多時，所得得解變得毫無意義。02dxddxd2212345x例題例題1：設(shè)有 ，（0）= 0，（5）= 1。試用二階截差與四階截差的格式將該式離散，采用均勻布置的4個節(jié)點進行求

7、解并比較其結(jié)果.0) h2 ()h1 ( 4) h2 (1ii21i得:, i=2、3、4解：解：采用區(qū)域離散法A，x = h, 代入二階截差的表達式:21ii1i22x2xx2x1i1ix2x1i1i0)h121h121()h128h1216(h1230)h128h1216()h121h121(2i21i2i21i22i284x2xeeee微分方程的精確解為：對于節(jié)點對于節(jié)點2 2和和4 4，四階格式不適用，只能采用二階截差格式，四階格式不適用，只能采用二階截差格式. .則對節(jié)點3有:22i1ii1i2i22x12163016x如果采用四階得差分表達式不同格式計算結(jié)果的比較格式格式 1 23

8、 45 精確解精確解 00.0473 0.1350 0.3679 1二階格式二階格式00.0582 0.1552 0.3944 1四階格式四階格式00.0505 0.1348 0.3918 1區(qū)間數(shù) 48163264精確解 1 0.0582 0.0502 0.0480 0.0475 0.0473 0.0473 2 0.1552 0.1404 0.1364 0.1353 0.1350 0.1350 30.3944 0.3752 0.3697 0.3683 0.3679 0.3679 網(wǎng)格疏密的影響（二階格式） 目前一般認為擴散項采用二階精度截差，對流項采用二階到三階的離散格式比較合適. 例題2：

9、設(shè)一維非穩(wěn)態(tài)導(dǎo)熱問題： ，0 x 2 t 0, T = 100oC, , T (0, t) = 400oC， T(2，t) = 400oC22xTatT1xta22x00t 試取 ，作數(shù)值計算并分析其結(jié)果。解：解：由于有對稱性取半個厚度作為計算區(qū)域，采用區(qū)域離散法A，等分4個子區(qū)域，得到5個節(jié)點，則顯式差分方程化為：按此式對開始的4個時層進行計算，計算結(jié)果如表中所示：nin1in1i1niTTTT22xTtT0t15 . 0 x5 . 00 x48. 022xtxta)1 ( 2xT例題例題3:設(shè)有以下初值問題，0 x1，T = 2x,及0.52下進行數(shù)值計算并與精確試在解相比較。計算結(jié)果計算

10、結(jié)果解：5.02xta初值不穩(wěn)定性初值不穩(wěn)定性：初值問題中由于時間步長取值不當(dāng)而引起的不穩(wěn)定性。對于所研究的一維非穩(wěn)態(tài)導(dǎo)熱的顯式格式，是一個臨界值，稍一超過就會引起數(shù)值解的震蕩。在一定的空間步長下，時間步長超過一定的值就會引起數(shù)值解的振蕩。3.2分析初值問題穩(wěn)定性的Von Neumann方法 以第一類邊界條件的一維非穩(wěn)態(tài)導(dǎo)熱問題為例： 22TTatxLx 0)()0 ,(xFxT)(), 0(1tftT)(),(2tftLT3.2.1誤差矢量隨時間傳遞的規(guī)律21112xTTTatTTninininini)(0iixFT)(10tnfTn)(2tnfTnI，n = 1, 2, 3 ，i = 1,

11、 2, 3 I， i = 1, 2, 3 I-1采用空間中心差分的顯式格式，其離散方程為:)()21 (111ninininiTTrrTT令 ，則差分方程可化為：2xtar固定n，將此式對i = 1, 2, 3 I-1寫出，得:ninnInInnnInInnrTrTTTTTrrrrrrrrrrTTTT00)21 (0)21 ()21 (0)21 (0122112121211nnngTAT1FT0初始條件為:A為系數(shù)矩陣。(2b)(2a)3.2.1誤差矢量隨時間傳遞的規(guī)律上式可改寫成以下簡潔的形式:nnngTAT100 FTnnnnnTTATT)(11000TTnnA1由此，得:式（3）與（2）

12、相減，得:(3b)(3a)假設(shè)在給出初值時引入了誤差矢量，把與這一含誤差的初值相對應(yīng)的解記為 ，則T 對于線性初值問題，如果在邊值的計算中不引入誤差，則誤差矢量的傳遞規(guī)律與原差分方程完全一樣，可應(yīng)用原差分方程來分析誤差隨時間的推移而傳遞得情形。結(jié)論：3.2.3 Von neumann分析的基本思想 假定所計算得初值問題的邊界值是準(zhǔn)確無誤的，而在某時層的計算中引入了一個誤差矢量（即小擾動).如果這一擾動的強度隨時間的推移不斷增大，則這一格式是不穩(wěn)定的；反之，是穩(wěn)定的?；舅枷牖舅枷?)()(ttt為進行某一離散格式穩(wěn)定性的分析，可以把誤差矢量的一個諧波分量表達式代入離散方程，以得出相鄰兩個時層

13、間該諧波分量振幅之比。格式穩(wěn)定的條件要求：放大因子放大因子Iiett)()(21112xTTTatTTninininini2)()()() 1() 1(2iIIiiIIieeextaetttt得代入離散方程 解：將試分析一維非穩(wěn)態(tài)導(dǎo)熱顯式格式的穩(wěn)定性。例題例題1 112sin)(41122xta5 .02xta注意：注意：Von Neumann分析法是針對第一類邊界條分析法是針對第一類邊界條件問題而發(fā)展的，上述穩(wěn)定條件僅適合于內(nèi)部節(jié)點。件問題而發(fā)展的，上述穩(wěn)定條件僅適合于內(nèi)部節(jié)點。右端自動成立，左端成立需滿足條件穩(wěn)定得條件為：2sin)(4122xta)cos1)(21)()(2xtattt經(jīng)

14、整理，得例題例題2:試分析一維模型方程格式的穩(wěn)定性條件。假設(shè)一維模型方程中的源項為常數(shù)，在分析格式穩(wěn)定時不予考慮?？臻g導(dǎo)數(shù)采用中心差分。ninininininininiSxxut21111122解：一維模型方程)2()(2(112111nininininininixtaxtu/a Iiett)()( )2()(2(1)()(2IIIIeextaeextuttt代入，經(jīng)整理得把其中方程可寫為：令xtuc稱為Courant數(shù)，得sincos221)()(Icrrttt穩(wěn)定性要求：1sincos221Icrrsincos221Icrr112 r1c 代表了在復(fù)平面上的一個橢圓，其中心在實軸上的(1-

15、2r)處，長軸為2r，短軸為 c。為使 ，橢圓必須位于以原點為中心，半徑為1的圓內(nèi)，因此，必要條件為：120 rrc2212/2rc1)(2021rr12/2/0222121rcrc21xtar22ytarxtuc1ytuc2Von Neumann分析法的不足是不能考慮邊界條件的影響。分析法的不足是不能考慮邊界條件的影響。其中對于二維對流-擴散方程的顯式格式，其穩(wěn)定性的條件為：或或一維模型方程的顯式格式穩(wěn)定的充分必要條件為：3.2.4 關(guān)于顯式格式的穩(wěn)定性分析的進一步討論 DCBAnininini111的離散方程中,穩(wěn)定的條件為系數(shù)A、B、C均大于零。以一維問題為例,在形如：對于非穩(wěn)態(tài)擴散方程

16、的顯式格式，可采用正系數(shù)法來獲得穩(wěn)定性條件。在數(shù)學(xué)上滿足穩(wěn)定性條件的解，不能保證一定得出在數(shù)學(xué)上滿足穩(wěn)定性條件的解，不能保證一定得出具有物理意義得結(jié)果具有物理意義得結(jié)果.以上分析初值問題顯式格式穩(wěn)定性的方法只適用于以上分析初值問題顯式格式穩(wěn)定性的方法只適用于線性問題，不能直接用于非線性的初值問題差分格線性問題，不能直接用于非線性的初值問題差分格式的分析。式的分析。對此可用局部線性化的近似的處理方法對此可用局部線性化的近似的處理方法。對每一時層的計算，速度都取上一時層的值，然后對每一時層的計算，速度都取上一時層的值，然后按線性問題的分析法，找出該時層各節(jié)點上穩(wěn)定性按線性問題的分析法，找出該時層各

17、節(jié)點上穩(wěn)定性條件所允許的最大時間步長，并以其中的最小值作條件所允許的最大時間步長，并以其中的最小值作為向下一時層推進的時間間隔。為向下一時層推進的時間間隔。注意：注意：3.3離散方程的守恒性(Conservativeness) 相容性，相容性，收斂性，收斂性，穩(wěn)定性。穩(wěn)定性。數(shù)學(xué)特性：守恒性，守恒性，遷移性，遷移性，人工粘性；人工粘性；q三種主要的物理性能：在對流-擴散方程的離散過程中，擴散項二階導(dǎo)數(shù)的中心差分離散格式具有優(yōu)良的物理特性和計算精度，物理特性的不良表現(xiàn)都是由于對流項的離散方式不完善所致。定義：定義：如果對一個離散方程在定義域內(nèi)的任一有限空間內(nèi)作求和的運算（相當(dāng)于對微分方程做積分）

18、，所得的表達式滿足該區(qū)域上物理量守恒的關(guān)系時，稱該離散格式具有守恒性.3.3.1離散方程具有守恒性的定義3.3.2 用直接求和法分析對流項中心差分的守恒特性 )2(11111xuutninininininiI2xinoutl1l2 xI1（a）21212)()(111IIiiiIIininixuut2121)(lllldxxudxx任取一段有限區(qū)間l1, l2來分析，積分:有)()()()()()(1111221211IIIiIIIiiuuuuuuu為截面的流量（b）tuuxIIiiiIIinini21212/ )()()(111或該式表明在t時間間隔內(nèi)流入與流出某一區(qū)域中的通量之差等于該時間

19、間隔中該區(qū)域內(nèi)的增量。因而格式（因而格式（a a）具有守恒性具有守恒性。tuuxoutinIIinini)()()(211則式（b）可寫作:outinIIIiIIIiiuuuuuuuu)()()()()()(21)()(211111221211因為：3.3.3 保證離散方程具有守恒性的條件 1. 導(dǎo)出離散方程的控制方程是守恒型的。2. 在同一界面上各物理量(及有關(guān)物性)及的一階導(dǎo)數(shù)是連續(xù)的（即從界面兩側(cè)的兩個控制容積來寫出該界面上的值是相等的）。()()ePeE() () ePeExxx如果用F表示界面上的流量,用 表示一階導(dǎo)數(shù)的離散形式，則連續(xù)性條件可表示為：xutxutniniinini2

20、111其中心差分格式為：一維純對流問題，非守恒型格式為：舉例：舉例：對右端作展開，求和的結(jié)果并不能寫成求和區(qū)域的進出邊界上通量的差，因而離散格式不具有守恒性。21212/)(111IIiiiiIIininitux對其作相應(yīng)的求和計算，得:3.3.4 對離散方程守恒性的評述 1. 采用守恒性的離散方程計算的結(jié)果能與原物理問題在守恒特性上保持一致；2. 可以使對任意大小的體積的計算結(jié)果具有對原離散格式所估計的誤差。當(dāng)對某個離散格式作截斷誤差分析時，是對任意取出的一個控制容積來進行的。當(dāng)把這一離散格式應(yīng)用于由若干個相互鄰接的控制容積所組成的有限容積時，如果界面上的通量可以相互抵消，則對該體積進行該物

21、理通量的數(shù)值計算的總體誤差只有在邊界上存在。 3. 一般地說，具有守恒特性的離散方程能給出比 較準(zhǔn)確的計算結(jié)果。 3.4離散方程的遷移特性 3.4.1 對流與擴散現(xiàn)象在物理本質(zhì)上的區(qū)別擴散擴散是由于分子不規(guī)則熱運動所致，分子不規(guī)則熱運動對空間不同方向的幾率都是一樣的，因此擴散過程可以把發(fā)生在某一地點上的擾動的影響向各個方向傳遞；對流對流是流體微團宏觀的定向運動，帶有強烈的方向性。在對流的作用下，發(fā)生在某一地點的擾動只能向其下游方向傳遞而不會逆向傳播。3.4.2擴散項中心差分可以將擾動均勻地向四周傳遞 22xt21112xtninininini的顯式格式:一維非穩(wěn)態(tài)擴散方程:證明：假定初始物理量

22、場均勻，且值為零。從某一時刻開始（如第n時層）,在某一節(jié)點i上突然有一個擾動，而其余各點上的擾動為零，隨著時間的推移，擾動傳遞的情形可由差分方程來確定。21112xtninininini對節(jié)點i：011nini)21 ()21 (221xtxtnini2/12xt所以有根據(jù)穩(wěn)定性要求:其中其中對節(jié)點i+1：021nini2121112xtninininini因此12102xt211xtni211xtni同理，對節(jié)點i1:因此，顯然n時刻發(fā)生在節(jié)點的擾動已經(jīng)均勻地向兩側(cè)傳遞。3.4.3對流項離散格式的遷移特性如果對流項的某種離散格式僅能使擾動沿著流動方向傳遞，則稱此格式具有遷移特性。3.4.3.1 定義0 xutxutnininini2111有:將中心差分用于一維非穩(wěn)態(tài)純對流方程的非守恒形式：證明：3.4.3.2 對流項的中心差分不具有遷移特性xutnininini22111021nini11()2niu tx因此其中對于節(jié)點i+1在n+1時層有:021ninixtuni2