# An Efficient Numerical Algorithm for Solving Fractional Higher-Order Nonlinear Integrodifferential Equations

Authors: Muhammed I. Syam; Qasem M. Al-Mdallal; M. Naim Anwar
Journal: Abstract and Applied Analysis
 (2015)

---

## Abstract

This paper is devoted to both theoretical and numerical study of boundary value problems for higher-order nonlinear fractional integrodifferential equations. Existence and uniqueness results for the considered problem are provided and proved. The numerical method of solution for the problem is based on a conjugate collocation and spline approach combined with shooting method. Some numerical examples are discussed to demonstrate the efficiency and the accuracy of the proposed algorithm.


---

## Body

## 1. Introduction

Within the context of fractional calculus, it is argued that anticipated sort of memory is being carried out from past states to current states; see, for example, the recent work of Agarwal et al. [1]. Therefore, in recent years, several phenomena in physics, chemistry, life sciences, geophysics and earth sciences, and fluid dynamics have been extensively investigated through mathematical models involving fractional calculus; see, for example, the works by Podlubny [2], Mainardi [3], and Kilbas et al. [4].In this paper we consider a class of boundary value problems for nonlinear integrodifferential equations of the form(1)Ly≔Dαyx+∫0xKx,tft,ydt+hx=0,x∈I=0,1,subject to(2a)y0=a0,y′0=a1(2b)y1=b0,y′1=b1,where3<α≤4, f∈C[I×R,R], K∈C[I×I,R+] is a positive kernel, h(x)∈C[I,R], and a1,b0,b1∈R. Here, Dα denotes the fractional differential operator of order α in Caputo’s sense and is given by(3)Dαyx=1Γk-α∫0xx-tk-α-1yktdt,where k=⌈α⌉ is the smallest integer greater than or equal to k. It is well-known that many mathematical models in engineering and other disciplines in science involve integrodifferential equations of fractional order, for example, problems in modeling of turbulent aerodynamic phenomena, continuum under viscoelastic situations, certain population dynamics problems, and heat transfer in composite materials with certain properties. Account for similar issues is given in [5–8] and the references therein.A survey of the literature reveals that theoretical and numerical investigations related to corresponding boundary value problems for fractional differential and integrodifferential equations are still in their early stages; see [9–17] and the references therein. Moreover, it is extremely difficult to find exact solutions of such problems; therefore, most researchers have focused on the numerical methods to approximate exact solutions. Examples of such methods are the Adomian decomposition method [18, 19], collocation spline method [20] and Zhao et al. [21], variational iteration method and homotopy perturbation method [22, 23], fractional differential transform method [24, 25], CAS wavelets [26], discrete Galerkin method [27], Chebyshev wavelets method adopted by Biazar and Ebrahimi [28], and Taylor expansion method [29].The rest of the paper is organized as follows: some definitions and preliminary results are presented in Section2. In Section 3, some relevant theoretical results such as the existence and uniqueness of the considered problem are presented. The numerical method of solution is presented in Section 4. In Section 5, numerical examples are discussed to demonstrate the efficiency and the rapid convergence of the present algorithm.


## 2. Definitions and Preliminary Results

This section presents some definitions and preliminary results that will be extensively used in this study. We first introduce the Riemann–Liouville definition of fractional derivative operator.Definition 1.
The left sided Riemann-Liouville fractional integral operator of orderα is defined by(4)Jαyx=1Γα∫axx-tα-1ytdt,where y∈L1(a,b) and α∈R+.The properties of the operatorJ(α) are summarized in the following lemma.Lemma 2.
Letα,β,x>0 and γ>-1. Then, see [2](i)
JαJβyx=Jα+βyx=JβJαyx,(ii)
J(α)x-aγ=Γ(γ+1)/Γ(γ+1+α)(x-a)γ+α.Note that the left sided Caputo fractional derivative (3) is originally defined via the left sided Riemann-Liouville fractional integral (4), as follows: (5)Dαyx=Jn-αynx=yx,x>0,where α∈R+, n=⌈α⌉, and y∈L1[a,b].Lemma 3.
Forα∈R+, n=⌈α⌉, and y∈L1[a,b], one has (i)
(D(α)J(α)y)(x)=y(x).(ii)
(J(α)D(α)y)(x)=y(x)-∑m=0n-1yma+(x-a)m/m!.(iii)
Dα(x-a)r=Γ(r+1)/Γ(r+1-α)(x-a)r-α, for [α]<r.


## 3. Analytical Results

The existence and uniqueness of the exact solution to problem (1) subject to the boundary conditions (2a) and (2b) are discussed herein. Since we are using the shooting method which requires converting the boundary value problem to initial value problem, we will discuss in the next theorem the existence of the solution to (1) subject to(6)y0=a0,y′0=a1,y′′0=a2,y′′′0=a3.Theorem 4 (existence).
Assume thath∈C[0,1], k∈C[0,1]2, and f∈C[0,1]×R. Then for any ϵ>0 and (7)χ=min⁡1,ϵΓα+1h∞+K∞f∞1/αthere exists y:[0,χ]→R solving the initial value problem (1) and (6).Proof.
Taking the Riemann–Liouville functional operator of both sides of (1) and applying Lemma 3, one obtains (8)yx=∑i=03yi0i!xi-1Γα∫0xhτ+∫0τKτ,tft,ydtx-τ1-αdτ.
LetB={y∈C0,χ:y-∑i=03y(i)(0)/i!∞<ϵ}. Obviously, B is a closed subset of the Banach space of all continuous functions on 0,χ equipped with the Chebyshev norm. Moreover, since y(x)=∑i=03y(i)(0)/i!xi for x∈[0,χ] is in B then B≠ϕ. Define the operator L on B by(9)Lyx≔∑i=03yi0i!xi-1Γα∫0xhτ+∫0τKτ,tft,ydtx-τ1-αdτ.The equation under consideration can be written as(10)Ly=y.Our aim is to show that (10) has a fixed point in B. Since h∈C[0,1], k∈C[0,1]2, and f∈C[0,1]×R, then L[y] is a continuous function. Notice that, since 2<α≤4 ensures that α-1>2. To achieve our target, we need to show that L is a self-mapping on B. For any y∈B and x∈[0,χ], we have(11)Lyx-∑i=03yi0i!xi=1Γα∫0xhτ+∫0τKτ,tft,ydtx-τ1-αdτ≤1Γα∫0xh∞+K∞f∞τx-τ1-αdτ.Since 0≤τ≤χ≤1, then(12)Lyx-∑i=03yi0i!xi≤1Γαh∞+K∞f∞∫0xx-τα-1dτ=1Γαh∞+K∞f∞xαα≤1Γα+1h∞+K∞f∞χα≤1Γα+1h∞+K∞f∞ϵΓα+1h∞+K∞f∞=ϵ.Therefore, L[y]∈B if y∈B; that is, L maps y to itself. Based on Banach’s fixed point theorem, it follows that the proof is complete.It is easy to see that the solution produced by Theorem4 depends on x, a2, and a3. To prove the existence of solution to problem (1)–(2b), it is enough to force the solution to satisfy the condition (2b). In this case, we can determine values of a2 and a3 so that the solution of (1)–(2b) depends on x only.Theorem 5 (uniqueness).
Letf be a Lipschitz function in the variable y with Lipschitz constant L. Let K be a bounded function on [0,1]×[0,1] such that (13)Kt1,t2≤M,∀t1,t2∈0,1×0,1,where M>0. If 0<LM/αΓ(α)<1; then problem (1)–(2b) has a unique solution.Proof.
Lety1 and y2 be two solutions to problem (1)–(2b); then(14)Dαy1x+∫0xKx,tft,y1dt+hx=0,(15)Dαy2x+∫0xKx,tft,y2dt+hx=0.Subtracting (14) from (15) then applying the Riemann–Liouville fractional integration one obtains(16)y2-y1=1Γα∫0x∫0τKτ,tft,y2-ft,y1dtx-τ1-αdτ≤LMΓαy2-y1∫0x1x-τ1-αdt≤LMαΓαy2-y1.Obviously, (17)1x-τ1-αdt=xαα≤1α.Since LM/αΓ(α)<1, then y1=y2 which completes the proof.


## 4. Method of Solution

The following is a brief derivation of the numerical algorithm used to solve problem (1) subject to (2a) and (2b). It is based on conjugate collocation approach with multiple shooting method. It consists of three main steps:(1)
Collocation step.(2)
Spline step.(3)
Multiple shooting step.


### 4.1. Collocation and Spline Methods

For the sake of simplicity, we discuss the solution of (1) as initial value problem with (18)y0=a0,y′0=a1,y′′0=a2,y′′′0=a3,where a2 and a3 are unknown constants which will be determined later.The intervalI=[0,1] is partitioned into N uniform subintervals λi=[xi-1,xi] (for i=1,2,…,N) of width h=1/N. Let ZN=xi=ih:i=1,2,…,N-1. For a given m≥1, let Sm+3(3)(ZN) be the spline space of piecewise polynomials on ZN which are 3 times continuously differentiable on the interval I, given by (19)Sm+33ZN=u∈C3I:uλi=uit∈Pm+3onλi,where Pm+3 represents the set of all real polynomials of degree not exceeding m+3. Notice that m represents the number of collocation points in each subinterval λi. Those points are defined as (20)Xi=xi,j=xi+hcj:i=0:N-1,j=1:m,with 0≤c1≤c2≤⋯≤cm≤1. Based on the collocation method, the exact solution y of problem (1) subject to (18) will be approximated by an element u∈Sm+3(3)(ZN) such that(21)Dαux+∫0xKx,tft,udt+hx=0,x∈X=⋃i=0N-1Xi,subject to(22)u0=a0,u′0=a1,u′′0=a2,u′′′0=a3.On each subinterval λi, the spline u can be expressed as a piecewise polynomial of degree m+3 of the form(23)ux=uixi+hχ=μ0i+μ1iχ+μ2iχ2+μ3iχ3+∑j=1mβjiχj+3,x∈λi,i:N-1,where χ∈I. Applying the results of Blank [30], we may evaluate the fractional differential operator of order q for the the collocation solution (23) at x=xi+cih as follows:(24)Dαuixi+cjh=NαΓ1-α∑k=0i∑r=03wj,ri-k,αμrk+∑k=0i∑s=0mwj,s+3i-k,αβsk,where (25)wi,kj,α=j+ci-α-δj,0∗j+ci-1-α,k=0j+ci-α+k∏l=1kll-α-δj,0∗j+ci-1ν-α∏l=1νk-ν+ll-α,k≥1and δj,0∗=0 if j=0, and 1 otherwise. Consequently, (21) for each Xi(i=0,1,…,N-1) can be expressed in the form(26)NαΓ1-α∑k=0i∑r=03wl,ri-k,αμrk+∑k=0i∑s=0mwl,s+3i-k,αβsk+∫0xi,lKxi,l,tft,∑q=03μqiclq+∑j=1mβjiclj+3dt+hxi,l=0,for l=1,2,…,N. Applying Simpson’s rule to approximate the integral in (26), one obtains(27)NαΓ1-α∑k=0i∑r=03wl,ri-k,αμrk+∑k=0i∑s=1mwl,s+3i-k,αβsk+xi,l3Kxi,l,0f0,∑q=03μqiclq+∑j=1mβjiclj+3+4xi,l3Kxi,l,xi,l2fxi,l2,∑q=03μqiclq+∑j=1mβjiclj+3+xi,l3Kxi,l,xi,lfxi,l,∑q=03μqiclq+∑j=1mβjiclj+3+hxi,l=0.Obviously, (27) can be written in the following matrix form:(28)giMi,Bi=Ri,i=0,1,…,N-1,where (29)gilMi,Bi=NαΓ1-α∑k=0i∑r=03wl,ri-k,αμrk+∑k=0i∑s=1mwl,s+3i-k,αβsk+xi,l3Kxi,l,0f0,∑q=03μqiclq+∑j=1mβjiclj+3+4xi,l3Kxi,l,xi,l2fxi,l2,∑q=03μqiclq+∑j=1mβjiclj+3+xi,l3Kxi,l,xi,lfxi,l,∑q=03μqiclq+∑j=1mβjiclj+3,l=1,2,…,N,M(i)=μ0(i),…,μ3(i)t, B(i)=β1(i),…,βm(i)t, and R(i)=-h(xi,1),…,-h(xi,N)t. Here ·t means the transpose of the vector. Based on the given definition of the approximate spline u, it can be easily verified that (30)M0=μ00μ10μ20μ30=a0a1a2a3.


### 4.2. Multiple Shooting Method

For simplicity, let us rewrite problem (1)–(2b) in the following form:(31)Ly≔Dαyx+∫0xKx,tft,ydt=-hx,subject to(32)y0=a0,y′0=a1,y1=b0,y′1=b1.The solution of problem (31)-(32) can be determined by solving the problem on the subintervals xi,xi+1, for i=0,1,…,N-1. In order to apply the shooting method, we introduce the following set of initial value problems:(33)Lyi=Dαyix+∫0xKx,tft,yidt=-hx,xi≤x≤xi+1subject to(34)yixi=a4i,yi′xi=a4i+1,yi′′xi=a4i+2,yi′′′xi=a4i+3,where a2,a3,…,a4N-1 are unknown real parameters. The solution of problem (33)-(34) will be obtained by the method described in the previous subsection, where the parameters a2,a3,…,a4N-1 are determined by solving the following algebraic system (35)y0nx1,a2:a3=y1nx1,a4:a7,⋮yi-1nxi,a4i-4:a4i-1=yinxi,a4i:a4i+3,⋮yN-11,a4N-4:a4N-1=b0,yN-1′1,a4N-4:a4N-1=b1,where n=0:3 and al:am represents al,al+1,…,am. It is worth mentioning that we use shooting method of order five (N=5) in our computations which are done using Matlab.


## 4.1. Collocation and Spline Methods

For the sake of simplicity, we discuss the solution of (1) as initial value problem with (18)y0=a0,y′0=a1,y′′0=a2,y′′′0=a3,where a2 and a3 are unknown constants which will be determined later.The intervalI=[0,1] is partitioned into N uniform subintervals λi=[xi-1,xi] (for i=1,2,…,N) of width h=1/N. Let ZN=xi=ih:i=1,2,…,N-1. For a given m≥1, let Sm+3(3)(ZN) be the spline space of piecewise polynomials on ZN which are 3 times continuously differentiable on the interval I, given by (19)Sm+33ZN=u∈C3I:uλi=uit∈Pm+3onλi,where Pm+3 represents the set of all real polynomials of degree not exceeding m+3. Notice that m represents the number of collocation points in each subinterval λi. Those points are defined as (20)Xi=xi,j=xi+hcj:i=0:N-1,j=1:m,with 0≤c1≤c2≤⋯≤cm≤1. Based on the collocation method, the exact solution y of problem (1) subject to (18) will be approximated by an element u∈Sm+3(3)(ZN) such that(21)Dαux+∫0xKx,tft,udt+hx=0,x∈X=⋃i=0N-1Xi,subject to(22)u0=a0,u′0=a1,u′′0=a2,u′′′0=a3.On each subinterval λi, the spline u can be expressed as a piecewise polynomial of degree m+3 of the form(23)ux=uixi+hχ=μ0i+μ1iχ+μ2iχ2+μ3iχ3+∑j=1mβjiχj+3,x∈λi,i:N-1,where χ∈I. Applying the results of Blank [30], we may evaluate the fractional differential operator of order q for the the collocation solution (23) at x=xi+cih as follows:(24)Dαuixi+cjh=NαΓ1-α∑k=0i∑r=03wj,ri-k,αμrk+∑k=0i∑s=0mwj,s+3i-k,αβsk,where (25)wi,kj,α=j+ci-α-δj,0∗j+ci-1-α,k=0j+ci-α+k∏l=1kll-α-δj,0∗j+ci-1ν-α∏l=1νk-ν+ll-α,k≥1and δj,0∗=0 if j=0, and 1 otherwise. Consequently, (21) for each Xi(i=0,1,…,N-1) can be expressed in the form(26)NαΓ1-α∑k=0i∑r=03wl,ri-k,αμrk+∑k=0i∑s=0mwl,s+3i-k,αβsk+∫0xi,lKxi,l,tft,∑q=03μqiclq+∑j=1mβjiclj+3dt+hxi,l=0,for l=1,2,…,N. Applying Simpson’s rule to approximate the integral in (26), one obtains(27)NαΓ1-α∑k=0i∑r=03wl,ri-k,αμrk+∑k=0i∑s=1mwl,s+3i-k,αβsk+xi,l3Kxi,l,0f0,∑q=03μqiclq+∑j=1mβjiclj+3+4xi,l3Kxi,l,xi,l2fxi,l2,∑q=03μqiclq+∑j=1mβjiclj+3+xi,l3Kxi,l,xi,lfxi,l,∑q=03μqiclq+∑j=1mβjiclj+3+hxi,l=0.Obviously, (27) can be written in the following matrix form:(28)giMi,Bi=Ri,i=0,1,…,N-1,where (29)gilMi,Bi=NαΓ1-α∑k=0i∑r=03wl,ri-k,αμrk+∑k=0i∑s=1mwl,s+3i-k,αβsk+xi,l3Kxi,l,0f0,∑q=03μqiclq+∑j=1mβjiclj+3+4xi,l3Kxi,l,xi,l2fxi,l2,∑q=03μqiclq+∑j=1mβjiclj+3+xi,l3Kxi,l,xi,lfxi,l,∑q=03μqiclq+∑j=1mβjiclj+3,l=1,2,…,N,M(i)=μ0(i),…,μ3(i)t, B(i)=β1(i),…,βm(i)t, and R(i)=-h(xi,1),…,-h(xi,N)t. Here ·t means the transpose of the vector. Based on the given definition of the approximate spline u, it can be easily verified that (30)M0=μ00μ10μ20μ30=a0a1a2a3.


## 4.2. Multiple Shooting Method

For simplicity, let us rewrite problem (1)–(2b) in the following form:(31)Ly≔Dαyx+∫0xKx,tft,ydt=-hx,subject to(32)y0=a0,y′0=a1,y1=b0,y′1=b1.The solution of problem (31)-(32) can be determined by solving the problem on the subintervals xi,xi+1, for i=0,1,…,N-1. In order to apply the shooting method, we introduce the following set of initial value problems:(33)Lyi=Dαyix+∫0xKx,tft,yidt=-hx,xi≤x≤xi+1subject to(34)yixi=a4i,yi′xi=a4i+1,yi′′xi=a4i+2,yi′′′xi=a4i+3,where a2,a3,…,a4N-1 are unknown real parameters. The solution of problem (33)-(34) will be obtained by the method described in the previous subsection, where the parameters a2,a3,…,a4N-1 are determined by solving the following algebraic system (35)y0nx1,a2:a3=y1nx1,a4:a7,⋮yi-1nxi,a4i-4:a4i-1=yinxi,a4i:a4i+3,⋮yN-11,a4N-4:a4N-1=b0,yN-1′1,a4N-4:a4N-1=b1,where n=0:3 and al:am represents al,al+1,…,am. It is worth mentioning that we use shooting method of order five (N=5) in our computations which are done using Matlab.


## 5. Numerical Results

Example 1.
Consider the following nonlinear fourth order fractional integrodifferential equations:(36)Dαyx+∫0xt+1y2tdt+hx=0,0<x<1subject to(37)y0=1,y′0=0,y1=2,y′1=3,where α=3.7, h(x)=-x8/8-x7/7-2x5/5-x4/2-x2/2-x, and the exact solution is y(x)=x3+1. Divide the domain [0,1] into five subinterval such that(38)⋃i=04xi,xi-1wherexi=0.2i,i=0,1,…,5,and then apply the multiple shooting method of order five yields to solve the following five initial value problems:(39)Ly0=-hx,0≤x≤0.2,y00,:=1,y0′0,:=3,y0′′0,:=a2,y0′′′0,:=a3,Ly1=-hx,0.2≤x≤0.4,y10.2,:=a4,y1′0.2,:=a5,y1′′0.2,:=a6,y1′′′0.2,:=a7,Ly2=-hx,0.4≤x≤0.6,y20.4,:=a8,y2′0.4,:=a9,y2′′0.4,:=a10,y2′′′0.4,:=a11,Ly3=-hx,0.6≤x≤0.8,y30.6,:=a12,y3′0.6,:=a13,y3′′0.6,:=a14,y4′′′0.6,:=a15,Ly4=-hx,0.8≤x≤1.0,y40.8,:=a16,y4′0.8,:=a17,y4′′0.8,:=a18,y4′′′0.8,:=a19,where(40)Ly=Dαyx+∫0xt+1y2tdt.Here y0x,:=y0x,a2,a3, y1x,:=y1x,a4:a7, y2x,:=y2x,a8:a11, y3x,:=y3x,a12:a15, and y4x,:=y4x,a16:a19. The above initial value problems (39) are solved using the collocation method with the following collocation points: (41)c1=1+cos⁡π/82,c2=1+cos⁡3π/82,c3=1+cos⁡5π/82,c4=1+cos⁡7π/82.Notice that these collocations points are generated from the roots of the Chebyshev polynomials of degree four. To find the parameters ai, i=2:19 we solve the following algebraic system:(42)y0n0.2,a2:a3=y1n0.2,a4:a7,y1n0.4,a4:a7=y2n0.4,a8:a11,y2n0.6,a8:a11=y3n0.6,a12:a15,y3n0.8,a12:a15=y4n0.8,a16:a19,y41,a16:a19=2,y4′1,a16:a19=3,for n=0,1,2,3, using Matlab software. Hence, we obtain (43)a2=5.999999998,a3=6.000000001,a4=1.008000002,a5=0.120000003,a6=1.199999998,a7=6.000000001,a8=1.063999996,a9=0.48000005,a10=2.40000002,a11=6.000000002,a12=1.21600007,a13=1.08000008,a14=3.599999995,a15=5.999999994,a16=1.51199989,a17=1.919999990,a18=4.800000012,a19=5.999999986.The graph of the exact and approximate solutions and the graph of the error function are, respectively, given in Figures 1 and 2. It is clearly seen that the two solutions are in excellent agreement. In addition, the computed L2 error norm is given by (44)ux-yx=∑i=04∫xixi+1ux-yi+1x2dx=2.1×10-14.Figure 1
The exact and approximate solutions for Example1.Figure 2
Computed absolute error between the present numerical solution and the exact one for Example1.Example 2.
Consider the following nonlinear fourth order fractional integrodifferential equations:(45)Dαyx+∫0xt2ytdt-x3α+3+1Eα,1xα=0,0<x<1subject to(46)y0=1,y′0=1αΓα,y1=Eα,11,y′1=Eα,α1α,where α=3.5 and the exact solution is y(x)=Eα,1(xα). Here Eα,β(z)=∑n=0∞zn/Γ(αn+β) represents the Mittag-Leffler function. Using the multiple shooting method of order five and then following the same steps in the previous example, the parameters ai, i=2,3,…,19 are found to be (47)a2=0.000396825,a3=0.0000000504,a4=1.0172022865,a5=0.086051128,a6=0.0003969266,a7=0.0000000503,a8=1.034420446,a9=1.49160490,a10=0.000397028,a11=0.0000000504,a12=1.051654498,a13=0.086209945,a14=0.000397123,a15=0.0000000505,a16=1.068904427,a17=0.086289345,a18=0.000397224,a19=0.0000000504.The exact solution is graphed along with the graph of the approximate solution for α=3.5 in Figure 3. Moreover, the computed absolute error between the numerical and exact solutions is given in Figure 4 which shows a satisfactory agreement. The computed L2 error norm is given by (48)ux-yx=∑i=04∫xixi+1ux-yi+1x2dx=1.87×10-15.Figure 3
The exact and approximate solutions for Example2.Figure 4
Computed absolute error between the present numerical solution and the exact one for Example2.Example 3.
Consider the following nonlinear fourth order fractional integrodifferential equations(49)Dαyx+∫0xy2tdt+Hx=0,0<x<1subject to(50)y0=0,y′0=1,y1=0.84,y′1=0.54,where α=3.25 and H(x)=0.005-0.13x-x2. By implementing the present algorithm, the parameters ai, i=2,3,…,19 are found to be (51)a2=0,a3=-1,a4=0.198669,a5=0.980067,a6=-0.198669,a7=-0.980067,a8=0.389418,a9=0.921061,a10=-0.389418,a11=-0.921061,a12=0.564642,a13=0.825336,a14=-0.564642a15=-0.825336,a16=0.717356,a17=0.696707a18=-0.717356,a19=-0.696707.Figure 5 shows the graph of the approximate solution. Since the actual solution for this problem is unknown, we may measure the error of the approximation by using the residual function, R(x), defined by (52)Rx≔Dαux+∫0xu2tdt+Hx.The graph of the residual function, R(x), is displayed in Figure 6. It is clearly evident from Figure 6 that MaxRx:0≤x≤1=0.0056.Figure 5
The exact and approximate solutions for Example3.Figure 6
Computed residual,R(x), for Example 3.


## 6. Conclusion and Future Work

In this paper, we have solved special class of higher-order nonlinear integrodifferential equations of order 4 subject to boundary conditions. The method of solution is based on conjugate collocation and spline technique with shooting method. The numerical results for given examples demonstrate the efficiency and accuracy of the present method. It should be noted that applying Theorem5 causes two difficulties. Firstly, the Lipschitz constant for f(x,y) with respect to y is not easy to find and, secondly, the necessary condition (0<LM/αΓ(α)<1) does not satisfy several cases of the problem (1)–(2b). Therefore, a further discussion on weaker necessary conditions should be followed in the future work to cover wider range of problem (1)–(2b).



---

*Source: 616438-2015-11-02.xml*