Pade Method for Construction of Numerical Algorithms for Fractional Initial Value Problem

In this paper, we propose an eﬃcient method for constructing numerical algorithms for solving the fractional initial value problem by using the Pade approximation of fractional derivative operators. We regard the Grunwald–Letnikov fractional derivative as a kind of Taylor series and get the approximation equation of the Taylor series by Pade approximation. Based on the approximation equation, we construct the corresponding numerical algorithms for the fractional initial value problem. Finally, we use some examples to illustrate the applicability and eﬃciency of the proposed technique.

Among all the numerical methods, the direct numerical method [16][17][18] is the basic one. Due to the special property of the fractional derivative, most of these methods have their inbuilt deficiencies mostly due to the calculation of Adomian polynomials, the Lagrange multiplier, divergent results, and some other huge computational works.
In this paper, we propose an efficient method for constructing numerical algorithms for solving the fractional initial value problem by using the Pade approximation of fractional derivative operators. e advantage of the proposed technique is that efficient numerical methods can be constructed without calculation of long historical terms of the fractional derivative.
Consider the following homogeneous fractional initial problem: where 0 D (α) is usually the Caputo derivative and m � [α]. Given the condition y (k) (0) � 0, k � 0, 1, . . . , m − 1, the Caputo derivative is equivalent to the Grunwald-Letnikov derivative and the Riemann-Liouville derivative, e.g., For the following linear nonhomogeneous fractional initial problem, the following transform can be employed to make it be a homogeneous one with regard to z(t): Based on the G algorithm [19], we have the following numerical scheme for problem (1): nj y j � f t n , y n , n � 0, 1, . . . , where c (α) nj and N are calculated according to different algorithms. Algorithm 5 is a technique which is easy to manipulate. However, this method involves huge computational work when T ≫ 0 because in general, these methods need to compute many terms to get the approximation to the fractional derivative. To some extent, the short memory principle [20] can be used to tackle the problem, but low accuracy will be the cost. So, it is significant to find an efficient approximant to the fractional derivative, which is of low computation cost on the one hand and highly accurate on the other hand. In this paper, we put forward a reliable method for the construction of numerical algorithms for solving the fractional differential initial value problem by using Pade approximation to fractional derivative operators. e rest of the paper is organized as follows. In Section 2, we briefly list some basics of Pade approximation. In Section 3, we get the approximant of the fractional derivative operators. In Section 4, we propose a method for constructing algorithms for solving the fractional initial value problem by using Pade approximation. In Section 5, we use some examples to illustrate the applicability and efficiency of the proposed technique. Section 6 is the conclusion.

Some Basics for Pade Approximation
In numerical mathematics, Pade approximation [21] is believed to be the best approximation of a function by rational functions of a given order. Under this technique, the approximant's power series agrees with the power series of the function it is approximating. e Pade approximant often gives better approximation of the function than truncating its Taylor series, and it may still work where the Taylor series does not converge. For these reasons, Pade approximants are often used in many fields of computations.
Given a function f(x) and two integers m > 0 and n > 1, the Pade approximant of order [m, n] is the rational function which agrees with f(x) to the highest possible order, which amounts to Equivalently, if R(x) is expanded in a Maclaurin series, its first m + n terms would cancel the first m + n terms of f(x), and as such e Pade approximant is unique for given m and n, that is, the coefficients a 0 , a 1 , . . . , a m , b 1 , . . . , b n can be uniquely determined. It is for reasons of uniqueness that the zeroth-order term at the denominator of R(x) was chosen to be 1; otherwise, the numerator and denominator of R(x) would have been unique only up to multiplication by a constant. e Pade approximant defined above is also denoted as [m, n] f(x) .

Pade Approximant for the Fractional Derivative Operator
ere are many definitions for the fractional derivative. e following equation is called the reverse Grunwald-Letnikov derivative: For any given α > 0, if we denote we have erefore, we can denote So, we have Let a 0 + a 1 x + a 2 x 2 + · · · + a n x n be the [n, n] Pade approximant to (1 − x) α ; we have a 0 + a 1 x + a 2 x 2 + · · · + a n x n erefore, a 0 + a 1 x + a 2 x 2 + · · · + a n x n where 1 can be understood as the identical operator. So, we get the following approximation: Noticing (13), we get So, we get en, we get the following approximation to fractional derivative operators: where h is the discretization step size. From (21), we have Due to the linearity and commutativity of the operators E − h and D α , we have a 0 y n + a 1 y n− 1 + a 2 y n− 2 + · · · + a n y 0 � h α D α y n + b 1 y n− 1 + b 2 y n− 2 + · · · + b n y 0 .
Noticing D α y n � f n , we get a 0 y n + a 1 y n− 1 + a 2 y n− 2 + · · · + a n y 0 en, we get the following numerical algorithm: a 0 y n � − a 1 y n− 1 − a 2 y n− 2 − . . . − a n y 0 is is a multiple step algorithm for problem (1), where a i s and b i s can be determined by the Pade approximation of (1 − x) α for any given α. From the construction process, we can see the local truncation error of this algorithm is 2n + 1.

Some Implicit Multistep Algorithms for the Fractional Initial Value Problem
Let α � 0.5; the [2, 2] Pade approximant to (1 − x) 0.5 is and the [3,3] Pade approximant to (1 − x) 0.5 is en, from (28), we get the following two multistep algorithms: in which 64y n � 112y n− 1 − 56y n− 2 + 7y n− 3 in which In the same way, from the following (36) So, we get the following numerical algorithm: We can also get the following algorithm by using the [3,3] Pade approximant to (1 − x) 0.8 : and many other more complicated numerical schemes can be constructed by using the Pade approximant for any given α.

Numerical Test 1.
e exact solution to this problem is y � sin(x). We use algorithms (1) and (2) to solve this problem. Because both algorithms (1) and (2) are implicit, we use the iterative method y n � g y n , to get the value of y n , where g y n � 5 4 y n− 1 − 5 16 for algorithm 1 and for algorithm 2, respectively. e computational errors are listed in Table 1 (h � (pi/2)/10). We also compare algorithms (1) and (2) with the exact solution in Figure 1.
We also make the comparison with the corresponding G algorithm. e results are listed in Table 2. e comparison shows that algorithms (1) and (2) are more efficient than the corresponding Grunwald-Letnikov-based G algorithm.

Numerical Test 2.
D 0.8 y � y cos 2π 5 + ����� 1 − y 2 sin 2π 5 , e exact solution to this problem is y � sin x. We use the following algorithm (4) and algorithm (5 to get the numerical solution. Because both algorithms (4) and (5) are implicit schemes, we use the iterative method y n � g(y n ) to get the value of y n , where g y n � 7 5 y n− 1 − 21 50 in algorithm (4) in algorithm (5). e errors of the numerical results are listed in Table 3 (h � (pi/2)/10). e comparison with the corresponding G algorithm is listed in Table 4. e comparison shows that algorithms (4) and (5) are more efficient than the corresponding G algorithm.

Numerical Test 3.
e exact solution is y � x 4 ; we use algorithm 5 to solve this problem. e comparison between the exact and the numerical solution can be seen in Figure 2. One can see that the numerical solution is very accurate when T ≫ 1, and the result is obtained in less than 50 steps. is is much more efficient than the G algorithm.

Conclusion
We can construct numerical algorithms for solving the fractional initial value problem based on the Pade approximation of fractional derivative operators. Numerical tests show that this method is more efficient than the corresponding G algorithms. Generally speaking, we can calculate more terms to achieve more accuracy when Grunwald-Letnikov-based method is employed, but that would lead to more computational work, while the algorithms derived from Pade approximation are proved to be more efficient.

Data Availability
e data used to support the findings of this study are available from the corresponding author upon request.

Conflicts of Interest
e authors declare that they have no conflicts of interest.