引言
现如今，随着社会的发展以及计算工具的更新迭代，分数阶微分方程作为一门新兴的建模工具正在运用于许多中日常生活领域，材料学，金融经济学，物理力学，生物遗传学等[1][7]。因为有了应用和需要的背景，以及运用分数阶微分方程建模分析问题的兴起，在一定程度上促进了人们对于分数阶微积分方程深入研究。由于分数阶微积分是数学建模中的重要工具，常见模型有分数阶扩散模型，麦克斯韦模型和本构模型等[10],随着模型的建立，也引发了人们对分数阶微分方程的求解问题的探索。
我们常常用拉普拉斯变换法和傅里叶变换法去求解分数阶微分方程，但是算出的解析解在日常领域中并没有什么实际作用。在许多情况下，我们也很难找到方程的解析解，所以在日常生活中我们对数值解的使用更为广泛。因为我们通常用有限差分，有线元和谱方法去求整数阶微积分方程的数值解，因此我们也可以用同样的方法求解出时间分数阶微分方程的数值解。求解整数阶微分方程的数值解关键在于如何处理方程的初边值问题，不同的是求解分数阶常微分方程的关键是方程的离散格式[12]。
截止到目前为止，对分数阶常微分方程的研究已经取得了许多成果，从中我们可以发现在现实生活中运用的最多的就是分数阶扩散方程[23]。而在研究分数阶微分方程的模型中，运用最多的方法就是有限差分法。相比较而言我们所要介绍的时间分数阶微分方程的数值解法，在研究方面可借鉴的方法就相对而言比较少。时间分数阶微分方程在本质上是经典扩散方程的一种，将经典扩散或者波动方程中的时间导数替换为阶分数阶导数即为时间分数阶微分方程。当时，时间分数阶扩散方程同时具有扩散和波动的性质，当时，该方程则为时间分数阶次扩散方程。人们通常用有限差分或者谱方法来求解这类方程[12]。
这篇文章主要介绍的是时间分数阶波动方程的数值解法[5][13]，为了更好的描述这个算法，我将整篇论文分成了三个部分，第一部分主要就是描述微积分方程相关的定义定理，Caputo导数,Riemann-Liouville积分和微积分算子的定义定理。第二部分就是时间分数阶方程数值解法的具体实现，主要分为两个步骤，首先用微积分算子的相关定理得到等价方程，在运用数值逼近中心离散的方法离散方程。第三部分为数值实验，主要就是目的就是求解数值算例，来验证此数值算法的有效性，并进行误差分析，最后得出结论。
1.相关预备知识
1.1分数阶微积分的定义与性质
定义1.1.1设对于实函数若存在实数使得，那么称函数属于空间；当且仅当时，属于空间，。
定义1.1.2若，则阶Riemann-Liouville分数阶积分定义为：
，.
特别地，当时，简写为.
引理1：若，则有
（1），（2）
（3），（4）
（5），（6）
定义1.1.3函数，定义阶Caputo分数阶导数为：
特别地，当时，简写为.
定义1.1.4 Riemann-Liouville分数阶导数定义
若函数，定义其阶Riemann-Liouville分数阶导数为：
.
引理2若，函数，并设存在，那么
.
1.2分数阶微积分算子的定义及引理
定义1.2.1 Riemann-Liouville积分算子的定义
设，算子定义在上，，
其中，称为Riemann-Liouville分数阶积分算子。
定义1.2.2 Riemann-Liouville微分算子的定义
设对算子定义为：
称为阶的分数阶微分算子。当时，是单位算子。
定义1.2.3 Caputo分数阶积分算子
设对算子定义为：
称为阶的Caputo分数阶微分算子。当时，是单位算子。
引理3 Riemann-Liouville分数阶积分的性质
算子集合关于复合运算构成可交换的半群。算子构成半群的单位元。由此可得下列等式成立，
，，
引理4 Riemann-Liouville分数阶导数性质
（1）设且，则.
（2）对于任给的设和为定义在上的两个函数且满足和几乎处处存在，则几乎处处存在，且.
（3）Riemann-Liouville算子的Leibniz’公式
设，对给定的和上解析，得到
.
引理5 Caputo分数阶导数性质
（1）对于任给的，设和为定义在上的两个函数且满足和几乎处处存在，且.
（2）Caputo算子的Leibniz’公式
设，对给定的，和在上解析，则
引理6 Riemann-Liouville分数阶积分与导数的关系
（1）设，且，则
（2）设，则对任一给定的，几乎处处成立。
（3）设，如果存在使得，则几乎处处成立。
引理7 Caputo分数阶导数与积分的关系
（1）若连续且，则.
（2）设，则.
引理8 Riemann-Liouville分数阶导数与Caputo分数阶导数之间关系与性质
（1）设，并假设存在且在点有阶导数。则
几乎处处存在，代表在点的阶展开多项式。
（2）设.假设和都存在。则当且仅当
1.3微积分数值计算的相关知识
1.3.1数值微分的条件和方法
当出现函数用离散点列给出要求导数值或函数较为复杂的情况时，我们就可以用数值的方法来求导原函数。取极限的近似值是我们最常用也最简单的方法，就是我们常说的差商。差商又可分为向前差商，向后差商和中心差商，插值型数值微分也可近似为原函数的导数，最后用Taylor展开分析[5]。
1.3.2数值积分的条件和方法
当函数由离散数据组成，求不出或过于复杂时，这种情况下我们通常使用数值积分求解这类问题。数值积分的定义：我们把离散点上函数值的线性组合称为数值积分。
其中为积分系数与无关，与积分区间和积分点有关。
我们常用的求解数值积分的方法有矩形法，梯形法和Simpson公式法。前两种，矩形法和梯形法是以面积法的思想来近似求解定积分的方法，矩形法是将函数划分的面积看作若干个矩形面积的总和从而得出积分的近似值，梯形法则是将函数划分的面积看作若干个梯形面积的总和从而得出积分的近似值。由此我们可以得到：
假设被积函数为，将[a,b]区间n等分，区间的长度为。
用矩形法所求积分结果为：
用梯形法所求积分结果为：
与矩形法和梯形法直线线段拟合函数曲线不同的是，Simpson公式法是以二次曲线逼近的方式取代矩形积分公式，来求定积分的数值近似解的[6][11]。将整个区间等分，子区间长度为.
若子区间是等距的，积分结果为：
2时间分数阶微分方程的算法设计
2.1时间分数阶微分方程的等价变换
以下是我们所要考虑的时间分数阶微分方程：
（2.1）
初始条件为
（2.2）
边界条件为
（2.3）
其中为扩散波动的常系数，表示空间方向变量，为源汇项，为关于时间的Caputo分数阶导数，定义如下
（2.4）
因为Caputo分数阶导数算子具有分裂性[21]，Caputo分数阶导数算子可以分裂成算子和算子表示的，即为
（2.5）
由等价于下面的偏微分方程
（2.6）
其中是阶关于时间的Riemann-Liouville分数阶积分算子
定义为：
查询文献资料可知[12]，构造有限差分格式有两个优点：一是降低函数在时间方向的光滑性要求，二是对积分的离散所得到的数值格式比离散高阶导数的格式更加稳定。
2.2时间分数阶微分方程的一种离散差分格式
有限差分离散的过程，设是任一给定的整数，令为步长，定义为区间的网格剖分，其中，用同样的方式定义为区间的网格剖分，其中，，且满足存在满足。
然后用向后差分来离散方程中的时间导数项，以及中心差分来离散项，由此可得下列等式
考虑含有奇异核的积分部分
记
则
（2.7）
对于和，
（2.8）
用数值积分法，我们可以得到公式（2.8）的数值逼近
由于积分（2.7）在趋于时收敛的，所以我们只需要考虑时的情况。
数值实验表明如果我们在在区间和上分别做出最优近似解，误差会更小。
（2.9）
（2.10）
当时，
（2.11）
因为（2.10）式中的第一部分积分，是奇异积分，我们用插值近似，
.（2.12）
在考虑公式（2.11）式中第二部分，
令
则
（2.13）
用梯形公式计算公式（2.13）中的积分部分得出
（2.14）
结合以上公式，可得一个求解方程的有限差分格式
（2.15）
其中
,
令
所以（2.15）可改写为
当时，
再利用近似公式（2.9）（2.10）来计算公式（15）里的后两部分积分
（2.16）
令
有
（2.17）
（2.18）
用梯形公式近似计算公式（2.17）（2.18）中积分的部分
,（2.19）
（2.20）
综合所述，可得方程公式（2.7）的有限差分格式
，.
3数值实验
3.1数值算例的条件
（3.1）
，
.
这个微分方程的数值解为：
其中为单函数的Mittag-Leffler函数，由以下的级数来定义
我们可以轻易得出方程（3.1）等价与下面的微积分方程
图1是,,时间分数阶方程的数值解和精确解的比较
图2是在,,时间分数阶方程的数值解和精确解的比较
通过比较图1和图2，我们发现求出的方程数值解的和方程精确解的曲线几乎重合，误差并不大，由此可以说明设计算法的准确性及有效性。
3.2误差分析结果
图3是在,,时间分数阶方程的数值解和精确解的比较
表1是在时，不同值下的误差及阶数
当时，在不同的值下的产生的误差
误差阶数误差阶数误差阶数
1/40 4.495e-4 5.213e-4 0.00469
1/80 2.320e-4 0.9541 2.811e-4 0.8910 0.00248 0.9192
1/160 1.169e-4 0.9888 1.450e-4 0.9550 0.00128 0.9541
表2是在时，不同值下的误差及阶数
当时，在不同的值下的产生的误差
误差阶数误差阶数误差阶数
1/40 6.113e-5 8.007e-8 1.300e-4
1/80 2.906e-5 2.0264 3.760e-5 2.0258 1.904e-4 2.0154
1/160 2.104e-5 2.0035 2.699e-5 1.9997 2.056e-4 2.0282
图3是将时间分数阶微分的数值解和精确解的局部放大图，从中我们可以看处我们设计的算法计算出的数值和精确解还有一定误差。
表1和表2是，时间分数阶微分方程的误差分析表，将计算得出的阶数分别与时间和空间阶数对比发现误差并不是很大，在有效范围内。因此进一步验证我们用算法的准确性和有效性。
查阅相关资料[20]，并通过对差分的格式的分析[21]，总结出算法的误差主要取决于公式（2.17）和公式（2.18）的近似和有限差分法梯形公式的整理。
结论
这篇论文主要描述了一些与微积分方程的相关基本定义与性质，介绍了一种计算时间分数阶微分方程数值解的差分格式。先用向后差分的方法将时间分数阶微分方程中的时间导数项离散化，再用中心差分的方法去离散空间项。因为方程离散后包含奇异核积分对计算量影响较大，为此我们利用了分级数逼近的方法构建了计算积分的递推公式，减少了计算量。最后，通过比较方程的精确解与数值解，证明了设计算法准确性与有效性并进行了误差分析。
参考文献
[1]刘发旺，庄平辉，刘青霞.分数阶偏微分方程数值方法及其应用[M].信息与计算科学丛书,2015.
[2]孙志忠,袁慰平,闻震初.数值分析（第二版）[M].东南大学出版社,2002.
[3]王自强，曹俊英.分数阶微分积分方程的数值解法及其误差分析[M].成都：西南交通大学出版社，2015.
[4]I.Podlubny.Fractional Differential Equations[M],New York:Academic Press,1999.
[5]Jianfei Huang,Yifa Tang,Luis Vazquez,Convergence analysis of a block-by-block method for fractional differential equations[J],Numerical Mathematics:Theory,Methods and Applications,5(2):229-241(2012).
[6]Sadia Arshad,Jianfei Huang,Abdul Q.M.Khaliq,Yifa Tang.Trapezoidal scheme for time-space fractional diffusion equation with Riesz derivative[J].Journal of Computational Physics,350:1-15(2017).
[7]蔡新,刘发旺.解空间Riesz分数阶扩散方程的一种数值方法[J].高等学校计算数学学报,27:243-244(2005).
[8]孙志忠,高广花.分数阶微分方程的有限差分方法[M].信息与计算科学丛书,2015.
[9]朱琳.求解空间分数阶扩散方程和对流扩散方程的有限差分格式研究[D].宁夏大学,2016.
[10]陈文.反常扩散的分数阶微分方程和统计模型[M].北京：科学出版社，2017.
[11]邓伟华.分数阶微分方程的理论分析与数值计算[D].上海大学,2007.
[12]孙志忠，高广花.分数阶微分方程的有限差分方法[M].科学出版社，2018.
[13]吴迎，黄健飞，赵维加.分数阶常微分方程的数值解法[J].青岛大学学报，28:159-160(2013).
[14]丁小丽.分数阶微分方程的解析研究方法[M].西北工业大学出版社，2019.
[15]黄健华，辛杰，沈天龙.分数阶偏微分方程的动力学[M].科学出版社，2019.
[16]刘发旺,庄平辉,刘青霞.分数阶偏微分方程数值方法及其应用[M].中国科学出版社，2013.
[17]龚晟.简明微积分[M].高等教育出版社，2017.
[18]覃平阳，张晓丹分数阶波动方程的一种数值解法[J].高等学校计算数学报，17:139-140(2016)
[19]白占兵.分数阶微分方程边值问题理论及应用[M].中国科学出版社，2013.
[20]张淑情.有限区间上的分数阶扩散-波方程的解[J].西北师范大学学报，5:171-174(2005)
[21]张明，张晓丹.解波动方程的精细积分法及其数值稳定性分析[J].石油大学学报，2004
[24]Liu,F.,Anh,V.and Turner,I.and P.Zhuang,Numerical simulation for solute transport in fractal porous media[J].ANZIAM J.,2004,45(E):461-473.
[25]于强，刘发旺.时间分数阶反应扩散方程的隐式差分近似[N].厦门大学学报，2006
[26]LAKSHMIKANTHAM V,VATSALA A S.Basic Theory of Fractional Differential Equations[J].Nonlinear Analysis,2008,69(8):2677-2682
[27]郭柏灵，蒲学科，黄凤辉.分数阶偏微分方程及其数值解[M].科学出版社，2013.
[28]Podlubny I.Fractional differential equations[M].San Diago,Academic Press,1999
致谢
这篇论文是在导师副教授的指导下完成的，首先在这里要对我的论文导师表达最诚挚的感谢。从一开始，导师对我们的论文进度和任务就有着明确的规划，在毕业论文的写作中遇到的种种困难，他也会悉心指导我们妥善解决。在论文写作过程中我也在相关网络上拜读了一些黄老师的论文期刊，可以说是受益匪浅。他在对分数阶微分方程的钻研与严谨精神也更值得我们学习。
在这大学四年本科的学习中，接触到很多优秀的老师和同学。这篇论文的顺利完成离不开学院勤勤恳恳的授课老师，也离不开同学们的鼓励与帮助，在这里我要感谢那些教育和帮助我的老师和同学们，他们教会了我很多东西。
最后要感谢我的家人，他们的鼓励和鞭策是我在学习和生活中最大的动力，感谢他们对我的支持和鼓励，让我成长。感谢在百忙之中评阅论文的各位专家和老师。
附录
MATLAB程序实现操作的具体代码：
main,m
clc;clear all;format long
%%%%%%%%%%%%%%%%%%%%%%%%%%%
global alpha coe alpha=1.6;%fractional index
coe=1;%diffusion coefficients
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
T=1;%time
L=1;%space tau=0.02;%time step size
h=0.05;%space step size
ga=tau^(alpha)/h^2;%步长比
TT=0:tau:T;
XX=0:h:L;
N=length(TT)-1;M=length(XX)-1;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
W=zeros(1,N+1);
W(1)=1;
for k=1:N W(k+1)=(1-(2-alpha)/k)W(k);%Lubich approximation coefficients end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
U=zeros(M-1,N+1);%初始化，解矩阵
A=zeros(M-1);B=zeros(M-1);C=zeros(M-1);
A=A+diag(ones(1,M-1))(1+2coega)+diag(ones(1,M-2),-1)(-coega)+diag(ones(1,M-2),1)(-coega);
B=B+diag(ones(1,M-1))(-2)+diag(ones(1,M-2),-1)+diag(ones(1,M-2),1);C=C+diag(ones(1,M-1))(1-2coegaW(2))+diag(ones(1,M-2),-1)coegaW(2)+diag(ones(1,M-2),1)coegaW(2);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
for n=0:N-1
tempMid=zeros(M-1,1);
for k=2:n+1
下载提示：
1、如文档侵犯商业秘密、侵犯著作权、侵犯人身权等，请点击“文章版权申述”（推荐），也可以打举报电话：18735597641(电话支持时间：9:00-18:30)。
2、网站文档一经付费（服务费），不意味着购买了该文档的版权，仅供个人/单位学习、研究之用，不得用于商业用途，未经授权，严禁复制、发行、汇编、翻译或者网络传播等，侵权必究。
3、本站所有内容均由合作方或网友投稿，本站不对文档的完整性、权威性及其观点立场正确性做任何保证或承诺！文档内容仅供研究参考，付费前请自行鉴别。如您付费，意味着您自己接受本站规则且自行承担风险，本站不退款、不进行额外附加服务。
原创文章，作者：写文章小能手，如若转载，请注明出处：https://www.447766.cn/chachong/14774.html,
