劉桃花，侯木舟

(1.中南大學(xué)數(shù)學(xué)與統(tǒng)計(jì)學(xué)院，中國 長沙 410083; 2.邵陽學(xué)院理學(xué)系，中國 邵陽 422004)

一類分?jǐn)?shù)階反應(yīng)擴(kuò)散方程的差分方法

劉桃花1,2，侯木舟1

(1.中南大學(xué)數(shù)學(xué)與統(tǒng)計(jì)學(xué)院，中國 長沙 410083; 2.邵陽學(xué)院理學(xué)系，中國 邵陽 422004)

分?jǐn)?shù)階反應(yīng)擴(kuò)散方程可以用來模擬反常擴(kuò)散運(yùn)動(dòng)，它是由傳統(tǒng)的反應(yīng)擴(kuò)散方程演變而來的.本文對(duì)帶變系數(shù)的空間分?jǐn)?shù)階反應(yīng)擴(kuò)散方程的初邊值問題進(jìn)行了數(shù)值研究，采用了移位的Grunwald公式對(duì)空間分?jǐn)?shù)階導(dǎo)數(shù)進(jìn)行離散，在此基礎(chǔ)上建立了經(jīng)典的隱性Euler差分格式.然后討論了該格式的解的存在唯一性，分析了該方法相容性、穩(wěn)定性及收斂性，得到了O(τ+h)收斂階.最后用數(shù)值實(shí)驗(yàn)證明了該格式的有效性.

分?jǐn)?shù)階反應(yīng)擴(kuò)散方程； 隱性Euler差分格式；相容性；無條件穩(wěn)定；收斂性

分?jǐn)?shù)階(空間分?jǐn)?shù)階、時(shí)間分?jǐn)?shù)階和空間-時(shí)間分?jǐn)?shù)階)反應(yīng)擴(kuò)散方程是傳統(tǒng)的整數(shù)階反應(yīng)擴(kuò)散方程的推廣,即用分?jǐn)?shù)階導(dǎo)數(shù)代替相應(yīng)的整數(shù)階導(dǎo)數(shù)[1].由于分?jǐn)?shù)階微分方程比整數(shù)階方程更精確地描述了客觀世界,因此它在描述自然界擴(kuò)散現(xiàn)象中起到非常重要的作用,已廣泛地應(yīng)用于物理[2-4]、化學(xué)[5]、環(huán)境[6]等領(lǐng)域.空間分?jǐn)?shù)階反應(yīng)擴(kuò)散方程被用來模擬反常擴(kuò)散運(yùn)動(dòng)，若粒子以連續(xù)方式擴(kuò)散，就是經(jīng)典的布朗模型.

本文將考慮如下帶變系數(shù)的空間分?jǐn)?shù)階反應(yīng)擴(kuò)散方程：

(1)

初邊值條件為：

c(x,0)=s(x),L≤x≤R,

(2)

c(L,t)=0,c(R,t)=b(t),t≥0.

(3)

其中擴(kuò)散系數(shù)d(x)≥0;v(x)≥0,耗散系數(shù) p(x)≥0，均為[L,R]上的連續(xù)函數(shù).

(4)

目前,已有很多學(xué)者研究其解的存在性和唯一性，取得了很好的結(jié)果[9-10].關(guān)于同類別的問題的數(shù)值方法有很多學(xué)者進(jìn)行了研究[11-14]，Xie等[15]在研究空氣中PM2.5的擴(kuò)散情況時(shí)，給出了此類分?jǐn)?shù)階擴(kuò)散方程的Crank-Nicholson格式.本文在采用移位Grunwald公式對(duì)空間分?jǐn)?shù)階導(dǎo)數(shù)進(jìn)行離散，構(gòu)造了經(jīng)典的隱性Eluer差分格式來求出此方程的數(shù)值解.

文章其它部分安排如下：第2部分，對(duì)帶變系數(shù)的空間分?jǐn)?shù)階的分?jǐn)?shù)階擴(kuò)散方程構(gòu)建了隱性Eluer差分方法并證明其解的存在唯一性，以及穩(wěn)定性、收斂性.第3部分，給出一個(gè)數(shù)值實(shí)驗(yàn)，來檢驗(yàn)此方法的有效性.

1 差分格式的建立

(5)

其中g(shù)k為Grunwald權(quán)系數(shù)

定義一階空間及時(shí)間導(dǎo)數(shù)：

其中，τ=Δt，因此，可得到

(6)

存在M為正常數(shù)，使得：

|Ri,n+1|≤M(τ+h).

由上式，對(duì)問題(1)～(3)可以建立隱性的Eluer格式如下：

(7)

(8)

(9)

由上可知式(7)～(9)與帶初邊值條件方程(1)～(3)相容.

方程(7)通過整理得如下形式

(10)

(Ei=viΔt/h,Bi=diΔt/hα).

可進(jìn)一步將分?jǐn)?shù)階方程改寫為矩陣形式

ACn+1=Cn+ΔtFn+1,

A=[Ai,j]是一個(gè)(K-1)×(K-1)系數(shù)矩陣，當(dāng)i=1,…,K-1和j=1,…,K-1時(shí)

(11)

2 差分格式的差分解的存在唯一性以及穩(wěn)定性及收斂性分析

引理1[10]當(dāng)1<α≤2時(shí)，Grunwald權(quán)系數(shù)gi(i=0,1,2,…)滿足：

定理2 差分格式(6)～(8)的解存在唯一，而且差分格式是無條件穩(wěn)定的.

證 假設(shè)λ為矩陣A的特征值，X為其對(duì)應(yīng)的特征向量，即AX=λX，設(shè)

(12)

將式(11)代入式(12)得

命題1 上面分析時(shí)，已知該方法是相容的，且得到該方法的截?cái)嗾`差為O(τ+h).根據(jù)Lax等價(jià)定理[16]，差分格式(6)～(8)的解cn以‖·‖∞收斂到初邊 值問題的解，且收斂階為O(τ+h).

3 數(shù)值實(shí)驗(yàn)

考慮如下分?jǐn)?shù)階反應(yīng)擴(kuò)散方程(α=1.5)

圖1 數(shù)值結(jié)果Fig.1 Numerical results

此方程的精確解為c(x,t)=x2e-t

圖1數(shù)值結(jié)果表明，在T=1時(shí)刻，取相同的空間、時(shí)間步長(τ=h=1/10)時(shí)，方程數(shù)值解與精確解相吻合，即所給出差分格式是有效的.

表1為T=1時(shí)刻，取相同的空間、時(shí)間步長時(shí)，數(shù)值解的最大誤差及誤差階.在表1中，當(dāng)空間步長和時(shí)間步長都減半時(shí)，誤差接近原來的一半，這表明差分格式的收斂階為O(τ+h).

表1 最大誤差及誤差階

Tab.1 Max error and error rate

τ=h最大誤差誤差階1/108.03228×10-31/204.21272×10-31.906671/402.17059×10-31.940821/801.22580×10-31.77075

[1] METZLER R,KLAFTER J.The random walk’s guide to anomalous diffusion:A fractional dynamics approach[J].Phys Rep,2000,339(1):1-77.

[2] SAICHEV A I,ZASLAVSKY G M.Fractional kinetic equations:solutions and applications[J].Chaos,1997,7(4):753-764.

[3] ROSSILIHIN Y,SHITIKOVA M.Applications of fractional calculus to dynamic problems of linear and nonlinear hereditary mechanics of solids[J].Appl Mech Rev,1997,50(2):15-67.

[4] METZLER R,KLAFTER J.Boundary value problems for fractional diffusion equations[J].Phys A:Stat Mech Appl,2000,278(1-2):107-125.

[5] YUSTE S B,LINDENBERG K.Subdiffusion-limited A+A reactions[J].Phys Rev Lett,2001,87(11):118-301.

[6] BENSON D A,WHEATCRAFT S W,MEERSCHAERT M M.Application of a fractional advection-dispersion equation[J].Water Resour Res,2000,36(6):1403-1412.

[7] TADJERAN C,MEERSCHAERT M M,SCHEFFLER H P.A second-order accurate numerical approximation for the fractional diffusion equation[J].J Comput Phys,2006,213(1):205-213.

[8] TUAN V K,GORENFLO R.Extrapolation to the limit for numerical fractional differentiation[J].Z Agnew Math Mech,1995,75:646-648.

[9] KEMPPAINEN J.Existence and uniqueness of the solution for a time-fractional diffusion equation with Robin boundary condition [J].Abstr Appl Anal Article,2011,(11):ID321903.

[10] PODLUBNY I.Fractional Differential Equations,Volume 198.An introduction to fractional derivatives,fractional differential equations,to methods of their solution and some of their applications[M].San Diego:Academic Press,1998.

[11] MEERSCHAERT M M.TADJERAN C.Finite difference approximations for fractional advection-dispersion flow equations[J].J Comput Appl Math,2004,172(1):65-77.

[12] LIU F,ANH V,TURNER I.Numerical solution of space fractional Fokker-Planck equation[J].Comp Appl Math,2004,(166):209-219.

[13] 尹修草，周 均，胡 兵.分?jǐn)?shù)階對(duì)流-彌散方程的有限差分方法[J].四川大學(xué)學(xué)報(bào)(自然科學(xué)版),2013,50(3):409-413.

[14] FIX G J,ROOP J P.Least squares finite-element solution of a fractional order two-point boundary value problem[J].Comput Math Appl,2004,48(7-8):1017-1033.

[15] XIE C P,LI L,HUANG Z Z,etal.Fractional difference approximations for fractional reaction-diffusion equations and the application in PM2.5[C].International Symposium on Energy Science and Chemical Engineering,2015-12,Guangzhou.

[16] CHAVES A.Fractional diffusion equation to describe Levy flights[J].Phys Lett A,1998,239(1-2):13-16.

(編輯 HWJ)

Finite Difference Approximations for Fractional Reaction-Diffusion Equations

LIUTao-hua1,2*,HOUMu-zhou1

(1.School of Mathematics and Statistics,Central South University,Changsha 410083,China；2.Department of Science,Shaoyang University,Shaoyang 422004,China)

Fractional reaction-diffusion equations are generalizations of classical reaction-diffusion equations,which are used in simulating the anomalous diffusion motion.In this paper,we examine a practical numerical method,which is called Euler method to solve a class of initial-boundary value a fractional reaction-diffusion equation with variable coefficients.Then we discuss the existence and uniqueness of solutions for the format.The stability,consistency and convergence of the method are established to get the convergence order ofO(τ+h).Finally,we use a numerical experiment to prove the effectiveness of the proposed format.

fractional reaction-diffusion equations; implicit Euler method; consistency; unconditional stability; convergence

10.7612/j.issn.1000-2537.2017.01.014

2016-09-04

國家自然科學(xué)基金資助項(xiàng)目(61672356,61375063,61271355,11271378,11301549)；中南大學(xué)2015年實(shí)驗(yàn)室研究項(xiàng)目

* 通訊作者,E-mail：liutaohua2005@163.com

O241.82

A

1000-2537(2017)01-0091-04