The spectral collocation method for solving a fractional integro-di ﬀ erential equation

: In this paper, we propose a high-precision numerical algorithm for a fractional integro-di ﬀ erential equation based on the shifted Legendre polynomials and the idea of Gauss-Legendre quadrature rule and spectral collocation method. The error analysis of this method is also given in detail. Some numerical examples are give to illustrate the exponential convergence of our method.


Introduction
All kinds of natural phenomena in nature can be explained by fractional calculus. For example, fractional integral and derivative theory are used to solve problems in nonlinear vibration of viscoelastic damping, control signal and processing, thermoelectric viscoelastic fluid, continuum and statistical mechanics, and other fields [1][2][3][4][5]. It is well known that generalized fractional operators are obtained by the extension of Caputo and Riemann-Liouville fractional derivatives. In convolution, because of the complexity of weight and kernel function, it is more difficult to design the higher-order numerical scheme of calculus equations with generalized fractional operators.
In the past decades, many numerical methods have been proposed for the differential and integral equations of different type problems. Such as, collocation method based on the Muntz-Legendre polynomials [6], Wavelet-Galerkin method and homotopy perturbation method [7], Legendre wavelets method [8], the rational second kind Chebyshev pseudospectral method [9], An hp-version spectral collocation method [10], and so on [11][12][13][14]. As one of the efficient numerical computing technique, the spectral collocation method provides a powerful tool for solving fractional calculus equations in recent years. Therefore, we also hope to provide a useful numerical method for fractional calculus equations.
In this paper, we considered the following fractional integro-differential equation: with an initial condition where 0 < α < 1, 0 < µ < 1, t ∈ I = [0, 1], D α t is the classical fractional derivative of order α, f (t) and K(s, t) are known functions, and y(t) is unknown function. This type of fractional integral-differential equation appear in several other places as in the modelling of particle motion in physics and mechanics in [15].
After careful planning, this article is developed as follows: the first section introduces the research background of fractional calculus equation model.
Then, some results of fractional integro-differential and the properties of the polynomials in Section 2. In Section 3, based on the polynomials, a new method of the spectral collocation and Gauss-Legendre quadrature rule is proposed for solving fractional integro-differential problems. In Section 4, we analyze the error truncation of the numerical method. It considers some examples to illustrate the high accuracy of the proposed approach in Section 5. Finally, the conclusion is given.

Fractional calculus
Several definitions of fractional calculus has been developed over time, it can be seen in [16,17]. Definition 1. [16] The Riemann-Liouville fractional integral operator of order α ≥ 0 is defined Definition 2. [16] The fractional derivative of order α in the Caputo sense is defined The properties of the operator I α and D α as follows x γ+α , x γ−α , otherwise.

The shifted Legendre polynomials
According to the definition in [16,18], we make it easier to solve Eq (1.1), we first introduce its definition about the shifted Legendre polynomials of degree i over the interval [0, 1], as follows and it is orthogonal for the legendre polynomial L i on [−1, 1].
The shifted Legendre polynomial ϕ i (t) satisies the following form and the complete basis is formed on [0, 1]. Then, the analytic form of ϕ i (t) may be expressed with t as

Function approximation
In reference [18], the following symbols are defined in L 2 (I) for weighted inner product and its norm as Property 1. [18] Suppose u ∈ L 2 (I) is arbitrary, there exists a unique q * m ∈ P m , the following relation holds forms an L 2 -orthogonal basis for P m .
By Property 1, for any y ∈ L 2 (I), there is a unique best approximation y N , such that The relationship between the function y(t) ∈ L 2 (I) and the shifted Legendre polynomials as follows Then, the shifted Legendre polynomials with first (N+1)-terms are considered as Similarly, for any function K(s, t) ∈ L 2 (I) with the variable s, it can be expanded by ϕ i (s) as which is used the Gauss-Legendre quadrature rule to compute the coefficients in this paper.

Fractional operational matrices of derivative and integral
An approximated of the derivative of the shifted Legendre polynomials is given as where D γ is the derivative operation matrix of the shifted Legendre polynomials and Φ(s) = (ϕ 0 (s), . . . , ϕ N (s)) T .

Lemma 1. [19]
Suppose D γ is an operation matrix of fractional derivatives of order γ > 0 for the generalized shifted Legendre functions of m × m in the Caputo sense, then the d (γ) i j are given for D γ as so called P υ is integral operator matrice of the shifted Legendre functions. Lemma 2. [16] Suppose P υ is the m × m the generalized shifted Legendre functions operation matrix of Riemann-Liouville fractional integral of order υ , then P (υ) i, j can be obtained for P υ as

Numerical algorithm
In this part, based on the shifted Legendre polynomials, we consider the idea of combination of Gauss-Legendre quadrature rule and spectral collocation method to solve the Eq (1.1). Then, using this numerical algorithm, we first expand y(t) and K(s, t) by the shifted Legendre polynomials of first (N+1)-terms as where the coefficients are obtained by (2.1) and (2.2), respectively. Substituting (3.1) and (3.2) in Eq (1.1), we obtain based on the orthogonality of polynomials on certain interval, have where To solving the coefficients Z, we substitute collocation points into Eq (3.3) lead to where the points of shifted Chebyshev-Gauss are t r = (1−cos((2r +1)π/(2N +2)))/2, r = 0, . . . , N −1. Thus, a system is formed by the N nonlinear equation. In addition, an algebraic equation with the initial condition is provided to obtain Z with y N (0) = y(0).
Then, by Eqs

Error analysis
The errors of the proposed method are analyzed, which is the work of this section. Then we introduce the following approximation operator, we can see in [20]. In the time direction, we refer to Π N as the traditional orthogonal projection operator of L 2 (I) and P N (I) is a space of polynomials of degree up to N with the time variable t. Π N : L 2 (I) → P N (I), that is, for any y ∈ L 2 (I), Π N y ∈ P N (I), it satisfies ((Π N y − y), φ) I = 0, ∀φ ∈ P N (I).
In the following lemma, we give an error estimate for approximating the operator Π N .
where C is a constant that is different from N. where C is a constant independent of N.
Similar to the Lemma 3.6 in [21], it is easy to proof the following Lemma.
Next, replace y(t) and K(s, t) with y N (t) and K N (s, t) in (1.1), rewriting the equation as follows we define the following operators In the following section, we estimate the error for L(y(t) − y N (t)).
Theorem 1. Let y and y N be the analytic solution and numerical solutions of (1.1) and (1.2). Also assume that y γ < ∞, γ > 1, ||||K(s, t)|| m || 0 ≤ K, m ≥ 0, then where K is a real numbers and C is a constant independent of N.
Proof. By Eqs (1.1) and (4.1) and triangle inequality, we obtain where R 1 , R 2 , R 3 are defined by the right term of (4.3), respectively.
For R 1 , it is easy to obtain based on the Lemma 3 as For R 2 , based on the Lemma 5 it is directly to obtain that where 0 < σ < min(1/2, µ).

Numerical examples
In this section, it is verified about the validity and applicability of the method by the following examples.
Example 2. The fractional integro-differential equation is considered as its exact solution is y(t) = e t − t − 1.
In Eqs (5.1) and (5.2) with the initial condition y(0) = 0. In the numerical implementation, we provide a value of µ and the nodes and weights of the Gauss-Legendre quadrature rule. With the increasing value of N, the error result of different α are shown in Figure 1. A logarithmic scale has been used for error-axis in figures. Clearly, the errors show an exponential decay, since in these semilog representations one observes that the error variations are essentially linear versus the polynomial degrees for α. The results show that the method with a high accuracy to solving fractional integrodifferential equations, and the error of y(t) − y N (t) 0 converge exponentially with the increase of N and reach machine precision. Meanwhile, the L 2 error has the advantage that it can better reflect the actual situation error of y(t) − y N (t) 0 in Tables 1 and 2.  Table 2. Error of y(t) − y N (t) 0 for Example 2.
In order to reflect the effectiveness of the proposed method, the results be compared with euler method in following example.
with the initial condition y(0) = 0. Its exact solution is y(t) = t 3 .
From Table 3, it is clear that the error of absolute with the present method by polynomials of degree N approximation. Clearly, in Example 3 the method gives more accurate results by using lower order polynomials and fewer points than the euler method, and our method is easier to implement and understand.  Table 3. e present (t) = |y(t) − y N (t)|, e euler (t) = |y(t) −ỹ(t)|,ỹ(t) is the computed solution by the euler method, h is the step size.

Conclusions
Based on the shifted Legendre polynomials, the idea of combining Gauss-Legendre quadrature rule and spectral collocation method, a new spectral collocation method is proposed for solving a fractional integro-differential equation in this paper. The error of this method is analyzed by Theorem 1. It is found that the error decreases exponentially from Figure 1, Tables 1 and 2. In Table 3, the method gives more accurate results than the euler method. The results show that this method is high accuracy and easy to be implemented. In the future, based on the idea of [22,23], we will use the shifted Legendre polynomials for stochastic fractional integro-differential problems and shape optimization.