孙春龙，李功胜，贾现正，杜殿虎

(山东理工大学 理学院， 山东 淄博255049)



含三个时间分数阶导数的反常扩散方程求解与微分阶数反演

孙春龙，李功胜，贾现正，杜殿虎

(山东理工大学 理学院， 山东 淄博255049)

摘要：对于一类带有三个时间分数阶的一维反常扩散问题,基于Caputo意义下时间分数阶导数的离散,给出了一个有限差分求解格式,并利用分离变量法及Laplace变换得到该问题的解析解.进一步应用同伦正则化算法,根据内点处的浓度观测数据对确定微分阶数的反问题进行数值反演,并讨论时间-空间步长及数据扰动等因素对反演算法的影响.

关键词：时间分数阶导数； 含三个时间分数阶的扩散； 反问题； 同伦正则化算法； 数值反演

在过去的几十年,分数阶偏微分方程在反常扩散现象模拟等领域开始发挥重要作用,其对复杂系统的描述具有建模简单、参数物理意义清楚、描述准确等优势.陈文案[1]较详细地描述了分数阶微积分理论在复杂力学行为建模和数值模拟方面的研究成果，强调分数阶微积分建模的物理和力学背景和概念.分数阶反常扩散方程是由经典扩散方程推广而得到的,按照其含有的时间-空间导数的形式可分为三类：分别是时间分数阶扩散、空间分数阶扩散和空间-时间分数阶扩散.当通常的整数阶扩散方程中关于时间变量的1阶导数换为α(0<α≤1)时,得到时间分数阶扩散方程；当整数阶扩散方程中关于空间变量的2阶导数换为α(0<α≤2)时,得到空间分数阶扩散方程；而当整数阶扩散方程中关于时间与空间变量的整数阶导数都换为相应的分数阶导数时,就得到空间-时间分数阶扩散方程.

对于时间分数阶扩散方程正问题的研究,刘发旺等人做了很多研究.于强与刘发旺[2]给出了时间分数阶反应-扩散方程的隐式差分近似；沈淑君与刘发旺[3]推导出了时间分数阶扩散方程初边值问题的数值解,并提出了方程的显式守恒差分格式,证明了此格式的稳定性和收敛性,然后将得到的结果推广到时间分数阶对流-扩散方程.刘法旺和Meerschaert[4]等人,给出了含有多个时间分数阶的扩散-波动方程的数值计算方法,并讨论了隐式差分求解格式.Luchko[5]与Liu等[6]分别研究了多个时间分数阶的对流-扩散方程,得到了广义Mittag-Leffler函数形式的解析解.

对于分数阶扩散中的反问题研究,已引起了大家的关注.程晋等[7]研究了一维时间分数阶扩散中由单个边界点处的连续测量值同时确定微分阶数与空间依赖扩散系数的反问题,应用特征展开方法证明了该反问题解的唯一性.刘继军和Yamamoto[8]应用准可逆性正则化(quasi-reversibilityregularization)方法探讨了一个时间分数阶扩散中由终值数据恢复初始状态的倒向反问题.郑光辉与魏婷等[9-10]应用Fourier谱正则化方法研究了带状区域一维时间分数阶对流扩散的Cauchy问题与空间Riesz分数阶扩散的倒向反问题.Sakamoto与Yamamoto[11]利用特征分解方法探讨了时间分数阶扩散-波动方程倒向反问题的适定性.徐翔、程晋与Yamamoto等[12-13]研究了具有1/2阶分数阶导数的时间分数阶扩散方程的Carleman估计方法,这是构建分数阶扩散相关反问题条件稳定性的新思路.熊向团等[14]构建了时间分数阶二维热传导方程的反向热传导问题的一种条件稳定性,李功胜等[15]利用最佳摄动量算法对于一维时间分数阶扩散方程的空间依赖的扩散系数进行了有效的数值反演,魏婷等[16]基于边界元离散提出一种正则化方法,研究了一维时间分数阶扩散方程中时间依赖源项的重建问题.

本文将探讨含有三个时间分数阶导数的一维反常扩散方程求解及确定微分阶数的反演问题.首先利用有限差分法给出正问题数值求解的差分格式,同时利用分离变量法及拉普拉斯变换得到正问题的解析解表达式,并且给出正问题的数值算例.进一步,应用同伦正则化算法对反常扩散方程中的三个时间分数阶微分阶数进行数值反演,给出数值反演算例,并讨论时间空间步长及数据扰动水平等因素对反演算法的影响.

1正问题及其数值求解

考虑区域Ω=(0,l)×(0,T)上含三个时间分数阶导数的一维反常扩散方程

(1)

其中(x,t)∈Ω;D>0是扩散系数,α∈(0,1)和αj(j=1,2)∈(0,1)是时间分数阶导数的微分阶数，且0<α2<α1<α<1.对于方程(1),给定非零初值条件和第一类零边值条件:

u(x,0)=f(x),u(0,t)=u(l,t)=0

(2)

这样,由方程(1)及初边值条件(2)就构成含三个时间分数阶导数的一维反常扩散方程的正问题.下面给出这个正问题求解的一个有限差分格式，并运用分离变量法推导其解析解的表达式.

1.1正问题求解的差分格式

(3)

(4)

对于方程中的整数阶导数项,按照通常的差分离散方法,即有

(5)

(6)

1.2正问题的解析解

本小节利用分离变量法及Laplace变换推导正问题的解析解.

设u(x,t)=X(x)T(t),代入方程可得

(7)

(8)

(9)

(10)

(11)

(12)

应用Laplace变换,得到

(13)

代入(12)式,可得

记c=Dλn为常数,可得

(14)

(15)

于是求得解析解为

(16)

1.3数值试验

利用上一节的有限差分法进行数值计算,取l=π,T=1,且离散点数M=100,N=100，微分阶数α=0.8,α1=0.5,α2=0.3,扩散系数D=1,系数r1=0.5,r2=0.5.另外,取初始函数u(x,0)=f(x)=sinx.此时解析解可表示为

(17)

对上述解析解,级数取前50项截断计算,记数值解u*(x,t),相对误差表示为Err=‖u(x,t)-u*(x,t)‖2/‖(u,t)‖2.数值结果分别列于表1～表3和图1.

表2 时间步长对解误差的影响(t=0.5)

表3　不同微分阶数对正问题求解的影响(t=0.5)

图1 　t=0.5时的数值解与解析解

从表1的计算结果可以看出,空间步长对误差的影响不大;由表2的计算结果可以看出,随着时间步长的减小,误差逐渐减小;由表1、表2及图1的结果可以看到，正问题的数值解和解析解吻合的较好.从表3看出,时间微分阶数对正问题求解具有一定的影响,微分阶数越接近,解误差有所减小.

2确定微分阶数的反问题与反演算法

2.1反问题的提出

假设方程(1)中的时间微分阶数未知,那么为了确定各个微分阶数的值,需要补充关于正问题解的部分条件信息,并联合正问题(1)～(2)形成一个确定微分阶数的反问题.本文给定内点x=x0,时刻tj(j=1,2,…，N)的观测值为附加数据,记为u(x0,tj)=θ(tj),则可定义附加数据向量V：

(18)

这样,由附加数据(18)联合正问题(1)～(2)构成了一个确定微分阶数(α,α1,α2)的数值反演问题.下面给出同伦正则化算法,并对上述确定微分阶数的反问题进行数值反演模拟.

2.2同伦正则化算法

记a=(α,α1,α2),利用上一节的差分方法求解正问题可得其解,并在x=x0处赋值,记之为u(x0,t;a),这可以看作对应于输入数据a=(α,α1,α2)的一个计算输出.

min‖U(a)-V‖2

(19)

为了克服病态性,应用同伦思想可将上式转化为如下极小值问题：

min{(1-λ)‖U(a-V‖2+λ‖a‖2}

(20)

(21)

根据同伦正则化思想,上述极小问题(20)的求解又转化为对于给定的an,通过求解最佳摄动量δan进而确定an+1的一种迭代算法：an+1=an+δan,n=0,1,2，…

(22)

(23)

将u(x0,tj;an+δan)在an处作泰勒展开,并略去高阶项,可以得到

则目标函数F(δan)可近似表为

(24)

(25)

式中ω为数值微分步长.不难验证泛函极小问题等价于求解规范方程：

((1-λ)GTG+λI)δan=(1-λ)GT(V-aU)，

(26)

则得到每一步迭代的最优摄动量δan的计算公式：

δan=((1-λ)GTG+λI)-1(1-λ)GT(V-U)

(27)

对于给定的精度eps,判断‖δan‖≤eps是否成立,若满足则δan即为所求;否则由式(22)得到an+1,再通过规范方程继续求解.

3数值反演

本节应用同伦正则化算法对确定微分阶数的反问题(1)～(2)及(18)进行数值反演.取初始函数f(x)=sinx,以下若无特殊说明,正问题计算中均取l=π,扩散系数D=1,终值时刻T=1,且M=100,N=100.附加数据取为在x0=0.5处的观测数据，数值微分步长ω=0.01,停止准则eps=1e-6,预估迭代次数N0=5,β=0.8,设式(1)中r1=r2=0.5,取α=0.9,α1=0.7,α2=0.5，则微分阶数真值a=(0.9，0.7，0.5)，ainv表示反演解,Err=‖ainv-a‖2/‖a‖2表示反演解与真解的误差.

3.1附加数据不带扰动的情形

(1)时间步长对反演结果的影响

令初始迭代值a0=(0.3,0.2,0.1),反演计算结果列于表4,n表示迭代次数.

表4　时间步长对反演结果的影响

表5　空间步长对反演结果的影响

由表4可以看出,时间步长对反演结果影响不大,只是随时间步长的变小迭代的次数稍微减少.

(2)空间步长对反演结果的影响

取初始迭代值a0=(0.3,0.2,0.1),反演计算结果列于表5,n表示迭代次数.

由表5可以看出,空间步长对反演结果的影响很小.

(3)微分阶数对反演结果的影响

令初始迭代值a0=(0.3,0.2,0.1),考察时间分数阶α,α1,α2的不同取值对反演结果的影响,计算结果见表6,n表示迭代次数.

表6　微分阶数对反演结果的影响

通过表6可以看出微分阶数对反演结果有一定的影响,随着微分阶数的接近误差有所增加,并且迭代步数也略微增加.

3.2附加数据有扰动的情形

(1)扰动水平对反演结果的影响仍取微分阶数真值a=(0.9,0.7,0.5),数值微分步长τ=0.001,分别取扰动水平δ=1%,0.1%,0.05%,0.01%,20次反演计算的平均值列于表7.

通过表7可以看出随着数据扰动水平的减小,反演结果精度越来越高.

(2)计算次数对反演结果的影响

表7　扰动水平对反演结果的影响

表8 给定扰动水平下计算次数对反演结果的影响

设扰动水平δ=1%,考察计算次数对反演结果的影响,结果见表8,其中p表示计算次数.

由表8可以看出,对于给定的扰动数据,计算次数对反演结果的影响不大.

4结束语

本文探讨了含有三个时间分数阶导数的一维反常扩散方程求解及确定微分阶数的反问题.正问题的数值模拟结果表明,数值解与解析解吻合较好,所建立的有限差分格式是有效的.对于确定三个微分阶数的反问题,虽然理论上还没有解决反演的稳定性问题,但应用同伦正则化算法可以实现对多个微分阶数的数值反演.通过反演计算结果发现,当数据扰动水平较大时,解误差变得较大.这说明含多个时间分数阶导数的反常扩散方程的微分阶数反问题具有较高的病态性,其他更有效的反演算法是需要进一步研究的课题.

参考文献:

[1] 陈文,孙洪广，力学与工程问题的分数阶导数建模[M].北京：科学出版社，2010.

[2] 于强,刘发旺.时间分数阶反应-扩散方程的隐式差分近似[J],厦门大学学报:自然科学版, 2006,45(3): 315-319.

[3] 沈淑君.分数阶对流-扩散方程的基本解和数值方法[D].厦门：厦门大学, 2008.

[4] Liu F, Meerschaert M M, McGough R J,etal. Numerical methods for solving the multi-term time-fractional wave-diffusion equations[J]. Fractional Calculus and Applied Analysis, 2013, 16: 9-25.

[5] Yury L.Initial-boundary-value problems for the generalized multi-term time-fractional diffusion equation[J].Journal of Mathematical Analysis and Applications, 2011, 374:538-548.

[6] Jiang H, Liu F, Turner I,etal. Burrage,Analytical solutions for the multi-term time-fractional diffusion-wave/diffusion equations in a finite domain [J]. Computers and Mathematics with Applications,2012, 64:3 377-3 388.

[7] Cheng J, Nakagawa J, Yamamoto M,etal. Uniqueness in an inverse problem for a one-dimensional fractional diffusion equation [J].Inverse Problems, 2009,25:115 002-115 017.

[8] Liu J J,Yamamoto M,A backward problem for the time-fractional diffusion equation [J], Applicable Analysis, 2010,89:1 769-1 788.

[9] Zheng G H,Wei T.Spetral regularization method for a Cauchy problem of the time fractional advection-dispersion equation [J]. Journal of Computational and Applied Mathematics, 2010, 233:2 631-2 640.

[10] Zheng G H, Wei T.Two regularization methods for solving a Riesz-Feller space-fractional backward diffusion equation [J]. Inverse Problems, 2010, 26:115 017-115 038.

[11] Sakamoto K,Yamamoto M. Initial value/boundary value problems for fractional diffusion-wave equations and applications to some inverse problems [J]. Journal of Mathematical Analysis and Applications, 2011,382:426-447.

[12] Xu X, Cheng J,Yamamoto M.Carleman estimate for a fractional diffusion equation with half order and application [J].Applicable Analysis, 2011,90:1 355-1 371.

[13] Yamamoto M,Zhang Y.Conditional stability in determining a zeroth-order coefficient in a half-order fractional diffusion equation by a Carleman estimate [J].Inverse Problems, 2012,28(10): 105 010-105 019.

[14] Xiong X T, Zhou Q,Hon Y C.An inverse problem for fractional diffusion equation in 2-dimensional case: Stability analysis and regularization [J].Journal of Mathematical Analysis and Applications, 2012,393,185-199.

[15] Li G S, Gu W J, Jia X Z. Numerical inversions for space-dependent diffusion coefficient in the time fractional diffusion equation [J].Journal of Inverse and Ill-Posed Problems, 2012, 20(3): 339-366.

[16] Wei T,Zhang Z Q.Reconstruction of a time-dependent source term in a time-fractional diffusion equation [J]. Engineering Analysis with Boundary Elements, 2013,37:23-31.

[17] Igor Podlubny, Fractional Differential Equations [M].San Diego:Academic press,1999.

(编辑：刘宝江)

Thesolutiontothreetime-fractionalanomalousdiffusion

equationandnumericalinversionofthefractionalorders

SUNChun-long,LIGong-sheng,JIAXian-zheng,DUDian-hu

(SchoolofScience,ShandongUniversityofTechnology,Zibo255049,China)

Abstract:Afinitedifferenceschemeisintroducedtosolvethe1-Dthree-termtime-fractionalanomalousdiffusionequationbasedonCaputo'sdiscretizationtothetimefractionalderivatives.UsingthemethodofseparationofvariablesandLaplacetransform,theanalyticalsolutionoftheforwardproblemisobtained.Furthermore,thehomotopyregularizationalgorithmisappliedtodeterminethethreetime-fractionalordersbytheadditionalmeasurementsataninteriorpointinthedomain.Numericalinversionsareperformedtodemonstrateeffectivenessoftheproposedalgorithm,andinfluencesofthetime-spacemeshgridsandthedatanoisesontheinversionalgorithmarediscussed.

Keywords:timefractionalderivative;threetime-fractionaldiffusion;inverseproblem;homotopyregularizationalgorithm;numericalinversion

中图分类号：O175

文献标志码：A

文章编号：1672-6197(2015)03-0001-07

作者简介：孙春龙，男，sunchunlong527@163.com; 通信作者： 李功胜，男，ligs@sdut.edu.cn

基金项目：国家自然科学基金资助项目( 11071148, 11371231); 山东省自然科学基金资助项目(ZR2011AQ014)

收稿日期：2014-09-20