A New Approach and Solution Technique to Solve Time Fractional Nonlinear Reaction-Diffusion Equations Inci

A new application of the hybrid generalized differential transform and finite difference method is proposed by solving time fractional nonlinear reaction-diffusion equations. This method is a combination of the multi-time-stepping temporal generalized differential transform and the spatial finite difference methods. The procedure first converts the time-evolutionary equations into Poisson equations which are then solved using the central difference method. The temporal differential transform method as used in the paper takes care of stability and the finite difference method on the resulting equation results in a system of diagonally dominant linear algebraic equations. The Gauss-Seidel iterative procedure then used to solve the linear system thus has assured convergence. To have optimized convergence rate, numerical experiments were done by using a combination of factors involving multi-time-stepping, spatial step size, and degree of the polynomial fit in time. It is shown that the hybrid technique is reliable, accurate, and easy to apply.


Introduction
The nonlinear reaction-diffusion equations have found numerous applications in pattern formation, in many branches of biology, chemistry, and physics [1][2][3][4].Reactiondiffusion (RD) equations have also been applied to other areas of science and can be successfully modelled by the use of fractional order derivatives.[5][6][7][8][9][10][11][12][13][14][15][16][17][18]; for example, the RD equations are employed to describe the CO oxidation on Pt (110) [5], the study of Ca 2 + waves on Xenopus oocytes [11], and the study of reentry in heart tissue [7,13].A great deal of effort has been expended over the last 10 years in attempting to find robust and stable numerical and analytical methods for solving fractional partial differential equations of physical interest.There has also been a wide variety of numerical methods, for example, finite difference techniques, finite element methods, spectral techniques, adaptive and nonadaptive algorithms, and so forth, which have been developed for RD's numerical solution [19,20].
The differential transform method was used first by Zhou [21] who solved linear and nonlinear initial value problems in electric circuit analysis.This method constructs an analytical solution in the form of a polynomial.It is different from the traditional higher order Taylor series method, which requires symbolic computation of the necessary derivatives of the data functions.The Taylor series method computationally takes long time for large orders.The differential transform is an iterative procedure for obtaining analytic Taylor series solution of ordinary or partial differential equations.The method is well addressed in [22][23][24][25][26]. Recently, the application of differential transform method is successfully extended to obtain analytical approximate solutions to ordinary differential equations of fractional order.Fractional differential transform method (FDTM) is a method that Arikoglu and Ozkol [23] developed for solving linear and nonlinear integrodifferential equations of fractional order.This method solves problems with high accuracy while constructing semianalytic solutions in the polynomial forms.FDTM is based on classical differential transform method, fractional power series, and Caputo fractional derivative [27].Arikoglu and Ozkol [23] tested their approach on several examples and the results obtained are in good agreement with the existing ones in the open literature.Momani et al. [28] presented a new generalization of the differential transform method that extended the application of the method to differential equations of fractional order.The new technique is named generalized differential transform method (GDTM) and is based on one-dimensional differential transform, generalized Taylor's formula [29], and Caputo fractional derivative [27].
Fractional partial differential equations (FPDEs) are also an interesting and an important topic.The fractional derivatives and integrals have been occurring in many physical and engineering problems with noninteger orders.Fractional calculus is based on the definition of the fractional derivatives and integrals.They play a major role in engineering, physics, and applied mathematics.FPDEs are used to model complex phenomena since the fractional order differential equations are naturally related to the systems with memory and nonlocal relations in space and time which exist in most physical phenomena.Fractional order differential equations are as stable as their integer order counterpart.One of the fundamental equations of physics is the Schrödinger equation which describes how the quantum state of physical system changes with time.The fractional Schrödinger equation provides us with a general point of view on the relationship between statistical properties of quantum mechanical path and structure of fundamental equations of quantum mechanics [30].In other words, the fractional quantum mechanics includes the standard quantum mechanics as a particular case.So FPDEs are obtained by generalizing differential equations to an arbitrary order.There are three popular methods for seeking approximate solutions for FPDEs which are the finite difference, finite element and spectral methods.In the literature there are many papers on these three methods.In these papers, the authors proposed the use of least squares finite element solution and fully discrete Galerkin method to solve nonlinear space fractional partial differential equations [20,31].
The differential transform is well suited to combine with other numerical techniques, as shown by Yu and Chen [32] who applied the hybrid method to solve the transient thermal stress distribution in a perfectly elastic isotropic annular fin.Kuo and Chen [33] employed the hybrid method to solve Burger's equation for flow systems with high Reynolds numbers.This method was also employed to analyze the dynamic response of an electrostatically actuated micro fixedfixed beam [34].
In the current study, the hybrid generalized differential transform/finite difference method is used for solving time fractional nonlinear RD equations.The validity of the proposed approach has been confirmed by comparing the results derived in the literature using the GDTM method [19], homotopy perturbation method (HPM) [35], and fractional variational method (FVIM) [36].
There are several approaches of definitions for the fractional derivative.Among them, one is called Riemann-Liouville fractional derivatives and defined by Here   is the -order Riemann-Liouville integral operator which is expressed as follows: If we use this definition, we must know the initial value of some fractional order derivative of the unknown function or we must have homogenous initial conditions.Unlike the Riemann-Liouville approach in the case of the Caputo derivative there are no restrictions on the initial conditions.Thus the following definition is used in this study: ( This operator is generally called "-order Caputo differential operator" [37,38].
For Caputo derivative we have

Generalized Differential Transform/Finite Difference Method
We firstly introduce the main features of GDTM [19,23,28,[39][40][41] according to the generalized differential transform of the th derivative of a function of one variable defined as follows: where 0 <  ≤ 1, (  )  =     ⋅ ⋅ ⋅    times.In (5), () is the original function,   () is the transformed function, and the differential inverse transform of   () is defined as follows: In case of  = 1, GDTM reduces to the classical DTM.From definitions ( 2) and (3) all fundamental properties of GDTM can be obtained easily [19,39].Since lim  → 1    =  has been proved from the definitions of fractional calculus, the fractional solutions   (, ) reduce to the standard solution (, ).
In this study, we use a hybrid method that is a combination of generalized, temporal differential transform and spatial finite difference methods to solve nonlinear fractional reaction diffusion equations.
We present a solution of a more general model of RD equation Mathematical Problems in Engineering 3 where  1 is the diffusion coefficient and () is a nonlinear function.We consider two different forms of () which are called time fractional Fisher equation and time fractional FitzHugh-Nagumo equation.
We apply GDTM to discretize fractional order time derivative and central difference method to discretize derivatives in  direction, respectively.After transforming (7) using the GDTM, we get the following equation: where   () is the transformed function of (  , ) and (  ) is the transformed function of ().The region 0 <  <  is divided into several equal intervals and each interval has a width ℎ.The time interval of interest is discretized using a time step Δ.After discretization of the equation we get solution at time Δ and these results are adopted as the initial values for the next time interval.This timestepping procedure assists in obtaining a converged solution to a desired accuracy [39].An appendix has been added to the paper to show that the error, if any, in the method is bounded.This implies stability of the scheme and by Lax equivalence theorem it thereby implies convergence.
The new algorithm has been developed to solve the nonlinear reaction diffusion equation and our aim of this approach is to combine the flexibility of differential transform and the efficiency of finite differences.This algorithm also provides an iterative procedure to calculate the numerical solutions; therefore, it is not necessary to carry out complicated symbolic computation.On applying the differential transform method with respect to time on the equation we are basically transforming the time-evolutionary equation to an elliptic type.In essence this means that the central finite difference approximation that is subsequently used on the transformed equation is a Poisson solver.The resulting system of linear algebraic equations is then diagonally dominant and hence the Gauss-Seidel iterative method used for solving the same has assured convergence as the coefficient matrix remains nonsingular throughout the computation.The algorithm used thus succeeds in segregating the time discretization from explicitly influencing the computation in the spatial domain and this presents a situation wherein the two can be handled independent of each other in the course of computation without having to bother about the stability of the solution if the differential transform part is properly handled.The latter is achieved deftly in the differential transform part of the algorithm by using the multistepping procedure as first enunciated by Yu and Chen [32] in their phenomenal work and used subsequently by Odibat et al. [29,39].In summary, this means that convergence is never in doubt in the algorithm but slow convergence can be if the time and spatial discretizations are badly handled.The convergence is optimized in the paper computationally by proper selection of time step in the differential transform part of the algorithm and then the spatial step size in the finite difference part of the algorithm.Next important step in the algorithm is the decision on the number of terms to be adopted in the inverse differential transform that gives us the solution of the problem as a power series in time.The aforementioned three vital components of the algorithm have been meticulously handled and a brief summary of the numerical experiment undertaken concerning the same is presented in a table.The numerical study recommends that the combination of 5, 10, and 50 spatial step size with a time step of 0.0005 or 0.001 assures the best rate of convergence if we take minimum ten terms in the time series.

Illustration of Generalized Differential Transform and Finite Difference Method
To show effectiveness of the proposed numerical solution using the temporal generalized differential transform and the spatial finite difference method and to give an understandable overview of the methodology, two examples of the reaction diffusion equations will be discussed in the following section.
Then our results will be compared with published work of Rida et al. [19] in which GDTM was used to solve the same equations.
Example 1.The time fractional Fisher equation is In this example, we have the nonlinear function () = 6(1 − ).The initial condition used is Operating the generalized differential transform on (9) gives us the following equation: where   (, ) is the generalized differential transform of (, ).
Now we apply the central finite difference method to the derivatives with respect to  and this gives us   The initial condition on discretization yields Equation ( 12) is a recurrence relation.The time series solution of the given equation is then obtained by using ( 12) and ( 13) with ℎ = 0.1 to obtain   ().Some of   () are recorded in Table 1.
The time series solutions of ( 12) with the initial condition (13) are obtained as follows: The numerical calculation results are shown in Figures 1 and 2, respectively.Our results are in agreement with the published work of Rida et al. [19], who considered the same equation.An exact solution of the standard form of Fisher equation for  → 1 is The comparison of our results with the exact solution is shown in Figure 1, for ℎ = 0.1, and quite clearly good agreement is found.Approximate solutions are shown in Figure 2 for  = 0.99 and  = 0.95.
The influence of  on the function (, ) is shown in Figure 3.This figure indicates a decrease in the fractional order  by choosing the fixed  = 5 that corresponds to an increase in the function and also indicates a slow diffusion for the values of  = 1 and  = 0.9 and a fast diffusion for the values of  = 0.8, 0.7, 0.6, respectively.It is clearly seen that (, ) increase for  = 1, 0.9, 0.8, 0.7, 0.6 with the increases in .
It is also found that the result is in complete agreement with the result of HPM [42,43] and ADM [44] for  = 1.
We investigate convergence criteria of our solutions for different values of ℎ and .To illustrate this, we compared our results with the analytical solution in case of  = 1.Here  is order of differential transformation method and denotes the number of terms to be calculated.
In Figures 4, 5, and 6 the difference between the results obtained in this study and the results of the analytical solution is of the order of 10 −5 .This is a pointer to the fact that there is convergence and is a restatement in numerical terms of what was shown in the Appendix.
One important observation made from the computation is that when the number of mesh points was increased, less number of terms was required in the time series solution to have convergence for a predetermined accuracy.The hybrid method of the present study gives faster convergence than other traditional methods; for example, if we take ℎ = 0.02 (mesh point is 50) then the solution converges for  = 3.We now consider another example.
Example 2. The time fractional FitzHugh-Nagumo equation is In this type of equation, the nonlinear function depends on  and it is () = (1 − )( − ).The initial condition is Using the hybrid method on the above initial boundary value problem (IBVP), as done in the previous example, we get Using second order finite difference method, the boundary values were obtained as follows: Table 3 presents some of the   ()'s.
Numerical solutions for the time fractional FitzHugh-Nagumo equation with various  values are shown in Figures 7 and 8.A comparison of the results in a limiting case wherein an analytical solution exists is shown in Figure 7.The results are in close agreement with those of Rida et al. [19] for the same equation.
For  → 1, it can easily be seen that the exact solution of FitzHugh-Nagumo equation is Figure 9 is prepared to show the influence of  on the function (, ).It is clearly seen that (, ) decrease for  = 1, 0.95, 0.85, 0.75, 0.65 with the decreases in .
As shown in the Table 4, our results show close agreement with the exact solution and agree with those of Rida et al. [19].From this table it is clear that the present work gives better approximation than GDTM as we increase .Numerical comparison between GDTM, FVIM, and hybrid method is shown in Table 5, which indicates hybrid method is more promising.

Conclusion
Many real physical problems can be best modelled with fractional differential equations but the fact is when the the fact that the error as a result of discretization and computation is bounded and hence implies stability of the method.Lax equivalence theorem further implies convergence of the scheme.Two time fractional nonlinear reaction-diffusion equations considered for illustration of the hybrid method highlight the usefulness of the method in obtaining the solution of IBVPs involving time fractional derivatives.The control of convergence through a judicious choice of time and spatial step sizes and also the number of terms in the time series solution spells assured convergence.The segregation of the time domain from the spatial domain in the solution method ensures the fact that problem of stability does not arise.Diagonal dominance of the coefficient matrix in the system of linear algebraic equations resulting from the use of the central difference approximation in the Poisson equation ensures the fact that the matrix remains nonsingular during iterations and hence has assured convergence.An appropriate computational decision on the number of terms to be taken in the time series solution results in a convergent solution with fast convergence.Excellent comparison of the present results with the previous works on generalized differential transform method [19] and homotopy perturbation method  [36].[35] and fractional variational iteration method [36] provides confidence in the methodology adopted for the solution of time fractional differential equations.

Figure 1 :
Figure 1: Numerical solution for the time fractional Fisher equation with  → 1 (a); comparison with the analytical solution (b).

Figure 4 :
Figure 4: Comparison of present results for ℎ = 0.2 and  = 10 with the analytical solution in case of  → 1 at  = 0.4.

Figure 9 :
Figure 9: Approximate solution for the time fractional FitzHugh-Nagumo equation with different  values at  = 5 (a) present and (b) [36].

Table 1 :
Some values of   () of Example 1.

Table 3 :
Some values of   () of Example 2.