朱 琳

(寧夏大學(xué)數(shù)學(xué)統(tǒng)計學(xué)院，寧夏銀川750021)

分?jǐn)?shù)階偏微分方程可用于描述一些反常擴(kuò)散現(xiàn)象，在物理、化學(xué)、水文地理學(xué)和經(jīng)濟(jì)學(xué)等諸多領(lǐng)域中應(yīng)用廣泛，近幾十年來，其理論研究得到了蓬勃的發(fā)展．解分?jǐn)?shù)階偏微分方程的數(shù)值方法主要有有限元法［1－2］、譜方法［3－4］、DG 方法［5］和有限差分方法［6－12］等．

由于分?jǐn)?shù)階算子具有非局部記憶性，使得分?jǐn)?shù)階偏微分方程的求解比整數(shù)階微分方程的求解需要花費(fèi)更多的計算時間和更高的存儲要求．文獻(xiàn)［13－14］利用非對稱迭代技術(shù)提出了一種半隱格式，既保留了隱格式穩(wěn)定性好的優(yōu)點(diǎn)，又不用求解線性代數(shù)方程組．本文擬采用這種非對稱迭代技術(shù)，結(jié)合加權(quán)移位 Grünwald －Letnikov算子［15］和中心差分算子構(gòu)造一個二階半隱式有限差分格式，求解一維分?jǐn)?shù)階對流擴(kuò)散方程．此格式形式上是隱式的，但是通過對奇數(shù)和偶數(shù)2個時間層分別選取不同的節(jié)點(diǎn)模板，從而可以進(jìn)行顯式計算，節(jié)約了計算時間和存儲空間，這種思想在本文中將進(jìn)行詳細(xì)介紹．

一些學(xué)者采用非對稱迭代技術(shù)構(gòu)造了許多不同類型的并行算法，比如在文獻(xiàn)［16－17］中，應(yīng)用對稱迭代技術(shù)解橢圓方程，僅給出了穩(wěn)定性，且只提到截斷誤差是 O(Δth－1+Δt+h);在文獻(xiàn)［18］中，只是在實(shí)際的計算中應(yīng)用了非對稱迭代技術(shù)．文獻(xiàn)［19］得到了空間分?jǐn)?shù)階對流方程解析解和離散解在 l2范數(shù)意義下誤差估計式為O(Δt2h－2(α－1)+Δt+h)．本文在此基礎(chǔ)上給出了空間分?jǐn)?shù)階對流擴(kuò)散方程的二階半隱有限差分格式，用Fourier分析法證明了此格式的穩(wěn)定性，給出了截斷誤差及解析解與離散解在l2范數(shù)意義下的誤差估計O(Δt2h－2(α－1)+Δt+h2)，其中 Δt和 h 分別為時間和空間方向的步長．

1 非對稱迭代格式

研究下列分?jǐn)?shù)階對流擴(kuò)散方程:其中，α 是分?jǐn)?shù)階空間項導(dǎo)數(shù)，且 1＜α≤2，d(x，t)≥0，c(x，t)≥0．

做出以下記號:

其中，n=0，1，…，≤T/Δt，l=0，1，…，K．

Riemann－Liouville［20］分?jǐn)?shù)階導(dǎo)數(shù)定義為其中，n＞α＞n－1≥0(n為整數(shù))，Γ 是伽馬函數(shù)．

計算Riemann－Liouville導(dǎo)數(shù)的經(jīng)典 Grünwald公式和移位 Grnüwald 公式［9－10］分別為:

其中，K 是正整數(shù)，gk是 Grünwald 系數(shù)，且

引理 1．1［21］當(dāng) 1 ＜ α ≤ 2，(4)式定義的Grünwald系數(shù){gk}時具有如下性質(zhì):g0=1，g1=－α ＜ 0，1 ≥ g2≥ g3≥ … ≥ 0，=0，≤0(m≥1)．

將(3)式中經(jīng)典和移位Grünwald算子加權(quán)平均用來近似計算Riemann－Liouville分?jǐn)?shù)階導(dǎo)數(shù)得

將(5)式簡化為

其中

由引理 1．1，經(jīng)過簡單推導(dǎo)，可得引理 1．2［15］．

用中心差分算子離散對流項得到

結(jié)合(6)和(7)式，分別采用向后或向前Euler公式離散時間導(dǎo)數(shù)得到方程(1)的2種隱式有限差分格式為:

或

接下來，結(jié)合非對稱迭代技術(shù)對(8a)和(8b)式進(jìn)行改造，設(shè) 1≤l≤K－1，n≤N=T/(2Δt)得到如下半隱式有限差分格式:

其中(9)式滿足的初邊值條件為:

下面介紹顯式解方程組(11)的過程．

1)求解方程(11a)，即從時間層t2n計算到時間層t2n+1．從給定的邊界值v20n+1開始，依次從左到右可以計算出 l=1，2，…，K－1的值 v2ln+1，在每次應(yīng)用(11a)式時，v2l－n1+1已經(jīng)知道，從而隱式格式進(jìn)行了顯式化的計算．

2)求解方程(11b)，即從時間層t2n+1計算到時間層t2n+2．從已知的邊界v2Kn+2開始，依次從右向左可以計算出 l=K－1，K－2，…，1 的值 v2ln+2，計算也是顯式的．

交替使用(11a)和(11b)式來進(jìn)行計算，就實(shí)現(xiàn)了對隱格式(11)式的顯式化計算，不用求解一個線性代數(shù)方程組．

2 差分格式的穩(wěn)定性分析

用Fourier方法討論非對稱迭代格式(9)的穩(wěn)定性．

其中，A=(Aij)(K－1)×(K－1)，B=(Bij)(K－1)×(K－1)，且 Aij、Bij定義如下:

一般地，假設(shè) s=0，bR=0，則(12)式變形為:

令v2ln+1=Z2n+1eiβlh，其中 β 是非負(fù)整數(shù)，代入(15)式得:

則由(16)式可得:Z2n+2=G2n+2Z2n，其中

有

引理 2．1 假設(shè) s=0，bR=0，對于(12)式，若α∈，則

．由于

所以

根據(jù)引理1．1和引理1．2得:

故

如果 ζα/2＞η 成立，同理可得 G1≤1，則，從而得到了結(jié)論，證畢．

令

將(21)式代入(13)式可得:

經(jīng)過簡單推導(dǎo)，可得:

其中m、n是非負(fù)整數(shù)．

上式等價于下列方程:其中，1≤l≤K－1，n≤N=T/(2Δt)．

欲證明‖(I－A)(I+A)－1‖≤1，只需證明(24)式的增長因子的模小于1，這個證明方法和引理2．1中類似，這里不再贅述．同理可得‖(I+B)－1‖≤1，證畢．

定理2．3 對已經(jīng)得到的滿足初邊值條件(10)式的非對稱迭代格式(9)式，假設(shè) c(x，t)=c，d(x，t)=d是常函數(shù)且bR=0，則當(dāng)α∈時，是一致穩(wěn)定的，且成立，其中m、n是正整數(shù)．

證明 由引理2．1、引理2．2 和(23)式，可得綜合(21)、(25)和(26)式，可得反復(fù)運(yùn)用(27)式可得結(jié)論，證畢．

3 差分格式的誤差分析

本節(jié)討論非對稱迭代格式(9)的誤差估計．首先假設(shè) c(x，t)=c，d(x，t)=d 是常函數(shù)．

設(shè)umi=u(xi，tm)是滿足方程組(1)的真實(shí)解，則有:

利用泰勒展開公式，可得:

由于 1＜α＜2，將(29)式代入(28)式得 t2n+1時刻的截斷誤差:

同理，可得t2n+2時刻的截斷誤差為:

設(shè) Um=(u1m，um2，…，umK－1)，下面給出方程組(1)的精確解和計算解的誤差估計．

定理 3．1 設(shè)u(xi，tn)是方程組(1)的精確解，vn是其計算解，如果 c(x，t)=c，d(x，t)=d 都是常函數(shù)，且空間步長h和時間步長Δt都足夠小時，存在與h和 Δt均無關(guān)的正常數(shù) C且 m≤T/Δt，使得‖Vm－Um‖≤C(Δt2h－2(α－1)+Δt+h2)．

注1 由(29)式可知，本文給出的半隱有限差分格式(9)的局部截斷誤差為而定理4．1給出的其精確解與解析解的誤差估計是應(yīng)用了非對稱迭代格式后具有的特點(diǎn)，該定理證明過程繁瑣且與文獻(xiàn)［19］中相關(guān)定理證明類似，此處略去．

4 數(shù)值算例

例1 考慮如下的空間分?jǐn)?shù)階對流－擴(kuò)散方程

其中h1、h2是空間步長．

由定理4．1可知，方程組(1)的精確解和計算值之間的誤差估計為 O(Δt2h－2(α－1)+Δt+h2)．為了驗(yàn)證這個結(jié)果，則:

1) 取 Δt=Chα，此時收斂階應(yīng)為O(hα);

2)取Δt=Ch2，此時收斂階應(yīng)為O(h2)．對應(yīng)不同的α選取不同的C從而更好的驗(yàn)證結(jié)果，其中C是非負(fù)整數(shù)．

表1～4分別給出了在t=1．0時刻的離散L2誤差和離散L!誤差，以及空間方向?qū)?yīng)的收斂階．分?jǐn)?shù)階導(dǎo)數(shù) α 分別取1．6和1．8．顯然，給出的對稱迭代格式(9)在空間方向達(dá)到了理論分析的精度，甚至更高．

表1 例1在t=1時刻的誤差和收斂階Tab．1 Error and convergence rates for example 1

表3 例1在t=1時刻的誤差和收斂階Tab．3 Error and convergence rates for example 1

表4 例1在t=1時刻的誤差和收斂階Tab．4 Error and convergence rates for example 1

5 結(jié)束語

本文結(jié)合非對稱迭代技術(shù)給出了一個半隱式的有限差分格式，形式上是隱格式，但是通過偶數(shù)和奇數(shù)時間層選取不同的節(jié)點(diǎn)模板進(jìn)行顯式化計算，并且通過分析，空間方向精確解和計算解的誤差估計為 O(Δt2h－2(α－1)+Δt+h2)，通過選取合適的時間步長，收斂階在l2范數(shù)意義下可以達(dá)到二階，甚至更高．