An Efﬁcient Numerical Scheme for Variable-Order Fractional Sub-Diffusion Equation

: The variable-order (VO) fractional calculus can be seen as a natural extension of the constant-order, which can be utilized in physical and biological applications. In this study, we derive a new numerical approximation for the VO fractional Riemann–Liouville integral formula and developed an implicit difference scheme (IDS) for the variable-order fractional sub-diffusion equation (VO-FSDE). The derived approximation used in the VO time fractional derivative with the central difference approximation for the space derivative. Investigated the unconditional stability by the van Neumann method, consistency, and convergence analysis of the proposed scheme. Finally, a numerical example is presented to verify the theoretical analysis and effectiveness of the proposed scheme.


Introduction
Fractional differential equations (FDEs) got the researchers attention and this fact states the phenomena in various areas such as biochemistry, wave propagation, medicine, anomalous diffusion, and electrical engineering, etc. [1][2][3][4]. The fractional diffusion equation is an important FDE, which has been hugely used in fractional diffusive systems, diffusion in transport processes, unification of diffusion phenomena and fractional random walk [5,6]. The fractional sub-diffusion equation (FSDE) is a subclass of the fractional diffusion equation and obtained by replacing the time derivative by fractional-order α lying between 0 and 1. The one-dimensional FSDE can be written as: where ∂ 1−α ∂t 1−α represent the Riemann-Liouville fractional derivative (RL-FD) of order 1 − α. Most of the FDEs are not easy even impossible to get an analytical solution, which contains complicated functions. Therefore, numerical methods have become the main way for the solution of such types of fractional order differential equations. Many researchers have constructed the various numerical methods to solve FDEs such as Yuste [7] considered the weighted average implicit difference scheme (IDS) for the solution of FSDE. They discussed the theoretical analysis by the well known Fourier series method. Ali et al. [8] constructed a new numerical scheme for the modified FSDE. They investigated the theoretical analysis and found that the scheme is unconditionally stable and highly accurate. In another survey Ali et al. [9] derived a new approximation for RL-FD and developed compact IDS for modified FSDE. They successfully discussed the theoretical analysis and accuracy. Mohebbi et al. [10] studied the compact order IDS and solved the modified FSDE. They approximated the fractional derivative by Grünwald-Letnikov and the space derivative by fourth order approximation. The theoretical analysis are discussed by Fourier method and the tested examples shown high accuracy. Zhuang et al. [11] formulated the physical-mathematical method for FSDE and investigated a new numerical scheme. They analyzed the theoretical analysis based on energy method and the numerical values fully support the theoretical analysis. Ali et al. [12] proposed the Crank-Nicolson method for 2D-FSDE. They replaced the fractional derivative with Grünwald-Letnikov approximation and the space derivatives with finite difference approximation. They discussed the theoretical analysis with unconditional stability and convergence. Khan and Rasheed [13] studied the mixed convection in an incompressible fractional Maxwell nanofluid with thermal transmission. They tackled the arising mathematical problem numerically via finite difference procedure in their reported study they have mentioned that escalating values of Prandtl number reduces the thermal transport. Hamid et al. [14] investigated the modeling of unsteady radiating fractional nanofluid problem in a channel. Fractional derivative in Caputo sense is engaged in the model problem. They computed the solution of the boundary layer equations with the help of Crank-Nicolson numerical scheme. They found that the magnetic parameter reduces the velocity, whereas, Grashof number boosts fluid velocity. Other related literature can be seen in [15][16][17][18]. The above cited studies are reported the fractional derivative of constant-order.
The constant-order fractional derivatives are sometimes not more suitable to explain the complex diffusion processes as in porous medium and medium structure because it changes with time [19,20]. In VO operator the order may vary as function of independent variables such as time or space or both time-space. Lorenzo and Hartly [21] discussed some physical phenomena which show that the fractional-order behavior may change with time and space. To solve such type of VO fractional model and developed some numerical techniques. Here, we include only the VO-FSDE, which can be written as: There are many researchers who have worked on a VO-FSDE and find the solution by different numerical techniques such as Chen et al [22] developed a numerical scheme for VO-FSDE with first-order temporal and fourth-order spatial accuracy. they investigated the stability, convergence, and solvability analysis by Fourier series method. Lin et al [23] consider the nonlinear VO-FSDE and construct a new numerical finite difference scheme, also discussed stability and convergence analysis. Sun et al. [24] discussed different complex diffusion processes that can be described by the VO model. They differentiate this model into four types such as time-dependent, space-dependent, concentration-dependent, and system parameter dependent models. Sweilam et al. [25] solved the VO Riesz space fractional wave equation by explicit difference method. They discussed the Riesz space fractional derivatives based on shifted Grünwald-Letnikov formula and also established the stability and convergence analysis. Chen et al. [26] formulated the numerical solution for the nonlinear VO-FSDE. They proposed the piecewise function for VOFD and then by operational matrices transformed into an algebraic equation. Bhrawy and Zaky [27] investigated a numerical method based on the collocation method combined with Jacobi operational matrix and successfully solved the one and 2D-VO fractional Cable equations. They analyzed the convergence analysis and the numerical results demonstrated that the algorithm is powerful with high accuracy. Chen [28] solved the 2D-VO modified FSDE by numerical method. They discussed the stability, convergence, and solvability by the Fourier method. The two new second order numerical approximations are derived for VO time-fractional operator by Zhao et al. [29]. They tested on VO-FSDE and super-diffusion problems, the results are different than standard diffusion and constant-order diffusion for the VO α(x, t).
Wang et al. [30] considered the VO-FDE with variable coefficient. They solved the problem by finite difference method and homotopy regularization method based on Legendre polynomials as the basis functions of approximation space. Ma et al. [31] discussed the numerical Adams-Bashforth-Moulton method for the VO fractional financial system of equations. The tested examples reported that the proposed method for VO-FDE is simple and effective. Shekari et al. [32] considered the moving least squares method to find the solution of 2D VO fractional diffusion-wave equation. They utilized the mesh free and collocation method based on moving least square approximation. Xu et al. [33] solved the multi-term VO space-time fractional diffusion equation on a finite domain by numerical method. Investigated the theoretical analysis via mathematical induction and numerical examples confirmed the efficiency and accuracy. Shen et al. [34] considered the Coimbra VO time derivative which is a Caputo type definition and used in the numerical scheme for variable -order fractional diffusion equation. The theoretical analyses are discussed by the well know Fourier series. Further, other related studies of fractional and VO-FDE can see in [35][36][37][38][39][40].
In this article, we developed a new implicit difference method to solve the VO-FSDE by applying the new discretized approximation for VO-RL fractional integral. The first order time derivative is replaced by backward difference approximation. Further, investigated the stability analysis by van Neumann method, consistency and convergence of the proposed method. The proposed scheme is very easy to implement, it decreases the computational complexity and increase efficiency. Another positive point of this scheme is, the easiness in finding the stability and convergence analysis.
This investigation is arranged as follows: Section 1, contains the comprehensive literature survey relevant to the reported material. In Section 2, briefly explain the mathematical preliminaries, we discussed the implicit scheme in Section 3. In Sections 3.1 and 3.2, we use the van Neumann method to prove the stability and consistency respectively. In Section 4, provided the numerical experiments to confirm the theoretical analysis and in the last section concluded the conclusion.

Mathematical Preliminaries
In this section, we discuss the VO fractional calculus [41] and deriving new approximation for VO-RL fractional integral, which is used in this paper: Definition 1. The VO RL-FD formula can be written as: Definition 2. The VO-RL fractional integral formula can be written as: Here, we formulate a new numerical approximation for VO-RL fractional integral of order α(x, t) and (0 < α(x, t) < 1) at the grid point (x i , t k ), as following: can be defined in discretized form at the grid point (x i , t k ) as:

The Implicit Scheme
To construct an IDS for VO-FSDE (2), first we use Equation (4), then applying lemma 1 for the discretization of RL fractional integral and for space derivative utilizing the central difference approximation. Here, the x variable steps as x i = i∆x, in the x-direction with i = 1, ..., M − 1, ∆x = L M and for t variable the step is t k = kτ,k = 1, ..., N where τ = T N . Let u k i be the numerical approximation to u(x i , t k ), put Equation (4) into Equation (2), we have by applying lemma 1 and backward difference approximation to Equation (7), we get Simplifying Equation (8), we get a new IDS for the one-dimensional VO-FSDE (2) with the conditions, as follows Here, and i = 1, 2, ..., M x − 1, and k = 1, 2, ..., N − 1, with

Stability
To analyze the stability of the proposed IDS using van Neumann method. Let U k i represent the approximated solution for (9), we obtained The error is defined as where e k i satisfies Equation (13), as The error initial and boundary conditions are given by By defining the following grid functions for k = 1, 2, . . . , N.
As 0 < α k i < 1, from (22) and lemma 1, we have From the above proof, it is clear that Proposition 1 and Equation (20), reported that the solution of Equation (9) satisfies this proved that the IDS in (9) is unconditionally stable.

Consistency
To analyse the consistency of the IDS, let U be the exact and u be the approximate solution and F(U) = 0 represent the approximated difference equation of the FDE (2) at point (x i , t k ), then F(u) = T k i represented the truncation error at grid point (x i , t k ). Theorem 1. The truncation error T(x, t) of the proposed difference scheme is : Proof.
By Taylor expansion we can write and After simplification, we get the following truncation error This theorem shows that the proposed scheme is consistent if ∆x → 0 and ∆t → 0, then the truncation error tends to zero.

Theorem 2.
According to Lax equivalence theorem, consistency and stability are both necessary and sufficient for convergence [39], the proposed VO-IDS is convergent.

Numerical Experiments
The numerical example of VO-FSDE are presented in this section to confirm the effectiveness of the new IDS. Here, we compute the E ∞ and E 2 norm error, is defined as follows: Example 1. Consider the VO-FSDE is as following [19]: , the initial and the boundary conditions are The exact solution is u(x, t) = e x t 2 .
Tables 1 and 2 shows the numerical results of the new VO IDS for Equation (26). In Figure 1, the numerical values are compared with the previous study in Cao et al. [19] for various values of time steps and fixed space step. The proposed VO scheme reported better level of accuracy and more efficiency. Table 2 show that the error is reduced as the number of space and time steps increased. , and y = 0.25, 0.083 and T = 1.0 respectively. The above discussion proved that new IDS for VO-FDEs is accurate and numerically more efficient.

Conclusions
In this article, first we derived a new approximation for VO-RL fractional integral and for space derivative used the implicit central difference approximation. Successfully utilized the new numerical approximation and the IDS for one-dimensional VO-FSDE. Investigated the theoretical analysis by the well-known von Neumann method. The scheme is unconditionally stable and convergent with order O(τ + (∆x) 2 ). The tested example result has compared with the previous study and with the exact solution for various values of α(x, t). We concluded that the scheme is a more effective, high level of accuracy and the solution is well-matched. We believe that this finding is another contribution to the literature. This numerical approximation can also be applied to other types of variable-order fractional differential equations. Acknowledgments: The authors are thankful to the reviewers for the useful suggestions to improve the quality of this manuscript.

Conflicts of Interest:
The authors declare no conflict of interest.