重心插值配点法求解分数阶Fredholm积分方程

虎晓燕, 韩惠丽

(宁夏大学 数学统计学院 宁夏 银川 750021)

重心插值配点法求解分数阶Fredholm积分方程

虎晓燕, 韩惠丽

(宁夏大学 数学统计学院 宁夏 银川 750021)

重心插值配点法是插值法和配点法的结合和推广，它具有稳定性好、高精度和计算效率高等优点.主要运用高精度无网格重心插值配点法求解分数阶Fredholm积分方程. 首先推导了基于分数阶Fredholm积分方程重心插值配点法的离散公式，然后通过理论分析得出其解的存在唯一性与误差分析，最后利用数值算例通过对等距节点与第二类Chebyshev节点的对比，验证了所用方法的高精度和可靠性，并得出影响精度的条件.

重心插值配点法； 高精度； 分数阶积分方程； 无网格

0 引言

分数阶微积分有着广阔的应用领域和前景，如生物传热、粘弹性力学、图像信号和风力发电等. 由于求分数阶微积分方程的解析解存在诸多困难，所以重点研究它的数值解. 目前微积分方程的数值解法较多，如有限差分法[1]、有限元法[2]、小波方法[3]、拉格朗日乘子法[4]、变分迭代法[5]、同伦摄动法[6]等，也有与统计力学相结合的随机行走方法[7]等. 上述方法较多研究整数阶微积分方程，而分数阶微积分方程更复杂，精确数值模拟分数阶微积分方程的困难更大，现有文献中研究分数阶微积分方程数值解的理论较少，而且很少考虑计算精度或收敛速度等问题. 高精度无网格重心插值配点法是一种新型求解分数阶积分方程的数值方法[8]，该方法相对其他数值方法而言应用前景广阔，克服了一般插值方法的振荡性和向前稳定性，具有计算精度高，不用划分网格，程序容易实现等诸多优点.重心插值配点法多用于数值模拟积分微分方程和偏微分方程. 文献[9]运用线性重心有理插值配置法求解了一阶线性双曲型初值边值问题.文献[10]运用两种重心型插值求解了带初值问题和边值问题的二阶线性微分方程，并将重心插值法与拟谱法以及微分求积法做了比较，得出重心插值法比其他数值方法精度更高.文献[11]使用直接线性重心有理插值(DRQ)求解了线性Fredholm积分函数，用间接线性重心有理插值法(IRQ)求解了线性Volterra积分函数，并给出了误差估计，从而说明了重心插值法的高效性.

与以往插值配点法不同，本文将运用重心插值配点法求解分数阶Fredholm积分方程，从而推导出分数阶Fredholm积分方程重心插值配点法的离散计算公式. 通过证明解的存在唯一性、误差估计以及数值算例验证本文方法的高精度和可靠性.

1 重心插值配点法

设xj(j=0,1,…,n+1)是函数f(x)的n+1个不同的插值节点，其对应的函数值为fj.如果要使用多项式插值，则可以在多项式空间中求插值多项式p(x)，使得p(xj)=fj,j=0,1,…,n，而多项式空间的次数不超过n.

定义1[12]重心Lagrange插值公式

(2)

由式(1)可知,分子和分母都包含重心权ωj,由式(2)可知,重心权只依赖于插值节点的分布，因此插值节点直接影响着重心权.本文主要采用等距节点和第二类Chebyshev节点来研究分数阶积分方程的数值解.

2 分数阶积分方程

定义2[13]设函数f(t)为定义在J0=(0,+)上分段连续，且在J=[0,+)的任意有限子区间上可积. 对∀t>0、Reν>0的任意复变量ν，函数f(t)的ν阶R-L分数阶定义,其中：ν>0；Γ(·)为Gamma函数.

3 格式推导

本文主要探讨如下形式的分数阶Fredholm积分方程，

(3)

其中：f(x)为未知函数；g(x)∈C[a,b]为已知函数；K(x,t)为区间[a,b]×[a,b]上的、关于变量x,t的二元连续积分核函数.将积分方程的积分区间[a,b]离散为a=x0

(4)

(5)

令方程(5)在节点x0,x1,x2,…,xn上均成立，可得

(6)

(7)

将式(7)写成矩阵形式为

(8)

其中:K=[Kij]:=[Kj(xi)]为未知函数f(x)在节点x0,x1,x2,…,xn的积分矩阵；I为n阶单位矩阵；F=[f(x0),f(x1),f(x2),…,f(xn)]T为f(x)在插值点处函数值；G=[g(x0),g(x1),g(x2),…,g(xn)]T为自由项g(x)所对应的插值点处的函数值. 这样式(8)就给出了积分方程(3)的重心插值配点法的离散计算公式.

由Gauss公式[14]，对给定的xi有

4 误差分析

4.1 解的存在唯一性

引理1[15](Fredholm选择定理)第二类Fredholm积分方程，要么对于所有函数f(x)都有解且解唯一，要么相应的齐次方程至少有一个非平凡解，即至少有一个不恒等于零的解.

引理2[16](压缩映射原理)设(X,d)是完备的距离空间，T：X→X是压缩映射.则T在X中存在唯一的不动点x0.

由引理1可知，本文所研究的方程必然有解. 在Banach空间中，设核K(x,t)在矩形区域D1=[a,b]×[a,b]上连续，自由项g(x)∈C[a,b]，下面利用压缩映射原理(引理2)证明其解的唯一性.

证明 定义映射T为

则T:C[a,b]→C[a,b],由于x-t>0,从而(x-t)α-1>0,且对于任意的f1,f2∈C[a,b]，有

4.2 误差估计

本文主要采用两种范数误差来刻画数值解与精确解的逼近程度问题，其中L误差定义为‖ek‖，L2误差定义为.

下面主要运用L误差研究本文方法的数值逼近问题.

(9)

若f(x)在[a,b]上存在n+1阶导数，则对∀x∈[a,b]，其重心插值有估计[17]

(10)

证毕.

5 数值算例

为了验证本文方法的高精度和可靠性，对2个具有精确解的数值算例进行数值模拟.算例中均分别选取等距节点和第二类Chebyshev节点，并对这两种节点下数值结果的L误差、L2误差和相对误差做比较. 下述算例的程序均在Matlab7.0中实现.

算例1

当α=1，取g(x)=11x/(12-x2)时,该分数阶方程退化为整数阶积分方程，精确解为f(x)=x(1-x). 表1为本文方法当α=1，n=10，q=8时两类节点和文献[18]方法的L误差、全局误差和相对误差对比，图1为α=1，n=10，q=8时两类节点误差与数值解、精确解比较. 当α=5/2时，适当选取g(x)使算例1精确解为f(x)=x(1-x).表2为α=5/2,n=10，q取不同值时两类节点的数值解、精确解、L误差、L2误差和相对误差. 其中全局误差定义为‖e‖2/n.

表1 算例1 本文方法当α=1，n=10，q=8时两类节点和文献[18]方法的L误差、全局误差和相对误差对比

Tab.1 The L errors, global errors and relative errors with α=1，n=10，q=8 for the present method and the paper [18] of two knots for example 1

表1 算例1 本文方法当α=1，n=10，q=8时两类节点和文献[18]方法的L误差、全局误差和相对误差对比

误差本文方法等距节点(10节点)Chebyshev节点(10节点)文献[18](64节点)L¥error1.387778780781446e-169.714451465470120e-173.113e-6全局误差6.081955556155792e-184.335447487907081e-185.505e-7相对误差1.053478278494604e-168.955261823477816e-17

图1 算例1当α=1，n=10，q=8时两类 节点的误差与数值解精确解比较

由表1、表2和图1可知，用重心插值配点法求解分数阶Fredholm积分方程时,α增大，q增大时，其计算精度随之增大. 当q由10变为100时，精度增加了3个数量级，当q由100变为500时，精度增加了2个数量级，这说明q影响数值解的计算精度为,q越大计算效果越好，但相应计算时间增大. 两类节点所得数值解的计算精度都比较高，其中第二类Chebyshev节点比等距节点的计算精度更高. 当α=1时，算例1退化为整数阶时，q取10时相对误差达到10-16，与文献[18]的方法相比，本文方法的计算精度更高.

算例2

当α=5/2时，选取适当的g(x)使算例2的精确解为f(x)=x3(1-x)ex. 具体计算中，分别取11个和21个插值节点，即分为n=10和n=20两种情况，q取10. 由表3可以看出，Gauss积分点数量q相同时，n=20比n=10计算效果更好，等距节点比第二类Chebyshev节点更好.

表2 算例1当α=5/2时两类节点的数值解、精确解、L误差、L2误差和相对误差

Tab.2 The computed solution,the exact solution,L errors,L2 errors and relative error with α=5/2 of two knots for example 1

表2 算例1当α=5/2时两类节点的数值解、精确解、L误差、L2误差和相对误差

x取值等距节点精确解数值解q=10数值解q=100数值解q=500Chebyshev节点精确解数值解q=10数值解q=100数值解q=5000.10.090.0900000.0899990.0899990.0238720.0238720.0238720.0238720.20.160.1600000.1599990.1599950.0863720.0863720.0863720.0863720.30.210.2099970.2099990.2099570.1636270.1636270.1636270.1636270.40.240.2400020.2399990.2398060.2261270.2261290.2261270.2261270.50.250.2500080.2500000.2493850.2500000.2500070.2500000.2500000.60.240.2399860.2399990.2384510.2261270.2261380.2261270.2261270.70.210.2100040.2100000.2066780.1636270.1636350.1636270.1636270.80.160.1600090.1599990.1536870.0863720.0863750.0863720.0863720.90.090.0900020.0899990.0791030.0238720.0238730.0238720.0238721.00.000.0000000.000000-0.0173600.0000000.0000000.000000L¥误差1.6427e-55.9319e-81.1255e-91.3723e-54.6133e-88.330e-10L2误差2.9879e-61.5451e-84.5119e-105.1144e-61.6984e-83.054e-10相对误差5.1755e-62.6764e-87.8153e-101.0564e-53.5083e-86.310e-10

表3 算例2当α=5/2时两类节点的L误差、L2误差和相对误差

Tab.3 The L errors, the L2 errors and the relative error with α=5/2 of two knots for example 2

表3 算例2当α=5/2时两类节点的L误差、L2误差和相对误差

误差数值解(q=10,n=10)等距节点Chebyshev节点数值解(q=10,n=20)等距节点Chebyshev节点L¥误差0.0010260.0010260.0010250.001026L2误差2.546e-43.416e-42.234e-43.257e-4相对误差6.242e-49.3661e-43.867e-46.314e-4

综合上述2个算例可知，应用重心插值配点法求解分数阶Fredholm积分方程时，可以得到很高的计算精度.计算精度随着阶数α的增大而增大，并且也受Gauss积分点数量q与插值节点个数n的影响. 当α=1，即积分方程由分数阶退化为整数阶，q与n取很小时就能得到很高的计算精度. 当α≠1，α越大，q越大，n越大时，分数阶积分方程的计算精度越高. 同时也可知，对于分数阶积分方程的精确解为基本初等函数时，第二类Chebyshev 节点比等距节点的精度更高，对于其精确解为更加复杂的函数时，等距节点比第二类Chebyshev 节点的精度更高一些.

6 小结

本文应用重心插值配点法求解了分数阶第二类Fredholm积分方程，推导出了分数阶Fredholm积分方程重心插值配点法的离散计算公式，证明了本文所讨论方程的解的存在唯一性，得出了误差估计，最后利用数值算例通过对比等距节点与第二类Chebyshev节点的数值结果，验证了本文方法的高精度和可靠性，并得出影响精度的条件，也给出两类节点的适用范围. 本文利用重心插值配点法对分数阶Fredholm积分方程进行了研究，得到了较好的结果.

[1] 丁恒飞. 分数阶偏微分方程的有限差分方法[D]. 上海：上海大学，2014.

[2] AGRAWAL O P. A general finite element formulation for fractional variational problems[J]. Journal of mathematical analysis and applications，2008，337(1)： 1-12.

[3] HEYDARIA M H，HOOSHMANDASL M R，MOHAMMADI F，et al. Wavelets method for solving systems of nonlinear singular fractional Volterra integro-differential equations[J]. Communication in nonlinear science and numerical simulation，2014，19(1)：37-48.

[4] 吴晓，黄志刚，杨立军. 用拉格朗日乘子法求解双模量静不定结构[J]. 力学与实践，2013，35(6)：79-81.

[5] HETMANIOK E，NOWAK I，STOA D,et al. A study of the convergence of and error estimation for the homotopy perturbation method for the Volterra-Fredholm integral equations[J]. Applied mathematics letters，2013，26(26)：165-169.

[6] ZHAN J，CHEN L，FU Y, et al． A new random walk simulation model for study of diffusion behavior of single particle within two-dimensional space[J]. Journal of electrochemistry，2012，18(5)：427-436.

[7] WU G C. A fractional variational iteration method for solving fractional nonlinear differential equations[J]. Computers and mathematics with applications，2011，61(61)：2186-2190.

[8] 王兆清，李淑萍，唐炳涛.一维重心型插值：公式、算法和应用[J].山东建筑大学学报，2007，22(5)：448-453.

[9] BALTENSPERGER R, BERRUT J P. The linear rational collocation method[J]. Journal of computational and applied mathematics，2001，134(1/2)：243-258.

[10]李淑萍. 基于重心型插值的数值计算方法[J]. 山东科学，2010，23(4)：13-16.

[11]KLEIN G，BERRUT J P. Linear barycentric rational quadrature[J]. Bit numerical mathematics，2012，52(2)：407-424.

[12]李树忱，王兆清. 高精度无网格重心插值配点法:算法、程序及工程应用[M]. 北京：科学出版社，2012.

[13]OLDHAM K B，SPANIER J. The fractional calculus theory and applications of differentiation and integration to arbitrary order[M]. New York：Dover Publications Inc，2002.

[14]李岳生，黄友谦. 数值逼近[M]. 北京：人民教育出版社，1978.

[15]吕涛，黄晋. 积分方程的高精度算法[M]. 北京：科学出版社，2013.

[16]匡继昌. 实分析与泛函分析[M]. 北京：高等教育出版社，2002.

[17]ZOU L,LI C. Barycentric Lagrange blending rational interpolation based on Padè approximation[C]//Proceeding of the 2011 International Conference on Computational and Information Sciences.Washington,2011：1124-1127.

[18]张倩，韩惠丽，张盼盼. 基于有理 Haar小波求解分数阶第2类Fredholm积分方程 [J]. 江西师范大学学报(自然科学报)，2014，38(1)：47-50.

(责任编辑：方惠敏)

Barycentric Interpolation Collocation Method for Solving Fredholm Integral Equation of Fractional Order

HU Xiaoyan, HAN Huili

(SchoolofMathematicsandStatistics,NingxiaUniversity,Yinchuan750021,China)

Barycentric interpolation collocation method was the combination and promotion of interpolation method and the collocation method. It was a new kind of numerical calculation method which had excellent numerical stability, high precision,and computational efficiency. This new method was used to solve fractional order Fredholm integral equation, and was deduced discrete calculate formula of the barycentric interpolation collocation method of fractional order Fredholm integral equation. Through theoretical analysis, it derived its existence and uniqueness of solutions and error analysis. Finally,some numerical examples were used to contrast equidistant nodes and the second Chebyshev nodes to demonstrate the effectiveness and precision of this method, and then the conditions that affect the precision was obtained.

barycentric interpolation collocation method; high precision; fractional order integral equation； mesh free

2016-08-15

国家自然科学基金项目(11261041, 11261045).

虎晓燕(1989—)，女，宁夏固原人，硕士研究生，主要从事积分方程数值解的研究，E-mail：huxiaoyancai@163.com；通讯作者：韩惠丽(1972—)，女，宁夏银川人，教授，主要从事积分方程数值解的研究，E-mail：nxhan@126.com.

O175.6

A

1671-6841(2017)01-0017-07

10.13705/j.issn.1671-6841.2016203