李玉山（甘肅政法學(xué)院 網(wǎng)絡(luò)空間安全學(xué)院，甘肅 蘭州，730070）

1 引 言

進(jìn)入21世紀(jì)以來(lái),隨著分?jǐn)?shù)階微積分在反常擴(kuò)散、非牛頓流體力學(xué)、多孔介質(zhì)力學(xué)、粘彈性材料力學(xué)、地球物理、生物醫(yī)學(xué)工程、經(jīng)濟(jì)學(xué)等諸多領(lǐng)域的建模及廣泛應(yīng)用[1-6],分?jǐn)?shù)階微積分的理論及應(yīng)用開(kāi)始受到越來(lái)越多的關(guān)注,相關(guān)的結(jié)果也陸續(xù)出現(xiàn),參見(jiàn)Podlubny, Kilbas等人的專著[7、8]以及中文專著[9、10]。分?jǐn)?shù)階反常擴(kuò)散方程是由經(jīng)典擴(kuò)散方程推廣而得到的。相比整數(shù)階擴(kuò)散方程，分?jǐn)?shù)階擴(kuò)散方程更適合描述反常擴(kuò)散現(xiàn)象，因?yàn)榉謹(jǐn)?shù)階導(dǎo)數(shù)能夠描述具有記憶和遺傳性質(zhì)的非均勻物質(zhì)[12、13]。

時(shí)間-空間分?jǐn)?shù)階擴(kuò)散方程

02α＜＜被運(yùn)用在反常擴(kuò)散模型中，通常,若考慮時(shí)間相關(guān)性或擴(kuò)散的記憶性，就得到時(shí)間分?jǐn)?shù)階擴(kuò)散方程；若考慮空間相關(guān)性或非局域性，則得到空間分?jǐn)?shù)階擴(kuò)散方程；如果既考慮時(shí)間相關(guān)性，又考慮空間整體性，則得到所謂的時(shí)間-空間分?jǐn)?shù)階擴(kuò)散方程。對(duì)于時(shí)間分?jǐn)?shù)階擴(kuò)散方程的數(shù)值求解，已經(jīng)有了很多的研究成果[15-20]。對(duì)于時(shí)間-空間擴(kuò)散方程的數(shù)值求解，相對(duì)來(lái)說(shuō)，研究的比較少，劉發(fā)旺等[21,22]給出了空間分?jǐn)?shù)階拉普拉斯算子的處理辦法，但是時(shí)間項(xiàng)依然是整數(shù)階的。

本文將探討一維時(shí)間-空間擴(kuò)散方程的數(shù)值求解，利用矩陣轉(zhuǎn)換技術(shù)處理空間分?jǐn)?shù)階拉普拉斯算子和有限差分法處理時(shí)間分?jǐn)?shù)階項(xiàng)，得到代數(shù)方程組，同時(shí)利用分離變量法得到解析解表達(dá)式，并且給出數(shù)值算例作比較。

2 問(wèn)題的提出

考慮如下的齊次Dirichlet邊界條件的時(shí)間-空間分?jǐn)?shù)階擴(kuò)散方程的初邊值問(wèn)題：

其中Γ為伽馬函數(shù)。

為階數(shù)為α(0 ＜α≤2)的空間分?jǐn)?shù)階拉普拉斯算子，由拉普拉斯算子-（）Δ的譜分解來(lái)定義，定義如下：

設(shè)

那么對(duì)于任意的fγ∈F，定義分?jǐn)?shù)階拉普拉斯算子為：

本文考慮給定f(x),p(t),g(x)和方程(1)-(3)，數(shù)值求解u(x,t)。

3 數(shù)值方法

利用矩陣轉(zhuǎn)換技巧求解齊次Dirichlet邊界條件的時(shí)間-空間分?jǐn)?shù)階擴(kuò)散方程的初邊值問(wèn)題。先考慮標(biāo)準(zhǔn)的一維擴(kuò)散方程：

取 正 整 數(shù) 1N＞ ,空 間 步 長(zhǎng) 1/hN= ,ixih=(0iN≤≤)，引進(jìn)標(biāo)準(zhǔn)的中心差分離散二階空間導(dǎo)數(shù)：

利用文獻(xiàn)[21,22]的結(jié)果，方程（14）-（16）可以離散為如下矩陣形式

取正整數(shù)1M＞,時(shí)間步長(zhǎng)(0jM≤≤)，用隱式差分格式離散時(shí)間分?jǐn)?shù)階項(xiàng)[20]：

將（19）寫為矩陣形式，代入（17）得：

利用Matlab軟件求解代數(shù)方程組（20），即可得到問(wèn)題（1）-（3）的近似解。

為了驗(yàn)證本文的方法，利用分離變量法，給出問(wèn)題（1）-（3）的解析解如下：

其中,Eβ β為廣義的Mittag-Leffler函數(shù)，定義如下：

在計(jì)算（21）的過(guò)程中，只計(jì)算無(wú)窮級(jí)數(shù)的前50項(xiàng)來(lái)近似得到u(x,t)，利用Matlab程序[23]計(jì)算Mittag-Leffler函數(shù)時(shí)取精度為10-6。

4 數(shù)值例子

例:取T=1,初始值u(x,0)=g(x)=x(1-x),f(x)=sin(x),p(t)=t,u(0,t)=u(1,t)=0,γ=1,M=N=100。

圖1 α =1.2, β =0.3時(shí)的數(shù)值解、解析解以及誤差圖

圖2 α =1.8, β =0.3時(shí)的數(shù)值解、解析解以及誤差圖

圖3 α =1.2, β =0.7時(shí)的數(shù)值解、解析解以及誤差圖

圖4 α =1.8, β =0.7時(shí)的數(shù)值解、解析解以及誤差圖

表1 數(shù)值解和精確解CPU占用時(shí)間對(duì)比 單位（s）

圖1-圖4給出了例1中α=1.2,1.8和β=0.3,0.7的數(shù)值結(jié)果，并與解析表達(dá)式（3.15）做了比較，可以看出，本文中的方法很有效，并且具有較高的精度。由表1可以看出，隨著時(shí)間和空間離散點(diǎn)的增加，解析解的CPU占用時(shí)間遠(yuǎn)高于本文提出的數(shù)值方法的CPU占用時(shí)間，這是由于利用解析表達(dá)式（3.15）計(jì)算時(shí)，牽涉到計(jì)算兩個(gè)無(wú)窮級(jí)數(shù)、數(shù)值積分以及Mittag-Leffler函數(shù)，在時(shí)效上有很大的缺陷，而本文提出的方法避免了這一缺陷，只需要在每一時(shí)間層解一個(gè)線性方程組，有較高的精度，且占用極少的CPU時(shí)間。

5 結(jié) 語(yǔ)

本文研究了一維具有齊次Dirichlet邊界條件的時(shí)間-空間分?jǐn)?shù)階擴(kuò)散方程的數(shù)值求解，采用矩陣轉(zhuǎn)換技巧求解且給出了數(shù)值例子，并和解析解做了比較。數(shù)值例子表明，1）本文的方法是很有效的，具有較高的計(jì)算精度，可以處理實(shí)際問(wèn)題；2）矩陣轉(zhuǎn)換技巧求解只需要求解線性方程組，而利用解析表達(dá)式求解計(jì)算比較復(fù)雜，且占用CPU時(shí)間較長(zhǎng)；3）利用矩陣轉(zhuǎn)換技巧可以推廣到第二類、第三類邊界條件以及非齊次邊界條件的問(wèn)題。