含兩個(gè)空間分?jǐn)?shù)階導(dǎo)數(shù)熱方程的非標(biāo)準(zhǔn)有限差分解法

劉子婷

(廣東工業(yè)大學(xué) 應(yīng)用數(shù)學(xué)學(xué)院， 廣東 廣州 510520)

分?jǐn)?shù)階導(dǎo)數(shù)由于具有非局部性，可以有效地描述具有記憶和遺傳特征的自然過程。近幾十年來，分?jǐn)?shù)階偏微分方程在控制工程、物理力學(xué)、生物學(xué)和計(jì)算機(jī)科學(xué)等學(xué)科中得到了廣泛的應(yīng)用[1-5]。尤其空間分?jǐn)?shù)階偏微分方程常被用來模擬超擴(kuò)散，是因?yàn)槠渲辛Ｗ佑鹆鞯臄U(kuò)散速度比經(jīng)典布朗運(yùn)動(dòng)模型預(yù)測(cè)的要快。Chen等[6]考慮了具有分?jǐn)?shù)階擴(kuò)散和整數(shù)階平流項(xiàng)的空間分?jǐn)?shù)偏微分方程，并且通過有效地結(jié)合Legendre spectral和Crank-Nicolson差分法得到了目標(biāo)方程完全離散的數(shù)值解。Zhou等[7]建立了Riesz空間分?jǐn)?shù)階偏微分方程的Crank-Nicolson有限差分格式，并且討論了該格式的穩(wěn)定性和收斂性，其中理論分析和數(shù)值實(shí)驗(yàn)的結(jié)果表明，該方法是有效可行的。Ali等[8]利用改進(jìn)的擴(kuò)展tanh方法，并使用數(shù)學(xué)軟件Mathematica對(duì)兩類時(shí)空分?jǐn)?shù)階偏微分方程進(jìn)行處理，給出了這兩類方程的新結(jié)構(gòu)。Mickens[9]提出了求解微分方程的數(shù)值方法：非標(biāo)準(zhǔn)有限差分法，在某種意義上，該方法克服了經(jīng)典有限差分法的一些不足。Agarwal等[10]結(jié)合非標(biāo)準(zhǔn)有限差分法和切比雪夫配置法，得到了一類空間分?jǐn)?shù)擴(kuò)散方程的非標(biāo)準(zhǔn)差分格式，通過數(shù)值模擬，發(fā)現(xiàn)該方法所得到的數(shù)值結(jié)果不僅優(yōu)于以前的方法，而且大大提高了收斂精度，進(jìn)而驗(yàn)證了該方法的有效性和準(zhǔn)確性。

在本文中，我們考慮以下空間分?jǐn)?shù)階熱方程：

(1)

初始條件為

u(x，0)=s(x)，0

(2)

邊界條件為

u(0，t)=u(L，t)=0，0

(3)

因?yàn)闊岱匠讨猩婕暗絻蓚€(gè)分?jǐn)?shù)階導(dǎo)數(shù)，這會(huì)使得數(shù)值計(jì)算有一定的困難，需借助Grünwald-Letnikov公式和帶位移的Grünwald-Letnikov公式來離散兩個(gè)空間分?jǐn)?shù)階導(dǎo)數(shù)，從而降低數(shù)值計(jì)算的難度。由Agarwal和El-Sayed[10]的研究受到啟發(fā)，本文將結(jié)合非標(biāo)準(zhǔn)有限差分法和Crank-Nicolson差分法得到含兩個(gè)空間分?jǐn)?shù)階導(dǎo)數(shù)熱方程的非標(biāo)準(zhǔn)有限差分格式。

1 非標(biāo)準(zhǔn)有限差分格式

f(xj,tn)，1≤j≤m-1，n≥1，

(4)

u(xj,0)=s(xj)，1≤j≤m-1,

(5)

u(0,tn)=u(m,tn)=0,n≥1，

(6)

將非標(biāo)準(zhǔn)有限差分法和Crank-Nicolson差分法結(jié)合起來，得到式(4)的非標(biāo)準(zhǔn)有限差分格式，如下所示：

(7)

(8)

將式(8)改寫成矩陣形式

(I+A)Un+1=(I-A)Un+φ1(τ)Fn+1/2，

其中

I是一個(gè)(m-1)×(m-1)的單位矩陣，A=(ajn)是一個(gè)(m-1)×(m-1)的系數(shù)矩陣，其中

2 穩(wěn)定性與收斂性分析

定理2.1 當(dāng)0<α<1，1<β<2時(shí)，非標(biāo)準(zhǔn)有限差分格式(4)是無條件von Neumann穩(wěn)定的。

(9)

定義

(10)

其中i是虛數(shù)單位，q是空間波數(shù)。將式(10)代入式(9)中得到

則

(11)

令

(12)

證明首先計(jì)算局部截?cái)嗾`差R的階數(shù)，R的表達(dá)式為

其中

(13)

(14)

u(xj-(k-1)h,tn)]+o(h)，

(15)

將式(13)—(15)代入R的表達(dá)式中得到

u(xj-(k-1)h,tn)]-f(xj,tn+1/2)+o(τ2+h),

所以局部截?cái)嗾`差R的階數(shù)為o(τ2+h)。

(16)

將式(16)寫成矩陣形式為

(I+A)En+1=(I-A)En+φ1(τ)R，

(17)

令K=(I+A)-1(I-A)，P=(I+A)-1，式(17)可改為

En+1=KEn+PY，

(18)

3 數(shù)值算例

考慮下述問題

f(xj,tn)，0

(19)

u(x，0)=s(x)，0

(20)

u(0，t)=u(1，t)=0，0

(21)

1)令φ1(τ)=τ，φ2(h)=h。圖1所示的是當(dāng)0

(a)解析解 (b)數(shù)值解圖1 問題(19)—(21)的解(0

τ=hα=0.2,β=1.2最大絕對(duì)誤差誤差率α=0.5,β=1.5最大絕對(duì)誤差誤差率α=0.8,β=1.8最大絕對(duì)誤差誤差率1/107.2899×10-3—5.6731×10-3—4.2503×10-3—1/204.0691×10-31.792.9616×10-31.922.2225×10-31.911/402.1501×10-31.891.5545×10-31.911.1449×10-31.941/601.4501×10-31.481.0534×10-31.487.6952×10-41.491/801.1066×10-31.317.9596×10-41.325.7923×10-41.331/1009.1476×10-41.216.3892×10-41.254.6480×10-41.25

2)改變時(shí)間分母函數(shù)，空間分母函數(shù)保持不變，即φ1(τ)=1-e-τ，φ2(h)=h。比較表1和表2中的數(shù)據(jù)，可以得出最大絕對(duì)誤差在減小，這表明改變時(shí)間分母函數(shù)可以提高精度。

表2 問題(19)—(21)的最大絕對(duì)誤差值以及誤差率(φ1(τ)=1-e-τ，φ2(h)=h，T=1)

3)在上一個(gè)實(shí)驗(yàn)的基礎(chǔ)上繼續(xù)改變空間分母函數(shù)，即φ1(τ)=1-e-τ，φ2(h)=eh-1。比較表2和表3的數(shù)據(jù)可以得出最大絕對(duì)誤差在繼續(xù)減小，這表明同時(shí)改變時(shí)間和空間分母函數(shù)會(huì)更好地提高精度。

表3 問題(19)—(21)的最大絕對(duì)誤差值以及誤差率(φ1(τ)=1-e-τ，φ2(h)=eh-1，T=1)

4 總 結(jié)

在本文中，我們使用了具有時(shí)間和空間步長的分母函數(shù)來構(gòu)造空間分?jǐn)?shù)階熱方程的非標(biāo)準(zhǔn)有限差分格式，并證明了該差分格式的穩(wěn)定性和收斂性。在數(shù)值算例中，通過構(gòu)造不同的分母函數(shù)來減小最大誤差，為理論證明提供了依據(jù)。結(jié)果表明，非標(biāo)準(zhǔn)有限差分法用于求解含兩個(gè)空間分?jǐn)?shù)階導(dǎo)數(shù)的熱方程是有效且可行的。