分数阶导数双边空间微分方程的显式差分解法 Finite Difference Methods for Space-Time Fractional Two-Sided Space Partial Differential Equations

doi:10.12677/orf.2012.21001

设为首页加入收藏期刊导航网站地图

期刊菜单

文章导航

Operations Research and Fuzziolgy 运筹与模糊学, 2012, 2, 1-7

http://dx.doi.org/10.12677/orf.2012.21001 Published Online February 2012 (http://www.hanspub.org/journal/orf)

Finite Difference Methods for Space-Time Fractional

Two-Sided Space Partial Differential Equations*

Yang Zhang, Ningping Li, Lu Chen

School of Mathematical Science, Nankai University, Tianjin

Email: zhangyang@nankai.edu.cn

Received: Dec. 4th, 2011; revised: Jan. 6th, 2012; accepted: Jan. 8th, 20 12

Abstract: Fractional order differential equations are generalizations of classical differential equations. Now,

they are widely used in the fields of physics, information; finance and others. In this paper, an explicit finite

difference method for space-time fractional two-sided space partial differential equations is established. Be-

sides, the stability and co nvergence order are analyzed.

Keywords: Fractional Derivative; Explicit Methods; Stability; Convergence; Error Estimates

分数阶导数双边空间微分方程的显式差分解法*

张阳，李宁平，陈璐

南开大学数学科学学院，天津

Email: zhangyang@nankai.edu.cn

收稿日期：2011 年12 月4日；修回日期：2012 年1月6日；录用日期：2012 年1月8日

摘要：分数阶微分方程作为广义的微分方程，被广泛地应用于物理，信息处理，金融等领域。本文

给出了数值求解时间空间分数阶导数的双边空间微分方程的一种显式差分格式，并对其稳定性和收敛

性进行了理论分析。

关键词：分数阶导数；显格式；稳定性；收敛性；误差估计

1. 引言

分数阶导数微分方程，是传统的整数阶微分方程的推广。相比整数阶导数微分方程，分数阶导数微分方程

能够更好拟合某些自然物理过程和动态系统过程，近年来被广泛应用于反常扩散，粘弹性本构建模，信号处理，

控制，流体力学，图像处理，软物质研究，金融等领域[1-8]。近年来，分数阶导数微积分的研究与兴起，伴随着

数字计算机计算技术的提高而迅速发展。双边空间分数阶对流–扩散方程是一种反常扩散，在地下水溶质运移

方面，它可描述含水层中溶质运移过程中的反常扩散现象。关于这类问题的数值解法 Meerschaert 等人分别对单

边对流–扩散方程和双边扩散方程，进行差分求解，文[9]针对变系数反应扩散方程

 









,,,

crtcrt crt

vrrf rt

trr





 



*资助信息：该工作由中国国家自然科学基金资助(11171168，11071132)和中国高校博士点基金资助(20100031110002)。

张阳等  分数阶导数双边空间微分方程的显式差分解法

建立显式和隐式两种差分解法，文[10,11]建立了双边空间分数阶导数微分方程

 









,,,

uxtuxtuxt

cxtcsxt

txx













数值解法并进行了稳定性分析，文[12]针对二维分数阶扩散方程











,,,, ,,

uxytuxyt uxyt

dxy exyqxyt

txx



 





建立了交替方向隐式差分解法；夏源、吴吉春在[13]中对一维时间分数阶对流–弥散方程











,,cxt cxtcxt





 



给出了差分算法；刘发旺等人用另一种方法对这个问题进行离散计算[14]。周璐莹等在[15]中针对二维分数阶对

流–弥散方程

  

,,,, ,,,,,,cxytcxyt cxytcxytcxyt

xytx







 









建立了数值解法。[16]则将时间空间分数阶导数微分方程

 







 

,,,

uxtuxt uxt

bxaxcxuxt qxt

txx





 



的隐式差分解法进行了研究。目前，关于带时间分数阶导数的双边空间分数阶导数微分方程的数值解法尚未见

到，本文正是针对该类问题建立了相应的显式数值解法并进行理论分析，涵盖了[10]和[13]的类型，且采用两种

理论方法分析，较已有的方法更为完善。采用显式差分解法可以减少计算量，但理论分析较为复杂，从而本文

得到的结论丰富了这一领域的研究成果。

本文结构如下：第二节针对研究的对象进行了描述，并给出了两个常见的引理；第三节建立了所描述方程

的显式差分解法并进行了方法的稳定性分析；第四节给出了方法的收敛性分析和收敛阶的估计；第五节总结了

结论与展望。

2. 问题描述

本文针对带分数阶时间空间导数的双边空间微分方程建立显式差分解法并进行理论分析，研究如下方程：

















,,,,

uxtuxtuxt uxt

kxtcxtcs xt

txxx







 



(1)

这里，。其中

LxR 0tT 01



，01





，12





，





,kxt，





,cxt

，均为非负有界函数，

并且假定初始条件：，以及边界条件



,cxt





,0t



f xuxLxR









,,Rt 0ux Ltux



 。

方程(1)右端的左边(+)以及右边(–)的分数阶导数为 Riemann-Liouville 分数阶导数。按如下方式定义[10]：

 







d1d d

fx f

















 

 (2)

和

 



dd.

nnR

fx f

nxx

















 

 (3)

这里是正整数，并且满足n1nn



 。空间分数阶导数





,uxt





是



次Riemman-Liouville分数阶导数，

张阳等  分数阶导数双边空间微分方程的显式差分解法

定义为：









(,) 1d, 01.

xut

uxt















 



 (4)

方程(1)左端的



,uxt





由Caputo 分数阶导数给出定义[17]：







1d,

uxt ux















  (5)

这里是 Gamma 函数。

()

根据 Grünwald 估计，为了得到稳定的数值方法[9]，定义左右两边分数阶导数分别按如下形式近似(Shifted

Grünwald)：





lim 1

fx kh







 













 (6)

和



d()1

lim 1

fx kh







 













 (7)

这里

，

均为正整数，





，Rx



，











 







0,1,2

。

引理 1. 上式定义的 k

满足如下条件：

1) ，2)

01g10g



 ，3)











k, >2

 



 

。

对于引理 1，结论 1)、2)是显然的，下面只说明结论 3)。

因为12



，所以 0



，10



 ，0k





，。从而有 2k

  















  

11 121

111

12 1

kkk

gkk



 

 





 

 

 



，



0

引理成立。对于  ，根据二项式定理，有











 1z









1z。令 ，有。由引理 1知









10, 1

, 1,2,,



 

m (8)

3. 显格式及其稳定性分析

先建立方程(1)在前文所述初边值条件下的显格式并通过圆盘定理证明该差分格式是无条件稳定的，并且还

将采用另一种方法来证明格式的稳定性及收敛性，并给出收敛阶的估计。

引入步长 h，时间步长



，则有 n

, >0; ,0

xLjhh tn



 。其中 ,0

TRL







。对于





uxt





，



uxt





用shifted Grünwald形式，有如下表达式[9]：



,11, ,

nki n

uxt

ux khth















(9)

张阳等  分数阶导数双边空间微分方程的显式差分解法



,11, ,

nki n

uxt

ux khth















(10)



,1,

ki n

uxt ,

ux khth













(11)

其中



  

(2)

1, 1,1,2,

gg k

 



 



易验证



具有类似 k

的性质：

且

 



22 2

1,0,0 2

ggg k



 











2,3, ,km。



对于时间分数阶导数有如下估计[17]：















 

in i

ij ij

ninj inj

uxt ux

uxt uxt

uxtuxt jj





















 

















 







 





 



















(12)

故可对方程(1)建立如下的显格式：

 



111

1,1 ,1

000

nj njn

niim

nnnnn

ii i

ikki ikki iki

jkkk

uu k

jjgugcu gcu









 













 





 







(13)

易知，该方法是相容的，其舍入误差为





。记



cj j





，则有

引理 2. 上面所定义的



cj j





满足以下条件

1) , 2) , 3)

0, 0,1,

cj 1,0,1,

cc j







1+1

22+1

jjk

cc c











，即

+ 1







 







对于 1)，3) 是显然的，2) 利用单调性也可以得到。

令











22 2

nn n

ii i

iii

cc b

rhhh

















 



则当时，格式(13)为 0n



100000

02111 20111

1+()

imi n

iiiiiiiiiiiikiikkiik

uggurg gurggugugus



















(14)

当时，格式(13)为 0 n













02111 20

111

331

nn n

iiii iiiiiii

imin

nnnjn

kiikkiikjinii

kk j

uggurg gurggu

gugudu cus



 

 



 





 



 

 



(15)

因此格式(14)、(15)有如下等价形式

100

1110

nnn n

UHU S

UHUdU dUcUSn









 









 (16)

张阳等  分数阶导数双边空间微分方程的显式差分解法

其中矩阵 ij

，有如下表达式 1,1;1, ,1;imj m 

111

iii

iji i

iji

ji i

rg gji

gji

















 

 







, 当

，当













其中，，；

00 1H00

H1, 2,,jm1



，0



，0,1, ,1jm



。

定理 1. 当01



，12





时，对于方程(1)所建立的显格式(16)在满足如下条件：













 

1,,

22maxmaxmax

maxmaxmax1 0

cc b











 

  







 

 

 



(17)

时是稳定的。

证明. 根据 Gerschgorin 定理，以及引理 1，知矩阵 H的特征值所在中心为









22 22

iiii iii i

rgr









 H



，

半径

0,0, 10, 1

() ()()

mimi

iikikikiiiiii

kkikkkk

Hggrrggr.









 





根据(8)式，有





22 2222 1,

ii i



 



 

由(17)式知，，即 H的特征值不超过1。



22 221,

ii i





 

对于时的情形，由于0n0

的特征值



满足





122 1

ii i

 





,r0

，又因为

故













122 221

i i





 2



的特征值亦不超过 1，从而稳定性得证。

下面用另一种方法给出关于最大模的稳定性分析。假设





1, 2,,1;1, 2,,

uimjK



，是 (15)所得的近似

值，误差(扰动)为





1, 2,,1;1, 2,,

jjj

iii

uuim jK



 



，并且满足：

当时， 0n



101000

021112011

22 ;

imi

iiiiiiii iiiikiikkiik

ggrgg rgggg



  









 





当时， 0n











0211120

111

331

nn n

iiii iiiiiii

imin

nnnj

ki ikki ikjini

kk j

ggrgg rgg

ggdc







 



 





 



 

 



令，。

12 1

,,,T

kkk









0,1,2, ,kK

定理 2. 在条件(16)下，有 0,0,1,2,

EEk



n；这里11

max

E1n







 

。

证明. 用数学归纳法。

1) 当时，令0n1

max



 

，由(16)知





22 0

ll l







 ，故有



22 0

ll l







 。再

张阳等  分数阶导数双边空间微分方程的显式差分解法

由引理 1，则有



10 0

021 1120

11 11

00 0

33 00

llllll llllll

lml lml

lk lklk lklklk

kk kk

ggrgg rgg

gg ggE



 

 



 

 

 





 







由引理 1-3)知，，故







1







10





，即 0k



。

2) 假设 0,0,1,2,

EEn



k成立，考虑 1nk



的情形。

令11

max

E1k





 





，类似于 1)，有



02111 201

111 111

1111

331 001

(2 2)

kkk k

llllllllllll

lmlklmlk

kkkj

ls lslslsjlkllklkjk

ss jkkj

ggrgg rgg

ggdc ggdc





 



 





 0



 









 





 

 

100









定理得证。

由上述定理知格式(16)条件稳定。

4. 显格式的收敛性分析与收敛阶估计

下面进行收敛性分析。令是方程(1)在相应初边值条件下于点





,1,2,,1;1,2,

uxt imkK









t的

精确解。定义



,1,2,,1;1,2,

iiki

euxt uimkK并且 12 1

,,, T

kkk k

Yeee













，则有如下表达式







02111 20

11 0

111

331

22( )

,1,

kk k

iiiiiii iiii i

imin

kkkj

siissiisjiki

ssj

eggerg gergge

gegedeck



 

 



 



 



 



 











K

1, 2,,1;0,1, 2,imk。记为局部截断误差，则存在常数，使得

R0C



1, 1,2,,1;0,1,2,.

RC himk





 K

定理 3. 当满足条件(15)时，





11kk

YCc h0,1,2, 1kK







，



。这里11

max

 

，C为正常数。

证明. 用数学归纳法。

1) 当时，令0k1

max

 

1

e，根据引理 1，则有





10 1

0ll

eRCc h











，即



YCc







h。

2) 假设





11,1,2,

YCc hn







k；成立，由(17)知





22 0

ll l







 ，考虑的情形。

令

1nk

max



 





，利用引理 1，引理2，可以得到









02111 20

11 1

33 1

111

001

kk k

lllllll l llll

lml k

kkjk

ls lsjlkll

ss j

lmlk

lkl kjkkk

kkj

eggerg gergge

geged eceR

lsl s

ggdcCchCc





 



 



 

 



 

 









 

 





 

















张阳等  分数阶导数双边空间微分方程的显式差分解法

故



YCc











h，定理得证。

根据定理 3，并注意到



lim limlim1

kk k

ckk

kkk



 









 





 







，因此存在常数 C，使得





11kk

YCc h









，因

为kT



，进而



YCT















，从而有下述定理

定理 4. 假设是方程(1)在相应初边值条件下在点





,1,2,,1;1,2,

uxt imkK



t处的精确解。

是格式(16)所得到的数值解，则在满足条件(17)时，存在正常数，使得







,,1,2,, 1;

ik i

uxtuCh imkK



 1,2,.

5. 结论与展望

本文将[10]中讨论的双边空间分数阶导数微分方程的时间导数推广到分数阶，建立了一种显式差分求解格式

并进行了相应的理论分析，得到了丰满的一阶时间和空间的收敛阶逼近。所得结论亦可推广到方程右端含项及

更为复杂的情形，所建立的方法也不难推广到高维情形，至于高维问题如何进行算子分裂尚有待进一步探讨，

对于所讨论方程隐格式的建立与理论分析作者将另文讨论。最后，作者对审稿人对稿件提出的修改意见表示感

谢。

参考文献 (References)

[1] R. Metzler, J. Klafter. The restaurant at the end of the random walk: recent developments in the description of anomalous transport by

fractional dynamics. Journal of Physics, 2004, A37: 161-208.

[2] Z. Deng, V. P. Singh and L. Bengtsson. Numerical solution of fractional advection-dispersion equation. Journal of Hydraulic Engineering,

2004, 130(5): 422-431.

[3] R. Gorenflo, F. Mainardi, E. Scalas and M. Raberto. Fractional calculus and continuous-time finance III: The diffusion limit. In: Kolhmann, S.

Tang, Eds., Mathematical Finance, Basel: Birkhauser Verlag, 2001: 171-180.

[4] 苏丽娟, 王文洽. 双边分数阶对流–扩散方程的一种有限差分解法[J]. 山东大学学报(理学版), 2009, 10: 29-32.

[5] F. Liu, P. Zhuang, V. Anh, I. Turner and K. Burrage. Stability and convergence of the difference methods for the space-time fractional

advection-diffusion equation. Applied Mathematics and Computation, 2007, 191(1): 12-21.

[6] 周激流, 蒲亦非, 廖科. 分数阶微积分原理及其在现代信号分析与处理中的应用[M]. 北京: 科学出版社, 2010.

[7] R. Hilfer. Application of fractional calculus in physics. Singapore, New Jersey, London and Hong Kong: World Scientific Publication

Company, 2000.

[8] A. A. Kibas, H. M. Srivastava and J. J. Trujillo. Theory and application of fractional differential equations. Amsterdam: Elser vier, 2006.

[9] M. M. Meerschaert, C. Tadjeran. Finite difference approximations for fractional advection-dispersion flow equations. Journal of Computational

& Applied Mathematics, 2004, 172(1): 65-77.

[10] M. M. Meerschaert, C. Tadjeran. Finite difference approximations for two-sided space-fractional partial differential equations. Applied Nu-

merical Mathematics, 2006, 56(1): 80-90.

[11] V. K. Tuan, R. Gorenflo. Extrapolation to the limit for numerical fractional differentiation. Zeitschrift für Angewandte Mathematik und Me-

chanik, 1995, 75: 646-648.

[12] M. M. Meerschaert, H. P. Scheffler and C. Tadjeran. Finite difference method for two dimensional fractional dispersion equations. Journal of

Computational Physics, 2006, 211: 249-261.

[13] 夏源, 吴吉春. 分数阶对流—弥散方程的数值求解[J]. 南京大学学报(自然科学版), 2007, 43(4): 44-446.

[14] F. Liu, V. Ahn and I. Turner. Numerical solution of the space fractional Fokker-Planck equation. Journal of Computational & Applied

Mathematics, 2004, 166(1): 209-219.

[15] 周璐莹, 吴吉春, 夏源. 二维分数阶对流–弥散方程的数值解[J]. 高校地质学报, 2009, 15(4): 569-575.

[16] Y. Zhang. A finite difference method for fractional partial differential equation. Applied Mathematics and Computation, 2009, 215(2): 524-

529.

[17] M. M. Meerschaert, D. A. Benson, H. P. Scheffler and B. Baeumer. Stochastic solution of space-time fractional diffusion equations. Physical

Review, 2002, E65: 1103-1106.