Cubic Spline Collocation Method for Fractional Differential Equations

In the recent past years, the use of fractional differential equations (FDEs) has gained considerable popularity in some fields (e.g., nonlinear oscillation of earthquake [1], fluid-dynamic traffic model [2], material viscoelastic theory and physics [3–6], etc.). Based on these requirements, the numerical approach has become very important to solve FDEs and analyze the experimental data which are described in a fractional way. The numerical approaches to FDEs have been recently studied by numerous authors [6–15]. However, the state of art is far less advanced for general fractional order differential equations. In the recent years, spline collocation methods including wavelet methods have been successfully applied to some initial value problems (IVPs) and initial-boundary value problems of FDEs [16–23]. Pedas and Tamme [16, 18] discussed spline collocation methods for some classes of IVPs of linear multiterm fractional differential equations and obtained the corresponding convergence results. Li et al. [19] used the higher-order piecewise interpolation for the fractional integral and fractional derivatives, proposed a higherorder algorithm based on the Simpsonmethod for FDEs, and gave the error and stability results. Li [21] applied the cubic Bspline wavelet collocation method to FDEs by the approach that such problems are converted into a system of algebraic equations which is suitable for computer programming. It is worth notice that the spline collocation methods for the IVPs of FDEs are often achieved by their solving of the equivalent IVPs of integral equations with weakly singular kernels. In this paper, the cubic spline collocation method is designed to solve directly the IVPs of general FDEs. And the corresponding theoretical results of the local truncation error, the convergence, and the stability of the cubic spline collocation method for the IVPs of general FDEs are given. This paper is organized as follows. In Section 2, we propose the spline cubic collocation method for solving the IVPs of FDEs. In Section 3, the corresponding theoretical results of the convergence and the stability are also given. In Section 4, the theoretical results are identified by some numerical examples. In the following text, we first recall the basic definitions of fractional calculus [6]. Usually, C a D α


Introduction
In the recent past years, the use of fractional differential equations (FDEs) has gained considerable popularity in some fields (e.g., nonlinear oscillation of earthquake [1], fluid-dynamic traffic model [2], material viscoelastic theory and physics [3][4][5][6], etc.).Based on these requirements, the numerical approach has become very important to solve FDEs and analyze the experimental data which are described in a fractional way.
The numerical approaches to FDEs have been recently studied by numerous authors [6][7][8][9][10][11][12][13][14][15].However, the state of art is far less advanced for general fractional order differential equations.In the recent years, spline collocation methods including wavelet methods have been successfully applied to some initial value problems (IVPs) and initial-boundary value problems of FDEs [16][17][18][19][20][21][22][23].Pedas and Tamme [16,18] discussed spline collocation methods for some classes of IVPs of linear multiterm fractional differential equations and obtained the corresponding convergence results.Li et al. [19] used the higher-order piecewise interpolation for the fractional integral and fractional derivatives, proposed a higherorder algorithm based on the Simpson method for FDEs, and gave the error and stability results.Li [21] applied the cubic Bspline wavelet collocation method to FDEs by the approach that such problems are converted into a system of algebraic equations which is suitable for computer programming.
It is worth notice that the spline collocation methods for the IVPs of FDEs are often achieved by their solving of the equivalent IVPs of integral equations with weakly singular kernels.In this paper, the cubic spline collocation method is designed to solve directly the IVPs of general FDEs.And the corresponding theoretical results of the local truncation error, the convergence, and the stability of the cubic spline collocation method for the IVPs of general FDEs are given.
This paper is organized as follows.In Section 2, we propose the spline cubic collocation method for solving the IVPs of FDEs.In Section 3, the corresponding theoretical results of the convergence and the stability are also given.In Section 4, the theoretical results are identified by some numerical examples.
In the following text, we first recall the basic definitions of fractional calculus [6].Usually,      denotes the Caputo fractional derivative of order  as where   is the classical differential operator of order , () is  times continuously differentiable, and   denotes the integral operator of order  as where  is the identity operator.
As is well known, there are some different definitions of fractional operator except the Caputo fractional derivative.From a theoretical point of view, the most natural approach is the Riemann-Liouville definition as (3) The relationship between the Caputo definition and the Riemann-Liouville definition can be given by the following Lemma.
Remark 2. In this algorithm, the initial value   (0) must be provided, but it does not need to be given in problem.That is, if we approximate   (0), then the approximation of   (0) is given by using the low-order methods.In this paper, we obtained approximation of   (0) by using one-order BDF method.That is, Then where ℎ is a tiny stepsize.
By using the collocation conditions in each subinterval   and (5), we have where  1 ,  2 are the collocation points, and 0 <  0 ≤  1 <  2 ≤ 1,  0 is a constant.And according to (10), we have −1 where Using the definition of the Caputo fractional derivative in the form (1) and using (13), we have Through calculation, we can get where In order to obtain the corresponding iteration formula, we denote ) , ))  (1)   −1 ))  (1)   −1 Obviously, it follows from the inequalities  0 ≤  1 <  2 ≤ 1 that the matrix Ã is invertible, and where Let Then, it follows that For (22), we can obtain their numerical solutions by using the Newton iterative method.
In this paper, for convenience, the cubic spline collocation method based on the direct discretization with ( 13) is called direct spline collocation method (DSCM for short); and the cubic spline collocation method based on the fractional order integral equations (7) and ( 27) is called indirect spline collocation method (ISCM for short).
Proof.Using the Taylor formula, we have Applying the integral mean value theorem, we can obtain (ii) If  0 ℎ ≤  − ( + ℎ), then Thus, By means of (30), Lemma 4 is established obviously.Proof.Applying the Taylor expansion and the definition of (, ℎ), we have From the definition of the Caputo fractional derivative, we obtain That is, where  = 1, 2.Moreover, Supposing (, ℎ) = (), for all  ∈ [0, −1 ], and substituting () into ( 13), we have where   is the local truncation error.
Combining (38), (37) with the Lipschitz condition yields where Ĉ is an appropriate positive constant.And the norm ‖ ⋅ ‖ is the 1-norm of the matrix (⋅) in this paper.
Proof.This result follows directly from Lemmas 6 and 9.
Now, we consider the stability of DSCM. , } are the numerical solutions of the problem (68), and { } are the numerical solutions of the perturbed problem (69), and where   > 0 is a constant, (, ℎ) is the numerical solution of the problem (68), and (, ℎ) is the numerical solution of the problem (69).
Proof.According to Lemma 9, if ISCM is stable, then DSCM is stable.In the following text, we give the proof of stability of ISCM.

Illustrative Examples
In order to demonstrate our theoretical results, we apply DSCM to the problem (5) and present some numerical examples in this section.Let The exact solution is () =  2 .
For different values of  ∈ (0, 1), the numerical solutions for problem (107) are obtained by using DSCM, fast wavelet collocation method (FWCM for short) reported in [21], and the method reported in [10].When  = 0.5, the time stepsize ℎ = 1/160; the absolute errors of DSCM, FWCM, and the method reported in [10] are shown in Table 1.When  = 0.75, ℎ = 1/160; the absolute errors of these methods are shown in Table 2.The numerical solutions and the exact solution are shown in Figure 3. From the numerical results, the results obtained by DSCM are better than by FWCM and the method reported in [10] in terms of accuracy if the exact solution is sufficiently smooth.Therefore, DSCM is a valid method in solving fractional differential equation.
Tables 3 and 4 list the errors and the error orders of DSCM with  = 0.5 and  = 0.75 and show that the errors between the numerical solutions and the exact solution are very small, respectively.Figures 4 and 5 illustrate high accuracy of DSCM with  1 = 0.9,  2 = 0.98 and show that the error is also very small.All of the numerical results show that DSCM for solving nonlinear FDEs is convergent and the method is robust.We use the cubic spline collocation method to solve the problems (110) and (111), respectively.Selecting () =  sin(),  = 20,  1 = 0.1,  2 = 0.2, we obtain the numerical results given in Figure 6.When () =  sin(),  = 20,  1 = 2.0,  2 = −2.0, the numerical results are given in Figure 7.When () =   ,  = 20,  1 = 1.0,  2 = 2.0, the numerical results are shown in Figure 8.

Example 3. Consider the initial value problem
From these figures, we can see that the absolute errors of the numerical solutions of the problems (110) and (111) decrease and finally tend to 0 as  increases.Thus, we can draw the conclusion that the cubic spline collocation method for nonlinear FDEs is stable.The numerical results verify our theoretical results.

Conclusion
In this paper, the cubic spline collocation method with two parameters is successfully applied to the IVPs of general FDEs.The result of the local truncation error of this method is given.And the convergence and stability results of the cubic  spline collocation method for the fractional order integral equations which is equivalent to the IVPs of general FDEs are obtained.By using the relationship between the numerical solutions from the cubic spline method for the IVPs of general FDEs and the numerical solutions obtained from the cubic spline method for the corresponding equivalent IVPs of fractional order integral equations, we also obtain some results of the convergence and the stability of the method for the IVPs of general FDEs.Some numerical examples successfully verify our theoretical results and show that the given method is efficient.

Figure 4 :
Figure 4: The error between the numerical solutions and the exact solution of the problem (109); ℎ = 0.05,  = 0.000001.