王 婉,張海湘,楊雪花
(湖南工業大學 理學院,湖南 株洲 412007)
分數微分方程現已被應用于越來越多的領域,如黏彈性材料[1-2]、控制理論[3]、金融市場中的期權定價模型[4-5]等。到目前為止,科研工作者們提出了很多解決具有弱奇異性分數階微分方程的方法,比如有限差分法[6-7]、有限元方法[8-9]、光譜方法[10],傅里葉變換法[11-12]等。對于分數階微分方程中的Caputo導數項,學者們采用不同方法近似處理。例如,Xu D.等[13]用L1離散Caputo導數項,用二階卷積求積公式近似Riemann-Liouville積分項,構造了緊差分格式,不僅給出了收斂性和穩定性的證明,且以數值算例驗證了理論分析結果。Guo J.等[14]用L1-2公式離散Caputo導數項,用二階有限差分法離散積分項,得知其空間收斂階和時間收斂階都達到二階,并證明了格式的穩定性與收斂性。Shen J.Y.等[15]提出H2N2數值微分公式(二次Hermite和Newton插值多項式的應用)來近似Caputo導數項,并建立了有限差分格式。為了增加計算效率,其利用指數和近似核t1-γ,并推導出一種快速差分格式。對于初始時間的弱奇異性也在等級網格中進行了討論。Jiang S.D.等[16]提出了分數擴散方程的快速L1差分格式,(0, 1)階的Caputo導數項被快速L1遞歸公式離散。Xu W.Y.等[17]考慮一種具有二階空間和時間精度的快速差分格式。用FL2-1σ公式近似時間卡普托導數,該公式使用了核指數和近似Caputo微分中出現的函數;通過離散能量方法證明了其無條件收斂性,最后通過數值例子驗證了該方案的數值精度和效率。
本文討論如下四階分數波方程的初邊值問題:
其初始條件和邊界條件分別如下:
式(1)~(3)中:Ω=(0,L)×(0,T];f(x,t)、u0(x)、ψ(x)為光滑函數;
本文中,C為一個正常數,在不同情況下可能具有不同的值。
對區間[0,L]作M等分,區間[0,T]作N等分,記h=L/M,τ=T/N,xj=jh, 0≤j≤M,tn=nτ, 0≤n≤N,其中h為空間步長,τ為時間步長。
引理1[16]對于給定的α∈(0, 1)、截止時間限制δ、誤差ε和最后時間T,有一個正整數Nexp,正點和相應的正權重,可得,其中
接下來介紹一種運用H2N2方法計算Caputo分數導數的快速算法。
令δ=τ/2,可得
其中,
引理2[15]設,γ∈(1, 2),則有,。
其中
引理3[17]若函數,則有
定解問題(1)~(3)可寫成如下形式:
定義網格函數
將離散緊算子A分別作用于式(6)和(7),由泰勒展開式、引理1及2,可得
式(8)(9)中:
注意:初邊值條件為
引理4由帶積分余項的泰勒展開式及不等式(10),可得
為了證明格式的收斂性,給出以下引理。
引理5[18]設u,v∈Vh,則
引理6[19]設u,v∈Vh,則有
引理7[20]設u∈Vh,可得
引理8[21]對于任意正整數p和函數G={G1,G2,…},有
其中,{aj}定義見式(4)。
引理9[19]對任意網格函數,有
引理10[22]設F(k)、g(k)為非負函數,且滿足
定理1設是問題(1)~(3)的解,是差分格式的解。為精確解,為數值解。則有
證明:令將式(8)(9)和(12)與(13)相減,得到如下誤差方程:
將算子A作用在式(13)中第二個式子,可得
利用引理2和引理3,則有
用式(16)減去式(15),可得
由引理6可知
將式(18)和(19)相加,可得
對于式(21)右端的第一項,利用引理9和young不等式,且考慮η0=0,可得
式(21)右端第二項,利用young不等式,得到
將式(22)~(25)代入式(21),有
令
因而式(26)可寫成
考慮式(10)和引理4可得
當m為任意值時,有
由式(14)和(17)知
引用引理7,有
定理1證畢。
應用快速緊差分格式計算下列定解問題:
其中精確解為
右端項為
不同步長時的最大誤差為
空間收斂階為
時間收斂階為
固定時間步長N=2 048,參數ε=10-12,γ取不同值時,差分格式的最大誤差以及相應的空間收斂階與CPU時間見表1。由表1中數據可知,上述問題的空間收斂階為4。

表1 N=2 048時最大誤差及相應的空間階、CPU時間Table 1 Maximum error, spatial order and CPU time with N=2 048
固定空間步數M=512,參數ε=10-12,表2中列出了不同γ取值和時間步長N下,計算得到的最大誤差,及其相應的時間收斂階與CPU時間。分析表2中的數據可以得知,時間收斂階為(3-γ)階,且所需CPU時間較少。

表2 M=512時最大誤差及時間階與CPU時間Table 2 Maximum error, time order and CPU time with M=512
固定γ值,圖1給出了在空間方向上的收斂階,圖2給出了在時間方向上的收斂階。由圖1和2可看出,實例的空間收斂階與時間收斂階與理論分析得到的結論是一致的。

圖1 當γ固定時的空間收斂階Fig.1 Spatial convergence orders with a fixed value of γ

圖2 當γ固定時的時間收斂階Fig.2 Temporal convergence orders with a fixed value of γ
本文研究了四階分數波動方程的快速緊差分格式,Caputo導數項用一種快速的H2N2方法來近似,通過使用降階法和離散能量法得到格式的收斂性。由數值算例可知,本文考慮的緊差分格式的收斂階為O(τ3-γ+ε+h4),且數值算例驗證了理論分析結果。該格式所需的CPU時間較短。