分数阶微分方程的数值解法概论
caputo分数阶微分方程求解 matlab 概述及解释说明

caputo分数阶微分方程求解matlab 概述及解释说明1. 引言1.1 概述在科学和工程领域中,微分方程是一种常见的数学模型,用于描述物质或现象之间的相互关系。
传统的微分方程主要基于整数阶导数进行建模和求解。
然而,许多现实中的问题不能仅用整数阶微分方程来完全描述,因此引入了分数阶微积分的概念。
Caputo分数阶微积分是世界上最早发表的一种分数阶导数定义方法之一,它在描述长尾动力学、非平衡统计物理、带记忆材料等领域具有广泛应用。
使用Caputo分数阶微积分可以更准确地对现实世界中各种复杂过程进行建模和仿真。
1.2 文章结构本文将首先介绍Caputo分数阶微积分的基本概念和定义,然后重点关注Caputo分数阶微分方程及其特性。
接下来,我们将详细探讨MATLAB在求解Caputo分数阶微分方程中所起到的关键作用,并提供实际示例以说明其应用方法和步骤。
随后,我们将选择一个具体的Caputo分数阶微分方程案例进行研究和求解,并通过结果及讨论来评估算法的效率。
最后,我们将对本文进行总结,并提出现有问题和未来工作方向的展望。
1.3 目的本文的主要目的是介绍Caputo分数阶微分方程在MATLAB中的求解方法,并通过案例研究和讨论来验证其有效性和实用性。
通过本文的阐述,读者将能够理解Caputo分数阶微积分的基本概念、MATLAB在求解Caputo分数阶微分方程中所采用的方法以及其应用领域。
此外,本文还旨在鼓励读者进一步研究该领域并探索新的解决方案。
2. Caputo分数阶微分方程概述:2.1 分数阶微积分简介分数阶微积分是传统整数阶微积分的推广,它引入了非整数阶导数和非整数阶积分的概念。
与整数阶微积分不同,分数阶导数和积分可以表现出一种记忆性的特点,使得在描述复杂自然现象、非线性动力学系统、多尺度问题等方面具有更好的适用性。
2.2 Caputo分数阶导数定义与性质Caputo导数是一种常用的描述物理过程中记忆效应的方法。
分数阶微分方程的高阶算法及理论分析

分数阶微分方程的高阶算法及理论分析分数阶微分方程的高阶算法及理论分析摘要:分数阶微分方程作为一种广义的微分方程,其在描述非平稳现象和复杂动力系统中具有重要的作用。
本文主要介绍了分数阶微分方程的高阶算法及其理论分析方法。
首先,对分数阶微分方程的定义和性质进行了简要的介绍,然后分别介绍了分数阶微分方程的常见解法及高阶算法,包括分数阶欧拉法、分数阶中点法、分数阶四阶Runge-Kutta法等。
随后,通过对这些高阶算法的理论分析,比较了它们在精度和稳定性方面的差异与优劣。
最后,通过数值实验验证了这些高阶算法的有效性和优越性。
关键词:分数阶微分方程、高阶算法、理论分析、数值实验1. 引言分数阶微分方程是一类常微分方程的广义形式,它与整数阶微分方程相比具有更广泛的应用领域和更丰富的数学性质。
尽管分数阶微分方程在理论研究和应用方面都具有重要意义,但由于其复杂的性质和数值方法的限制,研究分数阶微分方程的高阶算法成为一个热点问题。
2. 分数阶微分方程的定义和性质分数阶微分方程是指微分方程中导数阶数为分数的微分方程。
与整数阶微分方程不同,分数阶微分方程具有非局部性、记忆性和非马尔可夫性等特点。
这些特点使得分数阶微分方程对数学建模和实际问题的描述更加准确和全面。
3. 分数阶微分方程的常见解法为了描述和求解分数阶微分方程,已经发展了一系列的解法,包括变换法、数值解法、分数阶微积分和格林函数方法等。
这些解法在不同的问题和应用场景中具有不同的优势和适用性。
4. 分数阶微分方程的高阶算法4.1 分数阶欧拉法分数阶欧拉法是对经典的欧拉法进行推广得到的。
它采用离散化的方式逐步逼近分数阶导数,通过递推关系式计算逼近解,并能够保持一定的精度和稳定性。
4.2 分数阶中点法分数阶中点法是对经典的中点法进行推广得到的。
它通过对时间步长的二阶近似,能够更准确地逼近分数阶导数,提高了数值解的精度和稳定性。
4.3 分数阶四阶Runge-Kutta法分数阶四阶Runge-Kutta法是对经典的四阶Runge-Kutta法进行推广得到的。
第一讲分数阶微分方程

RL a
Dαx
f
(x)
Dn
(
a
Dαx −n
f
) (x)
=
1 Γ(n − α)
dn dxn
ˆx
a
f (t) (x − t)α−n+1
dt,
(1.6)
即先做 n − α 次分数阶积分, 然后再求 n 次导数. 我们注意到 0 < n − α ≤ 1.
RL a
Dαx
a
D−x α
f (x)
=
Dn
a
Dαx −n
a
D−x α
f (x)
=
Dn
a
D−x n
f (x)
=
f (x).
如果将次序反过来的话, 则有下面的复合公式.
定理 1.4 设 α > 0, 且 n − 1 ≤ α < n, 则
a
D−x α
RL a
Dαx f (x)
=
f (x)
∑n −
[a Dαx−i f (t)]t=a Γ(α − i + 1)
(
a
D−x 1
f
) (x)
=
dn dxn
f
(x),
因此, 当 α 是正整数时, R-L 分数阶导数与整数阶导数的定义是一致的. 所以, R-L 分数阶导数在整数阶导 数之间架起了 “桥”.
1.1.3 Caputo 分数阶导数
R-L 分数阶导数是最先提出来的, 理论分析也相对完善. 但与实际应用却存在一定的困难和障碍 [4, page 4]. 一个比较好的解决方法就是由 Caputo [1, 2] 提出来的 Caputo 分数阶导数.
随机分数阶微分方程

随机分数阶微分方程
随机分数阶微分方程是指带有随机项的分数阶微分方程,其形式如下:
$$D^\alpha y(t)=f(t,y(t))+g(t,y(t),\omega(t))$$
其中,$D^\alpha$表示分数阶导数,$f(t,y(t))$表示确定性项,
$g(t,y(t),\omega(t))$表示随机项,$\omega(t)$表示随机过程。
随机分数阶微分方程的解析解很难求得,因此需要采用数值方法进行
求解。
常用的方法有如下几种:
1. 分步法
分步法是将时间区间划分为若干个小段,然后在每个小段内用欧拉方
法或其他数值方法求解微分方程。
由于分数阶导数的计算比较困难,
通常采用格洛特曼公式或广义差分公式来逼近分数阶导数。
对于随机项,可以采用欧拉-玛鲁雅马方法或隐式方法进行求解。
2. 基于精确求导公式的方法
由于分数阶导数的计算比较困难,近年来出现了一些基于精确求导公
式的方法,可以准确计算分数阶导数。
这些方法包括Grünwald-Letnikov求导公式、Caputo求导公式等。
然后将求得的分数阶导数带入微分方程中,用离散化方法求解。
3. 基于模拟的方法
基于模拟的方法是通过模拟随机过程的轨迹来得到方程的数值解。
常见的方法包括蒙特卡洛方法、随机扰动法等。
这些方法可以通过增加模拟次数来提高精度,但计算量较大。
总之,随机分数阶微分方程的数值求解需要结合分数阶导数的逼近方法和随机项的数值模拟方法,应根据具体问题选择合适的数值方法。
同时,需要注意初始条件的选取和数值方法的稳定性,防止数值求解过程中出现震荡和发散。
分数阶微积分方程的一种数值解法

有限元法:适 用于求解复杂 的几何形状和
边界条件
谱方法:适用 于求解高维问 题和高阶微分
方程
迭代法:适用 于求解非线性
问题
算法优化策略
减少计算量: 采用高效的算 法和数值方法, 降低计算复杂
度
提高精度:采 用高阶数值方 法,减小误差, 提高解的精度
加速收敛:采 用加速收敛技 术,如松弛法、 预处理共轭梯 度法等,加快 迭代收敛速度
误差传播:误差在计算过程中 的积累和传播
误差估计:对数值解的精度进 行评估和预测
误差控制:采用适当的算法和 技术减小误差,提高数值解的 精度
05 数值解法的应用实例
实际问题建模
描述实际问题,建立数学模型 确定模型参数和变量 利用数值解法求解模型 分析结果,给出实际解决方案
数值模拟结果
展示了数值解法的实际效果 和精度
稳定性优化: 采用稳定性好 的算法和数值 方法,提高解
的稳定性
并行计算的应用
分数阶微积分方程数值解法的并行计算框架 并行计算在提高数值解法效率方面的优势 并行计算在优化数值解法精度方面的作用 并行计算在处理大规模分数阶微积分方程时的表现
误差分析
数值解法的误差来源:离散化 误差、舍入误差和截断误差
有限差分法
定义:有限差 分法是一种数 值求解偏微分 方程的方法, 通过将微分转 化为差分来近
似求解。
原理:基于泰 勒级数展开, 将微分算子近 似为离散的差 分算子,从而 将微分方程转 化为差分方程。
适用范围:适 用于规则区域, 如矩形、立方 体等,对于不 规则区域需要 进行适当的网
格划分。
优点:计算简 单、易于编程 实现、适合大
应用领域:广泛用于工程领域中的 各种微分方程的数值求解问题,如 结构分析、流体动力学、热传导等
分数阶微分方程课件

分数阶微分方程第三讲分数阶微分方程基本理论一、分数阶微分方程的出现背景及研究现状1、出现背景分数阶微积分是关于任意阶微分和积分的理论,它与整数阶微积分是统一的,是整数阶微积分的推广。
整数阶微积分作为描述经典物理及相关学科理论的解析数学工具已为人们普遍接受,很多问题的数学模型最终都可以归结为整数阶微分方程的定解问题,其无论在理论分析还是数值求解方面都已有较完善的理论。
但当人们进入到复杂系统和复杂现象的研究时,经典整数阶微积分方程对这些系统的描述将遇到以下问题:(1)需要构造非线性方程,并引入一些人为的经验参数和与实际不符的假设条件;(2)因材料或外界条件的微小改变就需要构造新的模型;(3)这些非线性模型无论是理论求解还是数值求解都非常繁琐。
基于以上原因,人们迫切期待着有一种可用的数学工具和可依据的基本原理来对这些复杂系统进行建模。
分数阶微积分方程非常适合于刻画具有记忆和遗传性质的材料和过程,其对复杂系统的描述具有建模简单、参数物理意义清楚、描述准确等优势,因而成为复杂力学与物理过程数学建模的重要工具之一。
2、研究现状在近三个世纪里,对分数阶微积分理论的研究主要在数学的纯理论领域里进行,似乎它只对数学家们有用。
然而在近几十年来,分数阶微分方程越来越多的被用来描述光学和热学系统、流变学及材料和力学系统、信号处理和系统识别、控制和机器人及其他应用领域中的问题。
分数阶微积分理论也受到越来越多的国内外学者的广泛关注,特别是从实际问题抽象出来的分数阶微分方程成为很多数学工作者的研究热点。
随着分数阶微分方程在越来越多的科学领域里出现,无论对分数阶微分方程的理论分析还是数值计算的研究都显得尤为迫切。
然而由于分数阶微分是拟微分算子,它的保记忆性(非局部性)对现实问题进行了优美刻画的同时,也给我们的分析和计算造成很大困难。
在理论研究方面,几乎所有结果全都假定了满足李氏条件,而且证明方法也和经典微积分方程一样,换句话说,这些工作基本上可以说只是经典微积分方程理论的一个延拓。
分数阶微积分方程

分数阶微积分是一个数学概念,它扩展了整数阶微积分到分数和复数领域。
分数阶微积分方程是分数阶微积分在数学建模中的应用,可以描述许多自然现象和工程问题。
分数阶微积分方程的一般形式为:
Df(x)=f(x)x∈[a,b]Df(x)=f(x) \quad x \in [a, b]Df(x)=f(x)x∈[a,b]
其中 D 是分数阶导数,f(x) 是待求解的函数,a 和 b 是定义域的上下限。
分数阶微积分方程的解法通常包括以下步骤:
1. 确定方程的形式和参数。
2. 确定初始条件和边界条件。
3. 使用数值方法求解方程,例如有限差分法、有限元法等。
4. 对解进行后处理,例如误差分析、可视化等。
分数阶微积分方程在许多领域都有应用,例如物理学、工程学、生物学等。
例如,在物理学中,它可以描述粘弹性材料的力学行为;在工程学中,它可以描述信号处理和控制系统;在生物学中,它可以
描述神经传导和扩散过程。
分数阶微分方程数值解的一种逼近方法讲解

分数阶微分方程数值解的一种逼近方法By:Pankaj Kumar, Om Prakash Agrawal摘要本文提出了一类分数阶微分方程(FDEs)的数值解方案.在这种方法中,FDEs 被Caputo型分数阶导数所表现. Caputo型分数阶导数的属性可以让一个分数阶微分方程减弱为一个Volterra型积分方程. 这样做了之后,许多研究Volterra 型积分方程的数值方法也可以应用于寻找FDEs的数值解. 本文总时间被划分为一组小区间,在两个连续区间中,用二次多项式逼近未知函数. 这些近似被替换成转化的Volterra型积分方程由此获得一组方程. 这些方程的解提供了FDE的解. 这种方法被应用于解决两种类型的FDEs,线性和非线性. 用这里给出的方法得到的解能与解析解和其他方法的数值解较好的吻合. 同时结果说明这种数值方法是稳定的.1.引言本文讨论分数阶微分方程的数值解. 分数阶导数和分数阶积分近年来收到了广泛的关注. 在许多实际应用中,分数阶导数和分数阶积分为考虑的系统提供了更加精确地模型. 比如,分数阶导数已经被成功地运用到模拟许多粘性材料的依赖频率的阻尼行为.1980年之前,Bagley 和Torvik提出了这个领域已经被研究的工作的一个回顾,并且说明了半阶导数模型可以非常好地描述阻尼材料的频率以来. 另一些学者说明了分数阶导数和分数阶积分在电化学过程,电解质极化,有色噪声,粘性材料和混沌领域的应用. Mainardi,Rossikhin和Shitikova 提出了分数阶导数和分数阶积分在一般固体力学,特定粘弹性阻尼模型中的应用的调查. Magin提出了分数阶微积分在生物工程的三个关键部分的回顾. 分数阶导数和分数阶积分在其他领域的应用以及相关的数学工具和技巧还可以在许多其他文献上找到.系统模型中分数阶导数的引进大多会导致分数阶微分方程的出现. 对某些特定的分数阶微分方程在通常系统条件下的解,已经有几种方法被找到. 这些方法包括,拉普拉斯变换,傅里叶变换,模态综合法和特征向量展开法,数值法以及基于Laguerre积分公式的方法. 然而,这些方法中大多数不能被应用到非线性分数阶微分方程. 更进一步的,正如Diethelm等人指出的,这些方法很多只能应用到特定类型的分数阶微分方程,并且人们并不知道他们能否被推广. 并且,在很多作者的研究成果中,并没有出现系统性的收敛性分析.最近,对于能被应用到线性和非线性分数阶微分方程的数值稳定数值逼近技巧,人们的兴趣愈发浓厚. 这些方法技巧大多利用了分数阶微分方程可以被减弱为Volterra型积分方程的特性. 因此,Volterra型积分方程的数值解法也可以应用到分数阶微分方程的解当中. Diethelm等人提出了分数阶微分方程数值解的一种PECE方法,其中P,C,E分别代表预测,校正和估计. 这样一来很多学者又推广了应用于常微分方程和分数阶微分方程的Adams–Bashforth–Moulton 型预测-校正格式. 这种方法的提出也是利用分数阶微分方程可以被转化为Volterra型积分方程的特点. 这些作者同时提出了误差分析和用Richardson外推法改善数值精度的延伸. Ford和Simpson提出了一种阶数大于1的分数阶微分方程的数值解法. 在该公式中,阶大于1的分数阶微分方程被减弱为阶小于1的分数阶微分方程,然后用相应的数值解法解由此导出的系统. 在所有这些方法当中,节点之间的未知函数用线性函数逼近. Kumar and Agrawal提出了阶数大于1的分数阶微分方程的数值解法. 这种方法要求就y(t)和它的导数在时间节点上连续.本文基于古典分数阶微分方程可以转化为Volterra型积分方程的特点也提出了一种数值方法来逼近分数阶微分方程的解. 特别地,我们用二次逼近函数来建立这种算法,结果说明这种方法可以被应用到寻求分数阶微分方程的数值解. 我们还通过两个例子,线性和非线性问题的解决,说明了这种方法的高效和准确,并且这种数值方法是稳定的.2.数值算法关于分数阶导数的定义已经出现有好几种,它们包括Riemann–Liouville, Grun-wald–Letnikov, Weyl, Caputo, Marchaud,和Riesz分数阶导数. 这里,我们规定使用Caputo导数.其中,Caputo导数的定义是, (n-1<α<n),(1)其中,α是导数的阶数,n是比α大的最小的整数.式(1)早在19世纪就在Liouville的论文中被提出,在Caputo的论文发表前一年它被Rabotnov所用. 然而,在文献中,被(1)式所定义的分数阶导数作为Caputo导数被广泛认知.在接下来的讨论中,我们考虑含有Caputo导数的初值问题:(2)在初始条件:, k=0,1,...,n-1,(3)下的解,其中,f是任意函数,是y的k阶导数,,k=0,1,…,n-1是指定初值条件. 假设这个函数关于参数和积分区间都是连续的,并且对于它的第二个参数满足Lipschitz条件.在纯数学中,Riemann–Liouville导数比Caputo导数应用更加广泛. 然而,这里考虑Caputo导数是因为以Riemann–Liouville导数为基础的分数阶微分方程要求y(t)在t=0点的导数和积分为0.一般来说,这些条件的物理意义不是已知的,并且在实际应用中,他们是不可用的. Lorenzo and Hartley讨论了寻找在更一般的情况满足下初始条件的正确格式的问题. 在Diethelm and Ford的文章中,方程(2)和(3)被证明可以等价描述为:,(4)其中g(t)为:. (5)为了解释以二次多项式为基础的数值方法,我们假设我要求的是由(2)式定义的分数阶微分方程从0到T的积分. 为了达到这个目的,我们把时间T等分成N 份,令h=T/N,作为时间区间的每一个部分. 时间在网格点上被表示为. 同时假设y(t)的数值逼近值被网格点所决定. 该方法的基本思想是在相邻的两个时间节点和上数值地获取函数y(t)的值,然后重复这个过程来接近所求积分直到取到终点T.为了便于接下来的讨论,我们规定如下记号:, ,这里的方法需要对方程(4)每一步求两次积分值. 这里有两种方法来达成目的.第一种,用一些近似函数逼近y(t),然后用一种数值方法确定式(4)的积分值.这里需要在未知积分的情况下对和作初始的估计. 第二种, 都用近似函数来显式地逼近y(t)和f(t,y(t))以及确定式(4)的积分. 注意在这种情况下,和会作为参数出现在函数f当中. 本文利用的是第二种逼近方法.现在,我们给出算法的详细思路. 首先我们需要确定y(t),,的值. 用二次插值函数可以在区间[0,]上逼近y(t)和f(t,y(t)):(6)以及(7)其中,是函数在第k个时间节点的值,,k=0,1,2是二次插值函数,其中下标(j,k)代表在第j+1,j=0,…,m步的第k个近似函数.我们首先确定y(t)在和处的值. 把(7)式带入(4)式,并积分得到:(8)其中,,(9)可以精确计算得到. 注意式(8)需要知道f在和的值(或者间接地说,y的值). 为了得到,在[,]上把f(t,y(t))近似为:,(10)其中,,是另一个二次插值函数. 函数,k=1,2,3由下式给出,函数由相似的办法定义.把(10)代到(4)中,积分得到:, (11)其中,, (12)可以(9)中一样被精确计算出来. 由(7)可以得到的值为:(13) 这里,我们充分利用了二次多项式的性质. 在非二次多项式的情况下,将会有不同的参数.把(13)代到(11)得到,+, (14)注意到(8)和(14)是关于两个未知量和的方程,可以用Newton–Raphson法,不动点迭代或者其他非线性方法求解. 这里,我们用Newton–Raphson法求解这些非线性方程. 这个方法需要对和作一个初始的估计. 当α大于1, 由t=0处的斜率可以得到关于和的更好的估计. 然而,在这里,对于α>1和α<1我们对把这些变量的初值估计为. 注意在每一次迭代式,时间步长取2h.现在我们假定在处,y的值是已知的,我们要求的是和处的值. 根据以上的逼近方法,和可以被表示为:+, (15)++(16)其中,,k=2m,2m+2,2m+2,,k=2m,2m+1/2,2m+1是和,用同样方法确定的系数. 注意(15),(16)的积分可以被数值确定因为y(t)在,处的值是已知的. 这些方程含有两个未知量和,而他们可以通过Newton–Raphson法得到. 本文中,我们把作为和的初值估计. 这样一来,方程(1)就可以在需要的区间上求积分.作为特殊情况,考虑如下非线性系统:.这种条件下,,式(16)和(15)减弱为:, (17) 其中=, (18)=, (19)=-, (20)=1+(21)=-, (22)=- . (23)(17)是一组线性联立方程,可以用线性方法求解.请注意以下两个补充说明. 第一,方程(1)只含有一个y(t). 当y(t)是一个向量函数时,算法同样成立.不过,所有的y(t)和f(t,y(t))必须换成相应的近似向量函数. 第二,算法需要保存所有算过的的y的值. 这是很多分数阶微分方程的典型特征. 这将会导致一些问题,特别是当y的维数和分数阶有限元系统一样大时. 这种情况下,系统有临近的短期记忆,y(t)在过去一段时间长度的值可以忽略不计,以此来改善对存储空间的需求和计算效率.3. 数值结果为了说明这种方法的效率,我们分别考虑一个线性和一个非线性的算例. 讨论这些例子是因为他们解的逼近格式是可靠的,并且可以用其他数值方法求解. 这样我们就可以把用这种方法得到的结果和解析解以及其他数值方法的结果相比较.3.1例1线性方程在第一个例子中,我们考虑如下给出的线性方程:,0<α<2,(24) 且. (25) 初始条件仅当α>1是成立. (24),(25)的解是:,(26) 其中,. (27) 是Mittag–Leffler函数的阶.图1.α=0.75时y(t)的比较(O:解析解,X:本文数值方法)表1.α=0.75时h在不同值下函数y(t)的误差对不同的α和h可以得到很多结果,这里给出其中一些. 在各种情况下,我们另T=6.4s.考虑这个区间是因为它接近α=2的系统的时间. 这里图(1)比较了α=0.75时的解析解和二次数值方法. 在这种情况下,我们令h=0.1s. 注意到这两个结果几乎完全重合. 为了强调收敛性,表(1)列出了当α=0.75,h分别等于0.4,0.2,0.1,0.05和0.025时的结果. 注意到随着步长的缩小,误差也如期望一样的缩小了. 在大部分节点误差比R=E(2h)/E(h)都非常接近3.37,这表明误差阶为1.75(或E(h)=O()).图2. α不同时y(t)的比较(O:α=0.25,X:α=0.5,+:α=0.75,Δ:α=0.95,*:α=1.)图3.α=1.5时y(t)的比较(O:解析解,X:本文数值方法)表2.α=1.5时h在不同值下函数y(t)的误差图(2)展示了h=0.025,α分别等于0.25,0.5,0.75,0.95和1时y(t)的数值结果.因为解析解和数值结果完全一致,因此图中没用画出解析解. 当α=1时,精确解为y(t)=. 注意到随着α越来越接近1,数值解越收敛到解析解y(t)=,即在极限情况下,分数阶微分方程的解接近整数阶导数的解.更进一步地,我们给出了α>1的一系列结果.α<1和α>1的结果是分离的,因为y(t)的斜率在α=1处有一个跃迁.图3比较了y(t)在α=1.5,h=0.4时的解析解和数值解. 两个结果又一次几乎完全重合.为了突出收敛性,表2给出了α=1.5,h分别等于0.4,0.2,0.1,0.05和0.025时的数值解. 正如之前观察到的一样,在这种情况下,随着步长的减小,误差也随之减小. 在这种情况下,误差比接近5.5,这表明误差阶为2.5(或E(h)=O()). 这样一来,通过观察α<1和α>1的收敛结果,可以知道误差的收敛阶为1+α(或E(h)=O()),即误差的阶不仅依赖于h,还依赖于导数的阶α.图4. α不同时y(t)的比较(O:α=1.25,X:α=1.5,+:α=1.75,Δ:α=1.95,*:α=2.)图5.α=0.25,0.75,1.25,1.75时y(t)的比较(O:解析解,X:本文数值方法)表3.本文数值方法和文献(35)中y(t)的误差的比较.图4展示了h=0.025,α分别等于1.25,1.5,1.75,1.95和2时y(t)的数值结果.又一次,因为解析解和数值结果完全一致,因此图中没用画出解析解. 当α=2时,精确解为y(t)=cos(t). 注意到随着α越来越接近2,数值解逐渐收敛到整数阶导数的解. 图2和图4展示的收敛结果十分重要,因为他们说明了在极限情况下,分数阶微分方程和他们的解逼近整数阶微分方程以及他们的解析解. 表3比较了t=1.0文献(35)的解的误差和用本文数值方法在α=0.5和1.25,h=0.1,0.05,0.025的解的误差. 注意到本文的方法得到了更低阶的误差. 这是因为,这里的方法是一种高阶方法. 当α和h取其他值时这种趋势也能显现出来.3.2 例2.非线性方程在第二个例子中,我们考虑一个如下定义的非线性方程:(28) 且. (29) 和之前一样,第二个初始条件仅适用于α>1. (28)(29)的精确解在文献(35)中已给出,(30)注意到当α<1, 方程的解在t=0处的斜率趋近于无穷. 因此,他可能导致在t=0附近出现一个巨大的数值误差.表4. α=0.75和1.5,h取不同值下y(t) 的误差.表5. 文献(35)中y(t)的误差和用本文数值方法得到的y(t)的误差的比较上面给出了一些在不同α和h下的数值结果. 图5表示了h=0.05,α分别等于0.25,0.75,1.25,1.75时解析解和数值解的结果. 由它可以得到(1).解析解和数值解基本重合,当α取其他值是,可以得到同样的结果. (2)尽管在t=0处斜率非常大,但是方法给出了非常精确的结果. (3)正如预期的那样, 在t=1处,对所有的α,y(t)的值收敛到0.25. 表4列出了当α=0.75和1.5,h=0.1,0.05和0.025的数值解的误差. 注意到误差随着步长的减小而减小. 同样的趋势在α取其他值时也能观察到. 在尝试过的α的值中,误差比R=E(2h)/E(h)表明没有特定的收敛阶. 然而,当α=0.75和1.5时,收敛的误差的平均值接近12和15.表5比较了文献(35)中在t=1.0处解的误差和用这种方法在α=0.25和1.25,h=0.05时的误差. 这里我们用的是文献(35)中用Richardson外推法得到的值. 观察得到,本文的方法又一次给出了更小的误差. 当α=0.25时,这种方法给出了比Richardson外推法小得多的误差. 这可能是因为,当α等于0.25时y(t)在t=0附近的斜率改变非常迅速,并且线性方法不能精确地获得结果. 需要指出的是,这种数值方法对于α和h取其他值时同样给出了更加精确的结果.4.结论本文给出了一种经典的分数阶微分方程的数值逼近方法. 这里的分数阶微分方程是依据Caputo型分数阶导数给出的, 这种导数的性质可以把分数阶微分方程减弱为Volterra型积分方程. 时间区间被分成一组网格,通过3个连续节点的二次插值多项式逼近未知的和已知的函数y(t)和f(t,y(t)). 把这些多项式带入Volterra型积分方程可以得到一组代数方程,这种数值方法的提出就是用来解这些方程以及获取y(t)的解. 通过一个线性一个非线性的例子的解决,说明了这种数值方法的作用. 用这种方法得到结果和解析解以及其他数值方法的结果是一致的. 还得到一个结论就是结果随着步长的减小而收敛. 在极限情况小,当α逼近整数值,数值方法会得到一个整数阶系统的解. 结果还表明这种方法是数值稳定的.。