劉建平,楊璐嘉,毛學志
(河北科技師范學院 數(shù)學與信息科技學院,河北 秦皇島 066004)
最近三十年,通過越來越多研究者的努力,確定了在動力學過程中的確具有分數(shù)階性態(tài)[1].分數(shù)階導數(shù)與整數(shù)階導數(shù)相比,它們能夠更加精確地模擬自然物理過程和動態(tài)系統(tǒng),它還能模擬物質(zhì)的記憶和遺傳性[3],因此分數(shù)階微分方程隨之成為一種解決復雜數(shù)學問題不可替代的數(shù)學工具,這就使許多學者開始在各個領(lǐng)域進行分數(shù)階微分方程的研究[2].目前,研究者在求解分數(shù)階微分方程的解析解的過程中,往往通過構(gòu)造特殊函數(shù),但該特殊函數(shù)構(gòu)造十分復雜,使得求解析解也變得很困難.所以,針對求解分數(shù)階微分方程的數(shù)值解的求解方法的研究就顯得更加重要了.
分數(shù)階擴散方程是在實際工程中有很廣泛的應用的一類方程,一些學者針對該類方程的求解方法展開了研究.文獻[4]與文獻[5]分別給出了變時間分數(shù)階擴散方程的隱式差分近似和時間分數(shù)階反應擴散方程的全解.文獻[6]與文獻[7]分別給出了變分數(shù)階非線性擴散方程的一個顯式有限差分格式和提出變階反常次擴散方程的兩種值格式;文獻[8]與文獻[9]分別給出了一個帶非線性源項的變分數(shù)階對流-擴散方程的顯式和一類變分數(shù)微分算子模型的擴散曲線.
筆者考慮到切比雪夫多項式具有良好的正交性質(zhì),展開公式中的冪函數(shù)數(shù)項更容易進行變分數(shù)階微分的計算.因此,從根本上具備了進行函數(shù)逼近處理的基礎,具有形成算子矩陣的條件.基于以上考慮,筆者探討利用移位切比雪夫多項式求解變時間分數(shù)階擴散方程的數(shù)值求法.考慮的變時間分數(shù)階擴散方程模型[10]如下:
初值和邊界條件為u(x,0)=g(x) 0≤x≤1;u(0,t)=h(t)0≤h≤1,其中u(x,t)區(qū)間[0,1]×[0,1]上的任意平方可積函數(shù),0<q(x,t)<1,f(x,t)已知,u(x,t)未知表示Caputo類型的變時間分數(shù)階導數(shù).
定義1[11]Caputo類型的變分數(shù)階微分定義
根據(jù)定義1,當u(t)=tn時,得到如下的計算結(jié)果
其中 0<q(x,t)<1.
本文利用了Caputo類型變分數(shù)階微分一些常用性質(zhì).
性質(zhì)1任何常數(shù)求α(t)階導都為零,即Cα(t)=0.
性質(zhì)2對于任意常數(shù)λ,μ,Caputo變分數(shù)階微分算子都有
性質(zhì)3可交換及疊加性:
切比雪夫多項式[12]定義區(qū)間為[-1,1],一般表現(xiàn)形式為
其在[-1,1]上遞推關(guān)系為
為了在x∈[0,1]上使用該多項式,作變量替換z=2x-1,得到移位切比雪夫多項式.用Hi(x)來表示,定義如下
此時移位切比雪夫多項式通項可表示為
函數(shù)y(x)在[0,1]上是可積的,通過移位切比雪夫多項式的變形,有
在實踐應用中只需要應用前(m+1)項
其中多項式的系數(shù)向量CT和向量函數(shù)Φ(x)如下:
式(4)是根據(jù)移位后的切比雪夫多項式(3)得到的y(x)函數(shù)近似表示.
對于任意二元函數(shù)u(x,t)∈L2([0,1]×[0,1]),都可以通過移位切比雪夫多項式進行函數(shù)的近似.一般地,也只對前n+1項作考慮
ui,j(i=1,2,…n;j=0,1,…,n)表示移位切比雪夫多項式逼近二元函數(shù)的待求系數(shù),如下
矩陣U用內(nèi)積表示為U=Q-1〈Φ(x),〈Φ(t),u(x,t)〉〉Q-1.Q,Φ(x),Φ(t)分別為
本節(jié)將推導變分數(shù)階矩陣算子.首先將式(5)表示成為矩陣形式(7)
其中 Tn(x)=[1,x,x2,…,xn]T,
稱為移位切比雪夫多項式的系數(shù)矩陣.矩陣A為一個上三角形且對角線沒有零元素,故可逆,從而有Tn(x)=A-1Φ(x).
對Φ(t)一階求導,有
其中D為(n+1)×(n+1)階的矩陣,稱為移位切比雪夫多項式的一階微分算子矩陣.根據(jù)式(7)可得
定義如下形式(n+1)×n階矩陣V(n+1)×n和n維列向量Tn*(t)
由于一元函數(shù)u(t)=CTΦ(t),則可將其導數(shù)u'(t)轉(zhuǎn)化為矩陣形式,即
對式(9)求二次導,可得Φ"(t)=D2Φ(t).針對任意的二元函數(shù)u(x,t),對x求二階偏導都可以得到如下式子:
下面將推導Φ(t)的q(x,t)階微分矩陣算子
其中
稱為移位切比雪夫多項式的q(x,t)階微分矩陣算子.
利用公式(14),將Dq(x,t)tu(x,t)表示為矩陣形式,可得
將式(13)和式(15)代入方程(1)中,得到如下形式:
例1 解如下變時間分數(shù)階擴散方程
該方程的精確解是u(x,t)=x2(1-x)(1+t3).
故n=2 時,數(shù)值解為u1(x,t)=ΦT(x)U1Φ(x),其中,U1如上所示.
故n=3時,該方程的數(shù)值解u2(x,t)=ΦT(x)U2Φ(x),其中U2如上所示.
分別取t=1/4,t=3/4時,在部分取值點處,精確解和數(shù)值解的絕對誤差見表1-2.
表1 t=1/4數(shù)值解與精確解的絕對誤差
表2 t=3/4數(shù)值解與精確解的絕對誤差
由上表可以觀察到當取n=2時,數(shù)值解與精確解的誤差大.而當n=3,4時,數(shù)值解與精確解的誤差接近于10-15.表1-2表明,當合理設置截斷項n,本文提出的方法是可行且有效.
本文在Caputo類型的變分數(shù)階微分定義下,根據(jù)變分數(shù)階微分的性質(zhì),推導出了移位切比雪夫多項式的變分數(shù)算子矩陣,結(jié)合配點法,將求解變時間分數(shù)階擴散方程問題轉(zhuǎn)化為求解代數(shù)方程組,并結(jié)合算例說明了該方法可行性及有效性.本文的研究為進一步探討變分數(shù)階微分方程的數(shù)值計算方法奠定了一定的理論基礎,具有一定的工程實用價值.另外,本文針對該方法的計算誤差、收斂性分析等問題還存在欠缺,這也是筆者下一步努力的方向.