奇異函數分數階導數的Hadamard 有限部分積分表示形式

婁汝馨，廉 歡，王同科

（天津師范大學數學科學學院，天津 300387）

利用分數階微積分可以更準確地描述一些實際問題，分數階微積分具有時間記憶性和全局相關性，在流體力學、圖像處理等許多領域有廣泛的應用，因此對分數階微積分的研究在不斷地深入.目前已有多種分數階導數的定義，其中Riemann-Liouville（RL）和Caputo 分數階導數的應用最為廣泛，這2 種導數的定義都是由RL 分數階積分算子

與經典的微分算子相結合導出的.

定義1[1-2]給定實數α >0，記[α]為α 的整數部分，n = [α] + 1.對于定義在x∈（a，b]上的函數f（x），其α 階RL 分數階導數和Caputo 分數階導數分別為

式中Tn-1[f]（x）表示f（x）在x=a 的n-1 階Taylor 級數展開式.

關于分數階微積分的計算方法研究一直是該領域的熱點問題，主要包括直接方法和間接方法.關于直接方法，文獻[3]通過直接計算瑕積分，給出了一種去掉積分瑕點的方法，在此基礎上設計了計算RL 分數階積分和導數的數值算法；文獻[4]利用移位Chebyshev 多項式逼近RL 導數，進而給出分數階邊值問題的數值算法；文獻[5]利用譜近似方法基于Legendre、Chebyshev和Jacobi 多項式的三項遞推公式逼近分數階積分和導數，并通過數值算例說明了方法的有效性；文獻[6]介紹了分數階導數的各類算法，如G 算法[7]和R 算法[8].間接方法主要通過先將分數階導數變形為Hadamard 有限部分（HFP）積分形式，然后再設計算法離散HFP 積分.

定義2[9]設f（x）∈Cn[s，r]，且n-λ>1，則HFP積分定義為（本文中用H ∫表示Hadamard 有限部分積分）

式中g（t）∈Cn[s，r]且使上式極限存在.

關于RL 分數階導數與HFP 積分的關系，有以下定理.

定理1[2]給定實數α>0，令n=[α]+1.設f（x）∈Cn[a，b]，則有

定理1 的結論非常重要，許多研究由此設計了計算RL 分數階導數的數值算法，文獻[10-11]利用定理1 的轉換設計算法求解分數階微分方程，得到了滿意的計算結果.定理1 要求f（x）∈Cn[s，r]，這是一個比較高的光滑性要求，而在求解分數階微分方程時，方程的解在初始點往往沒有如此高的光滑性，因此有必要對此結果進行推廣.

本文首先證明了當f（x）在x = a 處包含奇性時，HFP 積分形式仍然成立；然后基于該結論求得分數階導數在x=a 處的Psi 級數展開式；最后結合函數在奇點的Psi 級數展開式設計Chebyshev 譜逼近方法，用于給定區間上函數分數階導數的高精度計算.

1 分數階導數的Hadamard 有限部分積分表示形式

首先討論RL 分數階導數，證明當f（x）在初始點奇異時，其HFP 積分形式仍然成立.

定理2給定實數α>0，令n=[α]+1.設f（x）∈Cn（a，b]，則有

證明當f（x）在x = a 處奇異時，為處理在積分瑕點t=a 處的Hadamard 積分，引入一個充分小的正數δ，滿足a

式中：

由f（t）∈Cn[a+δ，b]，對f（t）在t=a+δ 處做Taylor 展開，得

其中余項使用積分表示形式，即

將式（6）代入式（5）中的I2，得

由HFP 積分定義可得

由Taylor 級數展開式的余項可得

將以上2 式代入式（7），得

再將式（8）代入式（5），得

下面說明RL 分數階導數定義（式（1））可以寫為式（9）的形式.由式（1）得

式中：

當x>a+δ，t∈[a，a+δ]時，（x-t）n-α-1關于x 充分光滑，對其求n 階導數，得

由定積分的可導性可得

對上式逐次進行分部積分和求導，可得

將式（11）—式（12）代入式（10），得

對比式（9）和式（13），可知式（4）成立.定理得證.

定理2 中f（x）∈Cn（a，b]，說明定理2 對f（x）在x =a 處奇異時依然成立.由定理2 可以導出f（x）的RL 分數階導數在x=a 處的級數展開式.

定理3給定實數α>0，令n=[α]+1.設f（x）∈Cn（a，b]，f（x）在x=a 處有Psi 級數展開式

式中：ui、μij為非負整數；βi為實數，且滿足-1<β0<β1<…→+∞.則f（x）的RL 分數階導數在x=a 處可以展開為

證明將式（14）代入式（4），得

對于上式中的HFP 積分，令t-a=（x-a）s，則有

將上式代入式（16），可得式（15）成立.定理得證.

注：式（14）和式（15）為函數在x=a 處的一般形式的級數展開式，包含實數次冪和對數多項式.這類級數展開式一般稱為Psi 級數[12]或Puiseux 級數[13]，是Taylor 級數和Laurent 級數的推廣.

式（14）和式（15）均為Psi 級數展開式，函數的Psi級數展開式的系數目前還沒有一般的公式表示，對于具體函數，可以使用數學軟件直接求出展開式中的冪指數和各個對數項的系數，具體算例可參見文獻[13].

利用定義1 中RL 分數階導數與Caputo 分數階導數之間的關系（式（2）），可以得到Caputo 分數階導數的HFP 積分表示形式及其在初始點的級數展開式.

推論給定實數α >0，令n=[α] + 1.設f（x）∈Cn-1[a，b]∩Cn（a，b]，則f（x）的Caputo 分數階導數的HFP 積分表示形式為

式中Tn-1[f]（t）為f（x）在x=a 處的n-1 階Taylor 級數展開式.進一步，若f（x）在x=a 處有Psi 級數展開式

式中：ui、μij為非負整數；βn為實數，且滿足n-1<βn<βn+1<…→+∞.則f（x）的Caputo 分數階導數在x=a 處可以展開為

2 分數階導數的Chebyshev 譜逼近與數值算例

當函數在區間上充分光滑時，譜逼近算法[5]有很高的計算精度.本節針對函數在左端點奇異的情形，基于RL 分數階導數的HFP 積分表示設計一種奇點分離的Chebyshev 譜逼近方法.

設函數f（x）在左端點x=a 處弱奇異，即f（x）的Psi 級數展開式（14）中β0>-1.設δ>0，滿足a

I1為弱奇異積分，可以使用文獻[13]設計的改進復合Gauss-Legendre 算法有效計算.I2為Hadamard 有限部分積分，當t≥a+δ 時，f（x）充分光滑，下面設計Chebyshev 譜逼近法離散該積分.

在區間[a+δ，b]上離散積分I2，做變換

則t、x∈[a+δ，b]轉化為s、y∈[-1，1].記g（s）=f（t），σ=，則式（19）轉化為

g（s）在區間[-1，1]上充分光滑，對其做Chebyshev插值，插值節點選為m 次Chebyshev 多項式Tm（s）=cos（m arccos s）的零點si=-cos θi，其中i=1，2，…，m，則g（s）的Chebyshev 插值形式為

式（22）中Ij（α，y）為超奇異積分.文獻[5]針對Ij（α，y）弱奇異的情形使用Chebyshev 多項式的三項遞推公式導出了計算Ij（α，y）的快速遞推算法.經檢驗，該遞推算法對于這里的超奇異積分仍成立，從而可以快速求得f（x）的RL 分數階導數在一些點的高精度計算值.

在求解分數階微分方程時，可以通過迭代的方法求出解在左端點的Psi 級數展開式，分離出解的奇性，在正則區間上利用Chebyshev 譜逼近方法則可以得到高精度的數值解.

算例當α=1/2 和3/2 時，求（fx）=K（0）的RL 分數階導數，其中K0（z）為零階第二類修正的Bessel 函數.

首先給出（fx）=K（0）在x=0 處的級數展開式，即

由展開式可知f（x）在x=0 處對數奇異，其RL 分數階導數可由數學軟件求出，其形式復雜，這里不再給出.利用定理3，當α=1/2、3/2 時，直接計算可得和在x = 0 處的有限項級數展開式S1（x）和S2（x）分別為

利用定義1 的式（1）直接求（fx）=K（0）的1/2 階和3/2 階RL 分數階導數，并對其在x=0 處做級數展開，通過對比可以發現與利用定理3 的結果完全一樣.由S1（x）和S2（x）的表達式可知，f（x）的1/2 階和3/2階RL 分數階導數在x =0 處奇異，當x>0 時，這2個分數階導數均存在.

利用Chebyshev 譜逼近法計算f（x）的RL 分數階導數，計算其與精確導數值的絕對誤差，并與級數展開式法的絕對誤差進行比較，結果見圖1.

圖1 對數尺度下級數展開式和Chebyshev 譜逼近法與精確導數值間的絕對誤差Fig.1 Absolute errors of series expansion and Chebyshev spectral approximation compared with precise derivative on logarithmic scale

圖1（a）和（b）分別給出了當α=1/2 和3/2 時2 種方法的誤差，其中：Eup 為S1（x）（S2（x））與精確導數值的絕對誤差對數圖形；CCM（m，δ）為Chebyshev 譜逼近法與精確導數值的絕對誤差對數圖形，取m=80，當δ=0時，CCM（80，0）為在整個區間[0，40]上做插值得到的結果，當δ=1 時，CCM（80，1）為在[0，1]上直接計算弱奇異積分，在[1，40]上做插值得到的結果.由圖1（a）可見，級數展開式S1（x）具有高精度的局部逼近性，其在零點附近的精度達到了10-14，但隨著x 的增大，其逼近精度逐漸降低.需要說明的是，S1（x）在x=0 處奇異，由于舍入誤差的影響，其在x=0 附近的誤差反倒比較小.利用Chebyshev 譜逼近法在全區間上做插值，其在x = 0 附近以及其他位置的精度都很低，只有10-1～10-4量級.這說明對于端點奇異的函數，全區間的Chebyshev 譜逼近法計算精度很低，不具有實用性.利用Chebyshev 譜逼近法在[1，40]上做插值，其在該區間的精度可以達到10-13，這說明級數展開式法和Chebyshev 譜逼近法聯合使用可以給出分數階導數在整個區間上的高精度近似值.由圖1（b）可見，α=3/2與α=1/2 的情況一致.

本算例還可以使用數值積分方法得到高精度的計算值.表1 給出當α=1/2 和α=3/2 時利用HFP 積分方法[15]（設置精度為10-16）計算的函數在一些點的分數階導數及其誤差，其中hα（x）為分數階導數的計算值，eα（x）= |0Dxαf（x）- hα（x）| 為絕對誤差.由表1 可見，數值積分方法的誤差在全區間上分布比較均勻.因此，使用分數階導數的Hadamard 有限部分積分表示形式，即使對于奇異函數也能夠得到高精度的計算結果.HFP 積分方法雖然計算精度很高，但在求解分數階微分方程時將面臨很大的困難，此時，奇點分離的Chebyshev 譜逼近法具有明顯的優勢.

表1 HFP 積分方法計算的函數在一些點的RL 導數及絕對誤差Tab.1 RL derivative values and absolute errors of the function at some points calculated by HFP integration method

3 結語

本文對于在初始點代數和對數奇異的函數，證明了其RL 和Caputo 分數階導數可以用Hadamard 有限部分積分表示，并由此導出了分數階導數在初始點的Psi 級數展開式，利用Psi 級數展開式在奇點附近的局部逼近性，結合Chebyshev 譜逼近方法，給出整個區間上分數階導數的高精度計算值.