Application of Bessel functions and Jacobian free Newton method to solve time-fractional Burger equation

Abstract In this paper, a novel method based on Bessel functions (BF), generalized Bessel functions (GBF), the collocation method and the Jacobian free Newton-Krylov sub-space (JFNK) will be introduced to solve the nonlinear time-fractional Burger equation. In this paper, an implicit formula is introduced to calculate Riemann–Liouville fractional derivative of GBFs, that can be very useful in spectral methods. In this work, the nonlinear time-fractional Burger equation is converted to a nonlinear system of algebraic equations via collocation algorithm based on BFs and GBFs without any linearization and descretization methods. Finally, by using JFNK, the solution of this nonlinear system will be achieved. To show the reliability and applicability of the proposed method, we solve some examples of time-fractional Burger equation and compare our results with others.

The spectral methods are one of the high accuracy numerical tools in computational and numerical algorithms [10,11,20,21]. But, spectral methods have not been used directly and without discretization and linearization to solve nonlinear partial differential equation with fractional derivative, widely. Because spectral methods convert the solving procedure of a nonlinear differential equation to solving of a nonlinear system of algebraic equations. In fact the major difficulty of spectral methods to nonlinear problems is solving this system of nonlinear equations. One of the best methods to solve a nonlinear system of equations is the Jacobian free Newton-Krylov sub-space methods (JFNK) [22,23]. Now, in this paper, we introduce a new method to solve nonlinear time-fractional Burger equation by using the spectral collocation method based on the Bessel functions and generalized Bessel function of the first kind and JFNK method with adaptive preconditioner.

Fractional Burger equation
The Burger equation has many application in traffic flow, shock waves in a viscous medium, gas dynamics, etc. [2,5]. The Burger equation is a simplified version of Navier-Stokes equations. In 1939 the dutch scientist J.M. Burgers simplified the Navier-Stokes equation by just dropping the pressure term and considered without external force, this equation can be investigated in one spatial dimension [24,25]. Now, here we consider the time-fractional Burger equation that is defined as follows [1,5] with initial condition u(x, 0) = u 0 (x) and boundary conditions u(0, t) = g 0 (t) and u(l, t) = g 1 (t),that u is the velocity, ν is the kinematic coefficient, f (x, t) is an external force and the D α t , 0 < α < 1 is the Riemann-Liouville fractional derivative operation: where: Solving and finding the solution of the time-fractional Burger equation (1) has been studied for the last half century and still it is an active area of research to develop some better numerical algorithms and methods to approximate its solution [3,5,14,15,18,26].
The rest of this paper is organized as follows: The basic definitions and properties of the first kind of Bessel functions and generalized Bessel function, are presented in Section2. The function approximation based on Bessel function and generalized Bessel function via spectral methods to solve a nonlinear FPDE are described in Section 3. Then, to show the advantages, applicability and reliability of the proposed method we solve some examples of time fractional Burger's equation and compare our results with others in Section 5. Finally, the paper concludes in Section 6.

Bessel function and generalized Bessel function of the first kind
In this section we explain the Bessel functions and generalized Bessel functions of the first kind and some useful relations of theirs to use in spectral methods. The Bessel function of the first kind Jn(x) is defined as follows: where Γ(λ) is the Gamma function: The series (2) is convergent for all −∞ < x < ∞. Actually, the Bessel function is a solution of the following Sturm-Liouville equation [27]: for x ∈ (−∞, ∞), (n ∈ R).
It is clear if n ∈ N, Jn(x) is linear independent and also, the integration of this function is infinite [27] Lemma: One of the useful recursion relations of Bessel function of the first kind is: Proof. By deriving of (2) and using expansions of J n−1 (x) and J n+1 (x), the result is desirable.

Remark:
The derivative operational matrix of the first kind Bessel functions can be obtained as follows: where D is derivative operational matrix and is obtained by (4): where the a 0 , a 1 , a 2 , ..., an will be obtained by an interpolation technique. Now we define generalized Bessel function (GBF) of the first kind as follows: Definition: Let n > −1 then generalized Bessel function (GBF) of the first kind is defined as:︀ Lemma: A recursive relation of derivative of GBF of the first kind is as follows, Proof. By using GBF of the first kind definition and expansion of Jn(x) the result will be achieved.
Remark: the derivative operational matrix of the first kind GBF can be obtained as follow: wherê︀ D is derivative operational matrix and is obtained using (7): where the b 0 , b 1 , b 2 , ..., bn will be obtained by an interpolation technique. Notice that:̂︀ x Theorem: If α > 0 and n > −1 then the Riemann-Liouville fractional integral of generalized Bessel function of the first kind is: and the Riemann-Liouville fractional derivative of GBF is: Proof. By using definition of Bessel function of the first kind (2) and Riemann-Liouville fractional integral of︀ Jn(x), we can write: t r+n 2 2r+n dt by calculating and relations of Gamma function for real values, we have: Which ultimately results x r+n+α 2 2r+n+α = 2 α̂︀ Jn+α(x). (12) Also, the Equation (10) can be concluded immediately in the same way.

Function approximation and spectral method
Let ω is a certain weight function and Λ = {x|a < x < b}, therefore: In other words, let y be an arbitrary element in L 2 (Λ) × L 2 (0, T), since J N ×̂︁ J M is finite dimensional subspace of L 2 (Λ) × L 2 (0, T), y has a unique best approximation y 0 ∈ J N ×̂︁ J M such that ∀v ∈ J N , ‖y − y 0 ‖ ≤ ‖y − v‖.
assume that the u ∈ H r (Λ) × L 2 (0, T), r > 1, is the solution of the equation (1) where K is a constant related to N, M and C is a positive constant depending only on the norms of u in the space mentioned. Now we use these principles in spectral methods to solve fractional PDEs

Collocation method to solve time-fractional Burger equation
Spectral methods, in the context of numerical schemes for solving differential equations. In spectral methods, we suppose that, the solution of differential equation belongs to Hilbert space H 1 × H 2 , and then we approximate it by (13) and finally via an appropriate method find the coefficients of expansion. For more clarification, consider the approximation of the time-fractional Burger (1) equation: with enough initial and boundary conditions. The starting point of the spectral method based on BFs and GBFs is to approximate the solution u(x, t) by a finite summation: x ∈ Λ , t ∈ [0, T].
In spectral methods, the coefficients must be determined.
Substituting u with u MN in (14) leads to the residual function: where The notion of the WRM is to force the residual function to zero by requiring: where {ψ k } are test functions, and ω is positive weight function. The choice of test functions results to a kind of the spectral methods. A method for forcing the residual function to zero, is the collocation algorithm [21,32]. In this method, by choosing the Lagrange basis polynomials constructed on the collocation points tp and xq as the test function, such that ψ 1 i (xq) = δ iq and ψ 2 j (tp) = δ jp and applying the Gauss quadrature rule can write:

Newton -Krylov algorithm
Solving a nonlinear differential equation by spectral method directly (without linearization or discritization) leads to solving a large nonlinear system of equations f 3 (x), ..., fn(x)) T and x ∈ R n is a vector. So speed and accuracy of solving this nonlinear system is very important. Many works have been done to improve solving the nonlinear systems [21][22][23]. One of the best methods to solve a nonlinear system is the classical Newton iterative method: presuppose x n+1 be root of F(x), F(x n+1 ) = 0 : is the n × n Jacobian matrix and is defined as follows: therefore: In fact in each iteration, a linear system must be solved: The obtained nonlinear systems of spectral methods are usually large and ill-condition and for large and complicated nonlinear systems, calculation and reordering J(xn) and solving obtained linear system in each iteration could be most time consuming. One good idea is to use the finite difference technique to approximate the Jacobianvector product: where ε is a very small value. Jacobian-vector product, can be useful to approximate the Jacobian matrix. Also, for large dimensions, iterative methods such as GMRes or BiCGSTAB are preferred over direct solvers [22,23]. In this paper, we use the Jacobian-free Newton-GMRes method to solve the obtained nonlinear systems of equations [21][22][23] of collocation method.  [5,33] to solve example 1.
Proposed method B-spline FE [5] Haar wavelet [33] and exact solution u(x, t) = t 2 e x . This example has been solved by finite element (FE) [5,26] and by using Haar wavelet [33]. Now, we apply the proposed method to solve this example. To approximate the solution and satisfying the initial condition, we use the following expansion: and boundary conditions are satisfied in nonlinear system of equations. The L 2 and L∞ errors of obtained approximated solution by proposed method with M = 6, N = 11 are compared with results of FE [5] and Haar wavelet [33] in Table 1 and the convergent rate of the proposed method to solve this example for α = 0.5 is shown in and exact solution u(x, t) = t 2 cos(πx). To approximate the solution and satisfying the initial condition, we use the fol- and boundary conditions are satisfied in nonlinear system of equations. The L 2 and L∞ errors of obtained solution via proposed method with M = 6, N = 13 are compared with results of FE [5] in Table 3 and the convergent rate of the proposed method to solve this example for α = 0.5 is shown in Table 4.
and exact solution u(x, t) = t 2 sin(2πx). To satisfy the conditions u(x, 0) and u(0, t) we approximate the solution as follows: and the another boundary condition is satisfied in nonlinear system of equations. The L 2 and L∞ errors of obtained solution via proposed method with M = 11, N = 14 are compared with results of FE [5] in Table 5 and the convergent rate of the proposed method to solve this example for α = 0.5 is shown in Table 6. Example 4: Consider the equation (1) with conditions: and exact solution u(x, t) = e −t sin(πx). Now we apply the proposed method to solve this example for α = 0.5. To satisfy the conditions u(x, 0) and u(0, t) we approximate the solution as follows:    and the another boundary condition be satisfied in nonlinear system of equations. The error of obtained solution via proposed method with M = 12, N = 9 is compared with results of [18] in Table 7 and the convergent rate of the proposed method is shown in Table 8.

Conclusions
In this paper, the generalized Bessel function is introduced to simple calculations in fractional equation then a fully spectral collocation method based on Bessel function and generalized Bessel functions is introduced to solve the nonlinear time-fractional Burger equation. In this method the nonlinear time-fractional Burger equation is converted  to a nonlinear system of algebraic equations by using a new spectral method based on BFs and GBFs, without any discritization and linearization method. Afterward, we use Jacobian free Newton-GMRes method to solve this nonlinear system of equation. As indicated in the presented examples, the solutions of the nonlinear systems are obtained in 3, 4 and 11 Newton iterations, also in all examples the initial guess of Newton method is simple vector [0, 0, ..., 0] T , that show the speed and the power of the proposed method. Also the shown L 2 and L∞ errors in the presented tables and comparison with other methods show efficiently, applicability and reliability of the proposed method.