求解擴散方程的幾種差分格式的比較

祖麗胡瑪爾·卡迪爾，阿米娜·沙比爾

(喀什大學(xué)數(shù)學(xué)與統(tǒng)計學(xué)院，新疆 喀什 844000)

0 引言

擴散方程是一類基本的運動方程，是拋物型偏微分方程一個很重要的分支，在眾多領(lǐng)域都有著廣泛的應(yīng)用，可用于環(huán)境科學(xué)、能源開發(fā)、流體力學(xué)和電子科學(xué)等許多領(lǐng)域[1].因此，數(shù)值求解拋物型偏微分方程問題具有重要的理論意義和實用價值.求解擴散方程的數(shù)值方法主要是有限差分法、有限元法、有限體積法等多種方法.但是對于對流占優(yōu)問題，通常是用差分法或有限元法進(jìn)行求解.本文利用三種差分格式求解了一維擴散方程的數(shù)值解，并給出了穩(wěn)定性比較，最后通過數(shù)值算例進(jìn)一步說明了數(shù)值解法的有效性.

1 差分格式的建立

考慮以下常系數(shù)擴散方程的初邊值問題：

其中，a>0為正常數(shù)，f(x,t),φ(x),α(t),β(t)都是已知函數(shù)，L為空間長度，T為時間長度.

首先，對求解區(qū)域D={(x,t)|0 ≤x≤L,0 ≤t≤T}作矩形網(wǎng)格剖分，將空間[0,L]等分m份，節(jié)點為xi=0+ih（0 ≤i≤m，h=）.再對時間區(qū)間[0,T]等分n份，節(jié)點為tj=0+jτ（0 ≤j≤n，.這樣得到(m+1)(n+1)個網(wǎng)格節(jié)點(xi,tj).

1.1 向前歐拉格式

在網(wǎng)格節(jié)點建立節(jié)點離散方程，即

關(guān)于時間的一階導(dǎo)數(shù)的向前差商形式為

上式誤差為O(τ).

類似地，關(guān)于空間的二階偏導(dǎo)數(shù)的中心差商形式為

上式誤差為O(h2)

將式（3）和（4）帶入（2）式中的第一式，再用數(shù)值解近似代替精確解u(xi,tj)，并忽略高階小項，即可得到向前歐拉格式：

經(jīng)過整理可得向前歐拉格式：

由（6）式知，第j+1 個時間層上的數(shù)值解可由第j層上的已知信息表示出來.如此一層一層計算，可得到所有網(wǎng)格點的數(shù)值解，計算簡單、便捷[2].

1.2 向后歐拉格式

對時間的偏導(dǎo)數(shù)用向后差商，關(guān)于空間的二階偏導(dǎo)數(shù)用中心差商代替，再用數(shù)值解近似代替精確解u(xi,tj)，并忽略高階小項，即可得到向后歐拉格式：

實際計算時，在每個時間層都需要求解一個線性方程組，計算比較復(fù)雜[3].

1.3 Crank-Nicolson 差分格式

對時間的一階偏導(dǎo)數(shù)用中心差商代替，則有

上式誤差為O(τ2).

接下來處理虛擬節(jié)點處關(guān)于空間的二階偏導(dǎo)數(shù)

將上面兩式代入式（10）得

上式和（9）式一起代入（8）式，用數(shù)值解近似代替精確解u(xi,tj)，并忽略高階小項可得Crank-Nicolson差分格式：

易見，此格式的階段誤差為O(τ2+h2).

整理上述格式，可得

上式的矩陣形式可以利用追趕法來求解.

2 穩(wěn)定性比較

向前歐拉格式[5]是條件穩(wěn)定的，其穩(wěn)定性條件為r=，如果網(wǎng)格比不滿足穩(wěn)定性條件，導(dǎo)致最后結(jié)果失效.

向后歐拉格式是無條件穩(wěn)定的，網(wǎng)格比很大情況下數(shù)值解仍然收斂.

Crank-Nicolson 差分格式[4]是二階無條件穩(wěn)定的，其中一個至關(guān)重要的因素就是將原項方程弱化，使之在相鄰時間層網(wǎng)格節(jié)點的終點處成立，而不是在網(wǎng)格節(jié)點處成立.

3 數(shù)值算例

例：計算拋物型方程的初值問題

已知其精確解為u(x,t)=tsin(x).

由表1 和圖1 可見，r=0.5 時數(shù)值結(jié)果完全吻合；當(dāng)網(wǎng)格步長不滿足穩(wěn)定性條件時，誤差會以指數(shù)增長，結(jié)果不收斂，數(shù)值結(jié)果無效.這表明向前歐拉格式是條件穩(wěn)定的.由表2 和圖2 可見，無論r的取值如何，也就是無論時間、空間步長如何選取，向后歐拉格式恒穩(wěn)定，即該格式是無條件穩(wěn)定的.由表3 和圖3 可見，Crank-Nicolson 格式是無條件穩(wěn)定的，差分格式是二階精度的.

表1 當(dāng)r=0.5時，向前歐拉格式的數(shù)值解和精確解

表2 向后歐拉格式在不同網(wǎng)格步長下的數(shù)值解

表3 Crank-Nicolson差分格式在不同網(wǎng)格步長下的數(shù)值解

圖1 向前歐拉格式

圖2 向后前歐拉格式

圖3 Crank-Nicolson格式

4 結(jié)論

本文討論了求解拋物型方程初邊值問題的幾種常用方法，如向前歐拉格式方法、向后歐拉格式方法、Crank-Nicolson 格式方法，并了解各種方法的優(yōu)缺點.結(jié)果表明，向前歐拉格式方法計算簡單、直接，但穩(wěn)定性差，只有在時間、空間步長的合理選取之下才能保證數(shù)值解的收斂性；向后歐拉格式方法則彌補了向前歐拉格式方法的缺點，不需要考慮時間、空間步長的選取，但每個時間層上都要求解一個線性方程組，計算量會加大.這兩個格式方法的精度都是時間一階、空間二階；Crank-Nicolson 格式方法是無條件穩(wěn)定的，時間和空間都是二階精度的.