Journal of Applied Mathematics

Volume 2013 (2013), Article ID 864025, 20 pages

http://dx.doi.org/10.1155/2013/864025

## Cubic Spline Collocation Method for Fractional Differential Equations

^{1}Department of Mathematics, Huizhou University, Guangdong 516007, China^{2}School of Mathematics and Computational Science, Hunan Key Laboratory for Computation and Simulation in Science and Engineering, Xiangtan University, Hunan 411105, China

Received 11 February 2013; Revised 17 April 2013; Accepted 7 May 2013

Academic Editor: Alain Miranville

Copyright © 2013 Shui-Ping Yang and Ai-Guo Xiao. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

#### Abstract

We discuss the cubic spline collocation method with two parameters for solving the initial value problems (IVPs) of fractional differential equations (FDEs). Some results of the local truncation error, the convergence, and the stability of this method for IVPs of FDEs are obtained. Some numerical examples verify our theoretical results.

#### 1. 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–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 higher-order algorithm based on the Simpson method for FDEs, and gave the error and stability results. Li [21] applied the cubic B-spline 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 The relationship between the Caputo definition and the Riemann-Liouville definition can be given by the following Lemma.

Lemma 1 (see [5, 6]). *Let and . Assume that is times continuously differentiable, and exist; then
*

#### 2. Cubic Spline Collocation Method for FDEs

Consider the initial value problem where is a given continuous mapping and satisfies the Lipschitz condition It is well known that the IVP (5) is equivalent to the problem [24]

Let be the grid points of the uniform partition of into the subintervals . On each subinterval , the cubic -spline function can be represented as where and , .

Let , . We can write (9) as where . Obviously, , and , for all .

*Remark 2. *In this algorithm, the initial value must be provided, but it does not need to be given in problem. That is, if we approximate , then the approximation of is given by using the low-order methods. In this paper, we obtained approximation of 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 are the collocation points, and , is a constant. And according to (10), we have 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

Obviously, it follows from the inequalities 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 addition, applying the collocation conditions in each subinterval and using (7), we have where are the collocation points, , and where ,,,, .

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).

#### 3. Theoretical Results

Lemma 3 (see [25]). *Suppose that , and satisfy
**
where is independent of . Then
**
where is a Mittag-Leffler function.*

Lemma 4. * Assume that the function . Then, for all , one has
**
where , is a constant and satisfies or , , is infinitesimal of higher order.*

*Proof. * Using the Taylor formula, we have
Applying the integral mean value theorem, we can obtain
where . Since , is bounded on . We have

Consider the following. (i) If , then
(ii) If , then

Thus,
By means of (30), Lemma 4 is established obviously.

Theorem 5 (the local truncation error). * If the analytical solution of the problem (5) , then the local truncation error of DSCM is . *

*Proof. * Applying the Taylor expansion and the definition of , we have
From the definition of the Caputo fractional derivative, we obtain
That is,
where . Moreover,

Supposing , for all , 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.

To prove the convergence of DSCM for the problem (5), we first give the following lemmas.

Lemma 6. *If the analytical solution of the problem (5) and the matrixes and in (14) satisfy that exists, and
**
then the numerical solutions of ISCM for problem (5) satisfy that . *

*Proof. *For (23), we have
It follows from (7) that the analytical solution of the problem (5) satisfies

Moreover, (44) and (45) yield

By using the Lipschitz condition, we have

It follows from the definition of that
where , . Obviously, for (14) and (46), the matrixes and are continuous on ; thus, we have
Therefore, there exist , such that

Obviously,
Applying the integral mean theorem, we have
It follows from the inequalities and that

Combining (45), (48), and (49) with (51) and defining
where , , we can obtain
where are some positive constants. If the matrix is invertible, then it follows from (47) that we have
Hence,
where is an appropriate positive constant. Take , . Then , , and

Since , we have
where , is an appropriate positive constant. According to (53), we obtain
where is an appropriate constant. Obviously, there exist and such that for . Thus,
Let ; then
where is an appropriate constant. Applying Lemma 3, we get
By using the convergence of the Mittag-Leffler function [6], we have

By means of (24), (47), we obtain
where is a positive constant. Consequently, we have

*Remark 7. * If , then
is invertible, and

Consider

When , the range of values can be obtained by using Algorithm 8, which can be shown in Figures 1 and 2, respectively.

*Algorithm 8. * For For * *if record ; end end end plot.

Lemma 9. *If the analytical solution of the problem (5) , the spline functions and defined as (9) and (24) on the same grids are the numerical solutions of the problem (5) obtained by DSCM and ISCM, respectively; then we have the following conclusions.*(i)*If , then , where is the analytical solution of the problem (5).*(ii)*Assume that and are the numerical solutions obtained by DSCM and ISCM for the problem
**and the perturbed problem
**
where satisfies the Lipschitz condition (6), respectively. Then there exists the constant such that, for all , when , there exists the constant such that ; here is a constant, .*

*Proof. *Firstly, we give the proof of the conclusion (i). Applying DSCM, we have
Since , we obtain
It follows from Lemma 4 that
where as .

By means of (70) and (72), we obtain
Applying Lemma 4, we get