分数阶微分方程数值解的一种逼近方法讲解

分数阶微分方程数值解的一种逼近方法摘要本文提出了一类分数阶微分方程（）的数值解方案在这种方法中被型分数阶导数所表现型分数阶导数的属性可以让一个分数阶微分方程减弱为一个型积分方程这样做了之后许多研究型积分方程的数值方法也可以应用于寻找的数值解本文总时间被划分为一组小区间,在两个连续区间中，用二次多项式逼近未知函数.这些近似被替换成转化的型积分方程由此获得一组方程这些方程的解提供了的解这种方法被应用于解决两种类型的线性和非线性用这里给出的方法得到的解能与解析解和其他方法的数值解较好的吻合.同时结果说明这种数值方法是稳定的..引言本文讨论分数阶微分方程的数值解.分数阶导数和分数阶积分近年来收到了广泛的关注.在许多实际应用中，分数阶导数和分数阶积分为考虑的系统提供了更加精确地模型.比如，分数阶导数已经被成功地运用到模拟许多粘性材料的依赖频率的阻尼行为年之前，和提出了这个领域已经被研究的工作的一个回顾，并且说明了半阶导数模型可以非常好地描述阻尼材料的频率以来.另一些学者说明了分数阶导数和分数阶积分在电化学过程，电解质极化，有色噪声，粘性材料和混沌领域的应用，和提出了分数阶导数和分数阶积分在一般固体力学，特定粘弹性阻尼模型中的应用的调查提出了分数阶微积分在生物工程的三个关键部分的回顾分数阶导数和分数阶积分在其他领域的应用以及相关的数学工具和技巧还可以在许多其他文献上找到.系统模型中分数阶导数的引进大多会导致分数阶微分方程的出现.对某些特定的分数阶微分方程在通常系统条件下的解，已经有几种方法被找到.这些方法包括，拉普拉斯变换，傅里叶变换，模态综合法和特征向量展开法，数值法以及基于积分公式的方法然而，这些方法中大多数不能被应用到非线性分数阶微分方程更进一步的，正如等人指出的，这些方法很多只能应用到特定类型的分数阶微分方程，并且人们并不知道他们能否被推广.并且，在很多作者的研究成果中，并没有出现系统性的收敛性分析.最近，对于能被应用到线性和非线性分数阶微分方程的数值稳定数值逼近技巧，人们的兴趣愈发浓厚.这些方法技巧大多利用了分数阶微分方程可以被减弱为型积分方程的特性因此，型积分方程的数值解法也可以应用到分数阶微分方程的解当中等人提出了分数阶微分方程数值解的一种方法，其中分别代表预测，校正和估计这样一来很多学者又推广了应用于常微分方程和分数阶微分方程的--型预测-校正格式.这种方法的提出也是利用分数阶微分方程可以被转化为型积分方程的特点这些作者同时提出了误差分析和用外推法改善数值精度的延伸和提出了一种阶数大于的分数阶微分方程的数值解法.在该公式中，阶大于1的分数阶微分方程被减弱为阶小于1的分数阶微分方程，然后用相应的数值解法解由此导出的系.统在所有这些方法当中，节点之间的未知函数用线性函数逼近提出了阶数大于的分数阶微分方程的数值解法这种方法要求就和它的导数在时间节点上连续.本文基于古典分数阶微分方程可以转化为型积分方程的特点也提出了一种数值方法来逼近分数阶微分方程的解.特别地，我们用二次逼近函数来建立这种算法，结果说明这种方法可以被应用到寻求分数阶微分方程的数值解.我们还通过两个例子，线性和非线性问题的解决，说明了这种方法的高效和准确，并且这种数值方法是稳定的..数值算法关于分数阶导数的定义已经出现有好几种，它们包括--和分数阶导数这里，我们规定使用导数其中，导数的定义是&，（。=^^:（"旧"""*伐门（工岫，n,1）其中，a是导数的阶数，n是比a大的最小的整数式早在世纪就在的论文中被提出，在的论文发表前一年它被所用然而，在文献中，被（）式所定义的分数阶导数作为导数被广泛认知在接下来的讨论中，我们考虑含有导数的初值问题在初始条件：产⑹=城？，，，n,下的解，其中，是任意函数，/g⑴是的阶导数，制%，…,n是指定初值条件.假设这个函数关于参数和积分区间都是连续的，并且对于它的第二个参数满足条件在纯数学中，-nn导数比导数应用更加广泛然而，这里考虑导数是因为以-nn导数为基础的分数阶微分方程要求在点的导数和积分为一般来说，这些条件的物理意义不是已知的，并且在实际应用中，他们是不可用的nI讨论了寻找在更一般的情况满足下初始条件的正确格式的问题在n的文章中，方程和3被证明可以等价描述为:v：.t；=gt：一士匚"一T'：c-~Lf-.T""：；d【，其中为fi（t）=比翁性

为了解释以二次多项式为基础的数值方法，我们假设我要求的是由(2式)定义的分数阶微分方程从到的积分为了达到这个目的，我们把时间等分成份，令，作为时间区间的每一个部分时间在网格点上被表示为号==04N同时假设的数值逼近值被网格点卬j=0,1，…,g<T所决定该方法的基本思想是在相邻的两个时间节点和t加+工上数值地获取函数的值，然后重复这个过程来接近所求积分直到取到终点为了便于接下来的讨论，我们规定如下记号:y(tj)=式也)々,目色)=虱也)%,f(=f(jKy(jh))々耳这里的方法需要对方程(4每值一步求两次积分值.这里有两种方法来达成目的.第一种用一些近似函数逼近，然后用一种数值方法确定式的积分值这里需要在未知积分的情况下对方m+t和%工作初始的估计第二种都用近似函数来显式地逼近和以及确定式的积分注意在这种情况下，□巧工会作为参数出现在函数当中本文利用的是第二种逼近方法现在，我们给出算法的详细思路首先我们需要确定，。，士工的值用二次插值函数可以在区间，］上逼近和以及是二2次f(t,v①尸理0是二2次其中，北七八％”%；；是函数在第个时间节点的值，m.’f：，插值函数，其中下标（,k代表在第,=0,…，步的第个近似函数我们首先确定在。和匕处的值把式带入式，并积分得到，住工）明—+友=d匚口/（士）九其中,匕二十.『：（「二-『厂中…工）"可以精确计算得到注意式需要知道在Q和七的值或者间接地说，的值为了得到外，在［，］上把，近似为£色¥（切君工=屯1/&1叫）曲（t）九，0其中；二二;♦1fl二二£「，『「］;,一:.”;・二二二二二是另一个二次插值函数函数中抽，=,，由下式给出，中。,」0=上铲，%人）=陪普，巾。,式0=分函数二；.Y=I23由相似的办法定义把（10代）到（4中），积分得到:丁：；三叫」-J二二；二小,其中，小三=士；；「.二］一［厂1中小切dj可以中%k一样被精确计算出来由可以得到后储的值为

(((l3-fi-5/2这里，我们充分利用了二次多项式的性质在非二次多项式的情况下，人将会有不同的参数.把代到)得到),巾］)=—d三0九一：3-7二)+】3注意到和是关于两个未知量先和先的方程，可以用-法，不动点迭代或者其他非线性方法求解这里，我们用-法求解这些非线性方程这个方法需要对外和先作一个初始的估计当a大于由处的斜率可以得到关于力和先的更好的估计然而，在这里，对于a和。我们对把这些变量的初值估计为加注意在每一次迭代式，时间步长取现在我们假定在"=：，•筋:处，的值是已知的，我们要求的是:二一1和光出+工处的值根据以上的逼近方法，比m+l和比m+工可以被表示为：V?m+2^.^->k=2ni^m,k+J=虱士加+3+a(亡.+2V?m+2^.^->k=2ni^m,k⑴九+乜&「4,一三(正,「江,一-江,一)+小,一公一储6)其中，」上：，k=2m,2m+2,2m+2,d;;1：.,k=2m,2m+1/2,2m+1是和二二=二，・二二用同样方法确定的系数.注意(15),(16)的积分可以被数值确定因为y(t)在「"：=【，.筋:处的值是已知的.这些方程含有两个

未知量；：二「和”:未知量；：二「和”:而他们可以通过法得到本文中，我们把”.作为F(t加+口和双〜十。的初值估计这样一来，方程就可以在需要的区间上求积分作为特殊情况，考虑如下非线性系统口，这种条件下，，=-丫『式和减弱为[alla12~\r^Zm+l'a21a2zJLvzm+Z-其中ana12~Q^2m+l/2aa12~Q^2m+l/2a21a22r+-2(22)(23)(22)(23)七=。垢-1)—{r’v&Li一T严一-'-Az-i~g二，，二nV：二=；:&—)一十一—『『丁⑺》C2m,2mVzm(17)是一组线性联立方程，可以用线性方法求解.请注意以下两个补充说明.第一，方程(1)只含有一个y(t).当y(t)是一个向量函数时，算法同样成立.不过，所有的y(t)和f(t,y(t))必须换成相应的近似向量函数.第二，算法需要保存所有算过的的y的值.这是很多分数阶微分方程的典型特征.这将会导致一些问题，特别是当y的维数和分数阶有限元系统一样大时.这种情况下，系统有临近的短期记忆，y(t)在过去一段时间长度的值可以忽略不计,以此来改善对存储空间的需求和计算效率.数值结果为了说明这种方法的效率,我们分别考虑一个线性和一个非线性的算例.讨论这些例子是因为他们解的逼近格式是可靠的，并且可以用其他数值方法求解.这样我们就可以把用这种方法得到的结果和解析解以及其他数值方法的结果相比较.例1线性方程在第一个例子中，我们考虑如下给出的线性方程:以『h=一『f：,O<a<2,且产0；=1丁：.0：二0初始条件仅当a是成立的解是：vt=4「吨，其中,EHf函数的阶

图a时的比较（解析解，本文数值方法）Errorinfunctiony(t)fora=0.75anddifferentvaluesofh1J..h=0.4h=0.2h=0Ah=0.05h=0.0250000000.89.832E-4-3.906E-4-1.509E-4-4.730E-5-1.423E-51.6-3.108E-4-2.113E-4-7.031E-5-2.139E-5-6.382E-62A-2.813E-4-L243E-4-3.955E-5-1.190E-5-3.537E-63.2-2.100E-4-8.003E-5-2.486E-5-7.433E-6-2.205E-64-1.564E-4-5.513E-5-L687E-5-5.022E-6-1.487E-64.8-1.191E-4-3.998E-5-L211E-5-3.592E-6-1.063E-65.6-9.287E-5一3.018E-5一9.069E-6-2.685E-6-7.939E-76.4-7.41IE-5-2.352E-5-7.029E-6-2.078E-6-6.139E-7表1a时在不同值下函数的误差对不同的a和可以得到很多结果，这里给出其中一些在各种情况下，我们另考虑这个区间是因为它接近a的系统的时间这里图比较了a时的解析解和二次数值方法在这种情况下，我们令注意到这两个结果几乎完全重合为了强调收敛性，表列出了当。分别等于0.，40.，20.，10.0和50.0时2的5结果.注意到随着步长的缩小，误差也如期望一样的缩小了在大部分节点误差比都非常接近3这表明误差阶为或5二「误差阶为或5二「0时的比较解析解，本文数值方法Table2Errorinfunctiony(t)fora=1.5differentvaluesofhtJh=0.4h=0.2A=0.1h=0.05h=0.0250000000.81.253E-31.664E-42.615E-54.415E-67.671E-7L64.999E-48.669E-51.598E-52.920E-65.267E-72A-1.346E-46.038E-63.683E-68.685E-71.716E-73.2-3.702E-4-3.438E-5-3.707E-6-4.766E-7-7.061E-84-2.915E-4-3.696E-5-5.525E-6-9.064E-7-1.554E-74.8-1.026E-4-2.108E-5-4.008E-6-7.327E-7-1.317E-75.64.302E-5-4.275E-6-1.647E-6-3.572E-7-6.824E-86.49.386E-55.095E-66.904E-8-4.760E-8-1.266E-8表a时在不同值下函数的误差图展示了a分别等于，，，和时的数值结果因为解析解和数值结果完全一致，因此图中没用画出解析解当a

时，精确解为注意到随着时，精确解为注意到随着a越来越接近，数值解越收敛到解析解已)即在极限情况下，分数阶微分方程的解接近整数阶导数的解更进一步地，我们给出了。的一系列结果a和。的结果是分离的，因为的斜率在a处有一个跃迁图3比较了在a,时的解析解和数值解两个结果又一次几乎完全重合为了突出收敛性，表给出了。,分别等于，，，和时的数值解正如之前观察到的一样，在这种情况下，随着步长的减小，误差也随之减小.在这种情况下，误差比接近，这表明误差阶为或)=这样一来，通过观察a和a的收敛结果，可以知道误差的收敛阶为a或h=o和a的收敛结果，可以知道误差的收敛阶为a或h=o即误差的阶不仅依赖于，还依赖于导数的阶a图a不同时的比较a.aa:图a.7,时的比较解析解，本文数值方法Table3Comparisonoferrorin>'(f)reportedin[35]andobtainedusingthisschemea=0.5a=1.5hRet[35]ThisschemeRef[35]Thisscheme0.1-1.3OE-O3-3.98E-04-5.46E-042.45E-050.05—3.93E—04-L43E-04-I.28E-044.22E-060.025-L26E-04-5.06E-05-3.04E-057.40E-07表本文数值方法和文献中的误差的比较图展示了a分别等于，，，和时的数值结果.又一次，因为解析解和数值结果完全一致，因此图中没用画出解析解.当a2时精确解为注意到随着a越来越接近，数值解逐渐收敛到整数阶导数的解.图2和图4展示的收敛结果十分重要，因为他们说明了在极限情况下，分数阶微分方程和他们的解逼近整数阶微分方程以及他们的解析解表比较了文献的解的误差和用本文数值方法在a和，的解的误差注意到本文的方法得到了更低阶的误差这是因为，这里的方法是一种高阶方法当a和取其他值时这种趋势也能显现出来.2例2非.线性方程在第二个例子中，我们考虑一个如下定义的非线性方程

⑷=SS—3著士T+汉。+1)+停或-t4)3-[y(t)]3/2,（。）=0,，（。）=。和之前一样，第二个初始条件仅适用于a和之前一样，第二个初始条件仅适用于a的（精2确9解）在文献（3中5已）给出注意到当a方程的解在处的斜率趋近于无穷因此，他可能导致在附近出现一个巨大的数值误差和7和7取不同值下的误差Table4Errorinfunctiony(t)fora=0.75and1.5,andfordifferentvaluesofhta=0.75x=1.5h=QAh=0.05h=0.02577=0.1h=0.05h=0.02500000000.1L967E-4-6.138E-6-5.584E-7L447E-44.013E-62.719E-70.2-1.374E-4-1.246E-5-1.037E-6L058E-47.124E-64.651E-70.33.060E-4-1.752E-5-1.401E-62.066E-49.188E-65.838E-70.4-2.519E-4-1.988E-5-1.531E-61.594E-49.576E-65.837E-70.51J17E-4-1.692E-5-1.213E-63.8O3E-57.308E-63.968E-70.6-1.220E-4-4.758E-6-1.381E-77.468E-51.152E-6-6.166E-80.7-7.829E-42.171E-52.093E-6-6.058E-4-1.033E-5-8.899E-70.86.911E-46.861E-55.953E-6-2.928E-4-2.879E-5-2.201E-60.9-3.243E-31.425E-41.194E-5-2.054E-3-5.630E-5-4.134E-61.02.867E-32.488E-42.043E-5-1.161E-3-9.566E-5-6.872E-6Table5Comparisonoferroriny(f)reportedin[35]andobtainedusingthisschemeha—0.25a=1.25Ref.[35]ThisschemeRef.[35]Thisscheme0.12.50E-012.51E-03一5.53E-03-1.14E-040.05-1.50E-0