Caputo導(dǎo)數(shù)的一種高精度算法及其在擴(kuò)散現(xiàn)象中的應(yīng)用

吳娟,張華,羅衛(wèi)華

(1. 湖南文理學(xué)院信息與現(xiàn)代教育技術(shù)中心,湖南常德,415000;2. 湖南文理學(xué)院數(shù)理學(xué)院,湖南常德,415000)

近年以來,分?jǐn)?shù)階導(dǎo)數(shù)作為一個(gè)新的數(shù)學(xué)工具,已經(jīng)在諸多科學(xué)研究領(lǐng)域得到了廣泛的應(yīng)用,特別是在高度復(fù)雜的介質(zhì)當(dāng)中,對(duì)相應(yīng)的擴(kuò)散行為作數(shù)值模擬時(shí),分?jǐn)?shù)階導(dǎo)數(shù)能夠提供更好的效果,相關(guān)結(jié)果可見文獻(xiàn)[1-4]。

作為分?jǐn)?shù)階導(dǎo)數(shù)的一個(gè)重要成員,Caputo 導(dǎo)數(shù)算子

在時(shí)間分?jǐn)?shù)階領(lǐng)域已經(jīng)得到了廣泛的應(yīng)用研究,這里 Γ(·)為Gamma 函數(shù)。而此類算子由于積分算子具有弱奇異性和非局部性,對(duì)其做數(shù)值離散是一個(gè)難題,這也是近年以來許多學(xué)者的重點(diǎn)研究對(duì)象。其中,高廣花等在文獻(xiàn)[5]中利用線性和二次的Lagrange 插值函數(shù)對(duì)式(1)進(jìn)行離散,得到的數(shù)值方法具有局部的3-β階精度,類似結(jié)果可見文獻(xiàn)[6-9]。

本文利用二次Lagrange 插值函數(shù)構(gòu)建式(1)的一種數(shù)值離散方法,給出所對(duì)應(yīng)的數(shù)值遞歸公式,理論分析顯示該方法可以達(dá)到全局的3-β階數(shù)值精度,且能保證良好的數(shù)值穩(wěn)定性。作為實(shí)際應(yīng)用,本文將該類數(shù)值方法應(yīng)用于時(shí)間分?jǐn)?shù)階擴(kuò)散方程的數(shù)值求解中,并用數(shù)值算例加以驗(yàn)證。

1 數(shù)值離散

在區(qū)間 [t1,t3]采用過t1,t2,t3的二次Lagrange 插值函數(shù),經(jīng)過計(jì)算和化簡可得t=t2處的離散格式為

其中A1,1=-2ρβ[1/(2-β)-1],A1,2=ρβ[1/(2-β)-1/2],δ1=ρβ[1/(2-β)-3/2],ρβ=1/(τβΓ(2-β))。

類似地,可以得到在t=t3處的離散格式為

其中A2,1=-22-βρβ[2/(2-β)-1],A2,2=2-βρβ[4/(2-β)-1],δ2=2-βρβ[4/(2-β)-3]。

對(duì)于區(qū)間 [t1,tk],(k=4,5,…,N+1):

(1)當(dāng)k=2p時(shí),將式(1)重寫為對(duì)于項(xiàng),采用過t1,t2,t3的二次Lagrange 插值; 對(duì)于項(xiàng),采用過t2j,t2j+1,t2j+2的二次Lagrange 插值。經(jīng)過計(jì)算和化簡可得

(2)當(dāng)k=2p+ 1時(shí),將式(1)重寫為對(duì)于項(xiàng),采用過t1,t2,t3的二次Lagrange 插值; 對(duì)于項(xiàng),采用過t2j+1,t2j+2,t2j+3的二次Lagrange 插值。經(jīng)過計(jì)算和化簡可得

綜上所述,Caputo 導(dǎo)數(shù)算子(1)在所有網(wǎng)格點(diǎn)處的近似值可以由式(2)～(5)逐個(gè)求出。此外,關(guān)于數(shù)值格式(2)～(5)的可解性和逼近精度,有如下結(jié)論。

定理1對(duì)于任意的右端函數(shù)f(t)和初始值u1=u(t1),線性方程組

有唯一的解(u2,u3,…,uN+1)T。

證明首先,從式(2)、(3)可得二元線性方程組此方程組的系數(shù)矩陣的行列式為從而,結(jié)合右端項(xiàng)f(2t),f(t3)和初始值1u,由該方程組可以唯一求出解23T(u,u)。

其次,注意到cj,s≠0,dj,s≠0,(j=4,5,…,N+1,s=1,2,3),由式(4)、(5)可知: 當(dāng)4≤k=2p≤N+1時(shí),uk可以由式算出。當(dāng)4k=2p+≤1≤N+1時(shí),ku可以由式算出。

因此,對(duì)于任意右端函數(shù)f(t)和初始值u1=u(t1),線性方程組(6)有唯一解。定理證畢。

定理2對(duì)于任意的β∈(0,1),u(t)∈C3[0,T],數(shù)值格式(2)～(5)的誤差總滿足

證明由Lagrange 插值的誤差理論知,對(duì)任意光滑函數(shù)u(t)和子區(qū)間[tk,tk+2],k=1,2,…,N-1,有u(t)-L[tk,tk+2](u(t))=u'''(ξk)(t-tk)(t-tk+1)(t-tk+2)/6,t∈[tk,tk+2],ξk∈(tk,tk+2)。

對(duì)于R(u(t2)),因?yàn)槎捎?/p>

對(duì)于R(u(tk)),4≤k=2p≤N+1,將誤差重寫為并且利用分部積分法,可得

對(duì)于R(u(tk)),3≤k=2p+1≤N+1,可將誤差重寫為類似于k=2p的情形,可得。定理證畢。

2 數(shù)值算法在時(shí)間分?jǐn)?shù)階擴(kuò)散方程中的應(yīng)用

作為應(yīng)用,本節(jié)利用上述數(shù)值離散格式(2)～(5)來求解如下時(shí)間分?jǐn)?shù)階擴(kuò)散方程。

該方程是工程領(lǐng)域中模擬高度復(fù)雜介質(zhì)下的流體擴(kuò)散所涉及的一個(gè)數(shù)學(xué)模型,時(shí)間分?jǐn)?shù)階導(dǎo)數(shù)的介入,可以有效緩解整數(shù)階導(dǎo)數(shù)算子所產(chǎn)生的重尾效應(yīng)。對(duì)于模型(8)的數(shù)值模擬已經(jīng)得到了一些計(jì)算數(shù)學(xué)研究者的廣泛關(guān)注,其中文獻(xiàn)[6,10-11]均取得了一些良好效果。

對(duì)于模型(8),時(shí)間方向用數(shù)值格式(2)～(5),空間方向用經(jīng)典的緊差分格式進(jìn)行數(shù)值求解。在下面的表1、2 當(dāng)中,N表示時(shí)間方向的網(wǎng)格數(shù),空間方向取定網(wǎng)格步長為1/128。為了體現(xiàn)收斂階,取通用的誤差計(jì)算公式這里E(N)表示網(wǎng)格數(shù)為N所對(duì)應(yīng)的誤差。

表1 例1 的數(shù)值結(jié)果

表2 例2 的數(shù)值結(jié)果

例1

例2

由表1、2 中的實(shí)驗(yàn)結(jié)果可知,數(shù)值格式(2)～(5)的收斂階均達(dá)到了O(τ3-β),與定理2 中的理論結(jié)果相吻合。

3 結(jié)論

本文利用二次Lagrange 插值函數(shù)構(gòu)造了Caputo 分?jǐn)?shù)階導(dǎo)數(shù)算子的一種數(shù)值離散方法,該方法的全局?jǐn)?shù)值精度為O(τ3-β)。為了調(diào)查其實(shí)際使用效果,此數(shù)值方法被應(yīng)用到了數(shù)值求解時(shí)間分?jǐn)?shù)階擴(kuò)散方程當(dāng)中,發(fā)現(xiàn)其數(shù)值效果與理論分析保持一致性。鑒于數(shù)值格式(2)～(5)的高精度性質(zhì),在其他類型的時(shí)間分?jǐn)?shù)階方程的數(shù)值模擬當(dāng)中,該方法具有很好的應(yīng)用潛力。