A new formulation of finite difference and finite volume methods for solving a space fractional convection–diffusion model with fewer error estimates

Convection and diffusion are two harmonious physical processes that transfer particles and physical quantities. This paper deals with a new aspect of solving the convection–diffusion equation in fractional order using the finite volume method and the finite difference method. In this context, we present an alternative way for estimating the space fractional derivative by utilizing the fractional Grünwald formula. The proposed methods are conditionally stable with second-order accuracy in space and first-order accuracy in time. Many comparisons are performed to display reliability and capability of the proposed methods. Furthermore, several results and conclusions are provided to indicate appropriateness of the finite volume method in solving the space fractional convection–diffusion equation compared with the finite difference method.


Introduction
The mathematical model is one of the most important logical aspects used to help us explain and study physical phenomena. Many physical situations include non-Markov processes, earthquakes, and terrestrial surface transport can be modeled through a class of diffusion equations or wave equations. In particular, the physical quantities are transferred inside physical systems by diffusion and the flow of heat and stochastically changing systems are described. In recent years, researchers have focused on finding effective numerical methods to solve such a type of partial differential equations (PDEs). Different PDEs have been studied and resolved including the space-time fractional diffusionwave equation [1], the space fractional convection-diffusion equation [2], the space fractional Cauchy equation [3], the fractional heat-and wave-like equations [4], the fractional Korteweg-de Vries equation [5], and so forth [6][7][8][9][10][11]. Different effective methods have been proposed by scholars such as the reproducing kernel method [12][13][14][15], the residual power series method [16][17][18], the homotopy perturbation method [19], the multi-step approach [20], and some others [21][22][23][24].
During the past years, certain numerical methods, such as the finite difference method, the finite volume method, and the finite element method, have been developed to solve some PDEs. The important goal of these methods is to convert the equation to a system of linear algebraic equations. For example, in [25], the time-fractional convection-diffusion equations have been reduced to a system of linear algebraic equations using the Gegenbauer polynomials and Lagrange polynomials space and time expansion, respectively. In [26], Saadatmandi et al. studied the sinc-Legendre collocation technique for solving a class of fractional convection-diffusion equations with variable coefficients. In [27], the finite volume element technique was used to solve the sub-diffusion problems in fractional order. In addition, Jia and Wang used the volume penalization approach to improve the finite volume technique for finding a solution for a space fractional diffusion equation [28]. In [29], Mahto et al. applied the continuation property to propose the approximate solution of the sub-diffusion equation via internal control. In [30], Zhao et al. applied the finite element technique for solving the time-fractional diffusion problem, they introduced an unconditionally stable scheme based on the spatial quasi-Wilson nonconforming; while in [31], the authors derived the optimal error estimates by studying the two mixed finite element methods, conforming and nonconforming, to solve the time-fractional diffusion problem. On the other hand as well, Stynes et al. used a different technique to solve the time-fractional diffusion equation and to find elaborate error boundedness based on the regularity of the exact solution [32]. In [33], Sayevand et al. proposed the finite volume element technique for analyzing the behavior of sub-diffusion equations and proved its stability analysis of obtained solution. In [34], the authors proposed the finite difference discretization based on the weighted trapezoidal rule to approximate the fractional Laplacian and proved the convergence in one dimension. More interesting applications can be found in [35][36][37][38][39][40][41][42][43][44][45].
In this paper, we deal with a new aspect of FVM and FDM to obtain approximate solution for the space fractional convection-diffusion equation (SFCDE), in which the space fractional derivative is represented in the Riemann-Liouville sense. More specifically, let us consider the SFCDE of the following form: subject to the initial condition where t ≥ 0, x ∈ [0, L], a parameter β describes the order of the space fractional derivative, 0 < β ≤ 1, so that the space fractional derivative is described in the Riemann-Liouville sense, ψ(x, t) is concentration, and ρ are positive parameters, ρ is a diffusion coefficient, is the velocity of water flow, and L is the length of the channel. For β = 1, the equation, indeed, reduces to the linear convection-diffusion equation. In particular, the exact solution of the linear convection-diffusion equation of the shape subject to the initial condition ψ(x, 0) = -sin(πx) is given in the form The structure of this analysis is as follows. In Sect. 2, some preliminary definitions and theorems related to Grünwald-Letnikov fractional derivative are introduced. Section 3 is devoted to constructing the FDM for FIVPs (1) and (2) along with its theoretical analysis. In Sect. 4, the FVM is proposed to solve FIVPs (1) and (2) along with its theoretical analysis. In Sect. 4, numerical examples are given to express validity and reliability of the methods and to make a comprehensive comparison for both FVM and FDM. Finally, some conclusions are given in Sect. 5.

Basic definitions
In the literature there are numerous integrals which are widely used in physics, astronomy, and engineering [46,47]. Consequently, different types of fractional integrals and derivatives have been found very reliable under some conditions. In this section, we introduce some preliminary definitions and theorems related to the Riemann-Liouville fractional derivative and Grünwald-Letnikov fractional derivative.

Definition 1 The Riemann-Liouville integral of fractional order
For α = 0, the Riemann-Liouville integral of fractional order J 0 a ψ(x) = ψ(x) is the identity operator.
Definition 2 Let n be the smallest integer that exceeds α, then the Riemann-Liouville fractional derivative of order α > 0 is defined as For α∈ N, the Riemann-Liouville fractional derivative of order α is given as Then the Grünwald-Letnikov fractional derivative of the function ψ of order α is given bỹ where h = x-a N .

For
The result is, therefore, positive. Part 3 is straightforward. Hence, we omit the details. This completes the proof of the lemma.

Finite difference scheme formulation
In this section, we propose a new finite difference scheme formulation (FDSF) for solving the SFCDE with constant coefficients subject to the initial condition is an analytical smooth function of a spatial x, ψ(x, t) is an unknown analytical function, and ρ are positive parameters, and β is a parameter that describes the order of the space fractional, where the space fractional derivative is described in the Riemann-Liouville sense. If β = 1, then the equation becomes a linear convection-diffusion equation. Now, by utilizing Definition 2 and Riemann-Liouville fractional derivative, we have where where τ is the time step and h is the spatial step. Then, we will be using the standard Grünwald formula to approximate the α order fractional Riemann-Liouville integral and the central difference formula to approximate the first and the second derivative as follows: ∂ψ(x, t) ∂x The starting point for applying the finite difference method by evaluating Eq. (10) at x = x i and using the above equations is that Hence, using the standard backward difference to approximate the temporal derivative in Eq. (14) reveals Letting ψ n i ≈ ψ(x i , t n ) denotes that the numerical solution implies By collecting like terms, Eq. (16) can be rewritten as Denoting the numerical solution vector by ψ n = [ψ n 0 , ψ n 1 , . . . , ψ n N ], the vector equation becomes where A is a matrix with elements a ij = b ij . Hence, the above given equation can be formatted as where M = (I + τ A) -1 is the iteration matrix. The next theorems show that the finite difference scheme FDS expressed in Eq. (20) is conditionally stable and consistent with first-order accurate in time and second-order accurate in space.
We commence to prove the strictly diagonally dominant property for the iteration matrix M.
where h is the spatial step. Then the coefficients a ij satisfy the inequality Proof Let i = 1, 2, . . . , N be given. Then consider the sum Therefore, the previous lemma and the hypotheses ρ > h α+1 2 guarantee that each term is negative. So, we have Replacing the finite sum with the infinite sum yields Therefore, the telescoping sum has the form (v α This completes the proof of the theorem. (20) is convergent, so the scheme itself is conditionally stable.

Corollary 1 The iteration matrix M in Scheme
Proof By the previous theorem, we infer that A is strictly diagonally dominant with positive diagonal elements, and so is I + τ A. Therefore, the iteration matrix M = (I + τ A) -1 exists, and its spectral radius satisfies The proof is therefore completed. (17) is consistent with second-order accuracy in space and first-order one in time.

Corollary 2 Numerical Scheme
By the fundamental theorem of numerical methods for differential equations, the given numerical scheme is convergent.

Finite volume scheme formulation
In this section, we propose a novel finite volume scheme formulation (FVSF) for solving the SFCDE with constant coefficients subject to the initial condition where t ≥ 0, 0 < β ≤ 1, g(x) is an analytical smooth function of spatial x, ψ(x, t) is an unknown analytical function, and ρ are positive parameters, and β is a parameter that describes the order of the space fractional, where the space fractional derivative is described in the Riemann-Liouville sense. If β = 1, it becomes a linear convection-diffusion equation. Now, to establish the finite volume scheme, we need to partition the finite domain = [a, b] to N + 1 uniformly spaced nodes x i = a + ih, i = 0, 1, . . . , N , where the spatial step is h = b-a N . Thus, by utilizing the Riemann-Liouville fractional derivative, we have where J 1-β a is the Riemann-Liouville integral with respect to x. Setting α = 1β gives 0 ≤ α < 1. Integrating Eq. (23) over the ith control volume [x i-1 2 , x i+ 1 2 ] suggests Dividing each side by h gives the standard finite volume discretization where is the control volume averages of ψ(x, t).
Using the fractionally-shift Grünwald formula to approximate J α a ψ(x, t) yields Hence, the Riemann-Liouville integral includes the order term Therefore, due to certain needs of the fractional shift p = 1 2 , we write where to write For Eqs. (31) and (32), use the standard averaging scheme ψ(x i+ 1 2 to construct approximations of first derivative in terms of function values at the nodes x j : So, Eq. (25) can be approximated as Notice that if ψ(x, t) is a smooth function, then the value of the control volume averages ψ i (t) will agree with the value of ψ(x, t) at the midpoint of the interval [x i-1 2 , x i+ 1 2 ] to O(h 2 ). So, we can rewrite Eq. (35) as Now, define a temporal partition t n = nτ , n = 0, 1, . . . , where τ is the time step, and use the standard backward difference to approximate the temporal derivative in Eq. (36) such that Let ψ n i ≈ ψ(x i , t n ) denote the numerical solution. Then we have By collecting like terms, we rewrite Eq. (37) in the form where k ij is given by The vector of the numerical solution ψ n = [ψ n 0 , ψ n 1 , . . . , ψ n N ] can be denoted by where the matrix A has the elements a ij = -k ij . The next two theorems explain that the finite volume scheme described in Eq. (37) is conditionally stable and consistent with first-order accuracy in time and second-order accuracy in space. However, we havê satisfies the von Neumann condition when |1r N j=0 k ij exp((ji)Iξ )| ≥ 1. Therefore, by using the reverse triangle inequality, we have Hence, the von Neumann condition will be satisfied whenever That is, either which is impossible to hold, or Thus, the symbol ρ(ξ ) of the numerical scheme satisfies the von Neumann condition if N j=0 k ij exp (ji)Iξ ≥ 2 r , ∀i = 0, 1, 2, . . . , N, which means that the numerical scheme is conditionally stable. This completes the proof of the theorem. (37) is consistent with second-order accuracy in space and first-order one in time.

Theorem 6 Numerical Scheme
Proof Insert ψ(x, t) into expression (37) to get By the control volume averages, combined with a temporal backward difference for the time derivative, we have Substitute (25), (27), and (30) above, it turns to T n i = O(τ + h 2 ). Hence, T n i → 0 as τ → 0 and h → 0. So, the given numerical scheme is consistent with second-order accuracy in space and first-order in time. This completes the proof of the theorem.

Corollary 3 Numerical Scheme (37) is consistent and stable.
By the fundamental theorem of numerical methods for differential equations, the given numerical scheme is convergent.

Numerical experiments
In this section, we introduce the numerical problems to verify the correctness and effectiveness of the proposed schemes. Both schemes can efficiently be solved, but the FVSF is more accurate than the FDSF. Two examples are introduced by comparing the numerical solutions with the exact solutions. The computations are performed by Mathematica Software 11.0.
Example 2 Consider the following SFCDE: subject to the initial condition where x ∈ [-2, 4], = 0.1, ρ = 0.02, t ≥ 0, 0 < β ≤ 1, and ζ = 1.1771243444677. For β = 1, the exact solution is u(x, t) = e ζ x-mt , m = 0.09. The FDSF and FVSF are tested numerically to solve the above problem for different values of β such that β ∈ {0.5, 0.7, 0.9, 1}, by  Table 3. When β = 1, t = 1, and x ∈ [-2, -1.25], the numerical results are reported in Table 4. These results demonstrate that the FVSF is more significantly accurate than the FDSF, and it can be observed through by the absolute errors between the numerical solution in both schemes and the exact solution. The 2D plot of the solutions behavior with varying β such that β ∈ {0.5, 0.7, 0.90} is drawn in Fig. 4 at t = 1 and x ∈ [-2, 1] using the FDSF and the FVSF, while Fig. 5 shows that the FVSF provides an excellent agreement with the exact solution rather than the FDSF. Moreover, the 3D plot of the exact solution behavior and the solutions behavior using the FVSF is drawn in Fig. 3 at β = 1, t ∈ [0, 1], and x ∈ [-2, 1]. Figure 6 shows the surface plot of the solution behavior at β = 1 by using the FVSF and the FDSF at t ∈ [0, 1] and x ∈ [-2, 4]. Table 3 and Table 4 depict that the solution obtained by the FVSF is more accurate and closer to the exact solution compared with the solution obtained by the FDSF in two different experiments. Also, it is easy to increase the nodes in the regular grid to get higherorder accuracy. Thus, it can be concluded that the FVM is effective for solving the SFCDE with a constant coefficient. This implicit scheme can be improved for implementation when tackling differential equations, in specific for problems with or without source terms that are deformed or nonlinear. Also, it can handle general condition without difficulty and can deal with unlimited kinds of the initial condition.

Conclusions
The SFCDE is one of the most important equations used to describe physical phenomena. In this paper, the FDM and the FVM have been developed for solving the SFCDE. The space fractional derivative in this equation was considered in the Riemann-Liouville sense. The Grünwald-Letnikov fractional formula has been used for the discretization of space fractional derivative, and in both formulations the implicit scheme has been introduced. It has been shown that both schemes are stable with first-order accuracy in time and second-order accuracy in space. Numerical results have been found using Mathematica software 11.0. Comparing illustrated examples justified the efficiency and accuracy of the proposed methods.