Solution of the Space Fractional Diffusion Equation Using Quadratic B-Splines and Collocation on Finite Elements

: In this paper, we consider the solution of the space fractional diffusion equation using orthogonal collocation on finite elements (OCFE) with quadratic B-spline basis functions. The main advantage of quadratic B-splines is that they have good interpolating properties and can be easily adapted for solving problems on non-uniform grids. The method is unconditionally stable and its convergence is also discussed. It is of order (3 − α ) for 1 < α < 2. We present various linear and nonlinear examples. The solutions compared favourably with previous results in the literature


Introduction
Fractional differential equations are extension of integer order ordinary and partial differential equations to that of a fractional order.This can be fractional ordinary differential equation (FODE) which involve an independent variable or fractional partial differential equations (FPDE) that consists of several independent variables.Many real life phenomena have been modelled by fractional differential equations.It has applications in long range spatial interactions, anomalous diffusion [1], finance, turbulent flow, physics, chemical kinetics, biology and so on.
The fact that many fractional PDEs do not possess exact solutions has led to proposition of various of numerical methods to solve FPDEs.A finite difference method based on implicit Euler method was used to solve the fractional advection-dispersion flow equation in [2].Weighted and shifted Grunwald difference operators were used to obtain the numerical solution of the space fractional diffusion equation in [3].Discrete Galerkin finite element method was used for solving time-fractional diffusion equations in [4].The fractional derivative considered were in the Caputo sense and the method was found to be of order (2 − α) where α ∈ (0, 1).Spectral method and finite volume methods for fractional partial differential equations were discussed in [5] and [6] respectively.
The authors of [1] discussed the stability and convergence of quadratic spline collocation method for timefractional diffusion equation.A combination of method of lines and spline approximation were used to obtain the numerical solution of space fractional diffusion equation by [7].A cubic B-spline collocation method was employed to solve time fractional wave equation in [8] and anomalous fractional diffusion equation in [9].A modified spline collocation was proposed and studied by [10] for linear fractional differential equations.
A time fractional anisotropic diffusion model and numerical schemes for the model were proposed to improve the quality of medical ultrasound images in [11] via meshless method and semi-implicit time integration.A combination of the meshless and Galerkin methods were applied in [12] to two-dimensional time-fractional Tricomi type equations.Collocation coupled with the meshless method was employed by [13] and [14] to solve problems on modified anomalous diffusion and reaction-diffusion processes as two-dimensional fractional differential equations.A kernel based collocation method with greedy algorithm was used to investigate anomalous mobile-immobile solute transport models in [15] as fractional differential equations.In [16], a generalized anomalous mobile-immobile solute transport model with variable coefficients was studied using the meshless technique.
In this paper, we demonstrate the capability of orthogonal collocation using quadratic B-spline basis functions to solve space fractional diffusion equation effectively.The main advantage of quadratic B-splines is that they have good interpolating properties and can be easily adapted for solving problems on non-uniform grids.The method is unconditionally stable and its convergence is also discussed.It is of order The paper is arranged as follows: section 1 is for introduction, section 2 discuss the derivation of fractional B-spline derivative and in section 3, we applied the method to the space fractional diffusion equation.Section 4 deals with the stability of OCFE while section 5 deals with convergence analysis.Numerical examples are provided in section 6 and section 7 concludes the paper.

Derivation of fractional B-spline derivative
Let a fractional partial differential equation with boundary conditions (2) , where B k (x) is a quadratic spline basis function defined by We then apply continuity at internal boundaries such that (5) where h is the uniform interval spacing.This gives (6) Using the first condition in (7), the solution u(x) in (3) can be simplified to The boundary conditions in (2) become ( 9) We shall make use of the Caputo fractional derivative together with the quadratic B-spline basis functions in (4) to solve equation (1).The Caputo derivative [17] is defined as ) where α ∈ R and n ∈ N. In our case, we are dealing with diffusion equation, so n = 2. Equation ( 10) is now We re-write equation ( 8) as to include time dependence for a PDE.Substituting ( 12) into (11) yields (13) We then evaluate (13) at the collocation point z = 0.5 to get ( ) ( ) , ( 2) ( ) ( ) We then incorporate equation ( 15) into (1) with x = x j + 0.5h substituted into the right hand side to get N equations.These equations together with the continuity equations and boundary conditions give a square linear system of size (2N + 1) given by ( 16) The solution to the fractional differential equation ( 1) is then obtained via equation (12).The rate of convergence is calculated as where the components of e are given by e i = u i (x) − U i (x), u i (x) and U i (x) are exact and approximate solution at the nodes respectively.

Application to fractional diffusion equation
We apply the trapezoid rule to the fractional diffusion equation ( ) where j is the index for discrete time t, i is index for discrete x, Nt = (t f − t 0 )/∆t is the number of time steps, t f is the final time, t 0 is the initial time and ( , ) Hence the system can be represented in the form

Stability
Consider now the fractional diffusion equation given by where d(x, t) and p(x, t) are coefficient of the fractional derivative and the source term respectively.Consider the homogeneous case of (21) i.e p(x, t) = 0 and assume that d where the superscript n on the coefficients b refer to its values at t n .Substitute the collocation point x m + h/2 into (24), let u n 1/2 = u(x m + h/2), and by using (25) it can be shown that ( ) , where 0.

Convergence analysis
Having shown that our method is stable in the previous section, we now show that our method is consistent.The determinant of the coefficient matrix 0 for all values of α.Therefore the coefficient matrix is nonsingular.This implies that our method is consistent.Hence our method is convergent.
We now consider (1) and illustrate the convergence analysis for the fractional partial differential equation for a fixed t.This is equivalent to the fractional ordinary differential equation where λ is a constant.This simplifies to The collocation solution in the i th interval is Substitute the collocation points x ¯i, i = 1, 2, …, N into (36) and (37) to get [ ] Hence, substitute the collocation points x ¯i into (34) gives (40)

(
) ( ) so that the length of this interval is given by w i h for x k e s e s Le x ds ds x s x s e ds e ds x s x s where we have used the integral mean value theorem, ξ k ∈ (x k , x k+1 ) and We assume that for h << 1, ξ k can be approximated by x ¯k so that e''(ξ k ) = e''(x ¯k) = u''(x ¯k) − u c '' (x ¯k).By Taylor's theorem, [ ] ( ) We have assumed that the collocation solution at the midpoint of [x k , x k+1 ] is approximately the average of the solutions at the endpoints.Hence b 2k ≈ u c (x ¯k). ( [ ] We assume that the approximation of the second derivative in (43) agrees closely with that in (46), thus for some constant M k .Equation (47) has been verified numerically.We shall approximate the integral using the midpoint rule. (49) . 24 ( ) ( ) Similarly, the second term in (41) can be simplified to = then from (42), (48), ( 50) and (51) it follows that (52) ( 1) .24 ( ) ( ) The function on the RHS of ( 52) is an upper bound for |Le(x)| ∀x.In particular, when x = x N+1 , i = N, we have Thus there exists constants k 1 and k 2 such that (54) Now assuming that L −1 is bounded above, we have 54).This shows that the fractional partial differential equation ( 1) is of order (3 − α) in space.
We then infer from the trapezoid rule that the local error is given by for some constants K.

Numerical examples
In this section we illustrate our method with some examples as follows:

Contemporary Mathematics Volume 5 Issue 2|2024| 1245
In Tables 1 and 2, we compared our method with that of [18] for different values of α.It is clear from the tables that our method is better than that of [18].Table 3 shows that as α is increasing, absolute error is increasing.Table 4 shows that the convergence order is (3 − α).Figures 1 and 2 display the efficiency of quadratic spline to handle fractional diffusion equation using a small number of intervals.
The exact solution is u(x, t) = x 4 e −t .The linearized form of ( 61) is , 2 2 Graphical solution for this problem is similar to that of example 6.1.Hence the graphs are not repeated here.Tables 5 and 6 show the comparison between the present work and that of [18] for different values of α.It is clear that the present results are better than that of [18] due to lower infinity norm of errors.The exact solution is u(x, t) = (x 2 + 1) sin(t + 1).The linearized form of ( 65) is where m(x) = Γ(3 − α) x α−1 and r(x, t) = (x 2 + 1)cos(t + 1) − 2xsin(t + 1).We compared the infinity norm of errors for different values of α in Tables 7 and 8.We found that B-spline collocation approximation is better than corresponding finite difference approximation of the same order for the fractional diffusion equation.Table 9 shows that absolute error decreases as α increases except when α = 1.8.The convergence order is shown in Table 10.Ord ∞ represents the global order computed using the L ∞ norm.The nodal order is better than the global order.Figure 3 showed the approximate solution overlaid with a mesh plot of the exact solution.There is a perfect match between the exact and approximate solutions.The error is represented by Figure 4.
The exact solution is u(x, t) = (x 2 − x 3 )e −t .The linearized form of ( 69) is where m(x) = Γ(4 − α) x α and r(x, t) = 2αx 2 e −t .Figures 5 and 6 show the approximate solution and the error respectively when α = 1.5.In Table 11, as α is increasing, the absolute error is increasing.The convergence order in Table 12 is approximately (3 − α).Therefore the approximate solution shown in for α = 1.5 is of higher accuracy than that of α = 1.8 but less accurate than that of α = 1.3.The exact solution is u(x, t) = tx 2 .The linearized form of ( 73) is (76)   The approximate solution in Figure 7 is closest to the exact that in the first half of x interval and Figure 8 shows that the error propagates towards the end of the x-interval.The absolute errors in Table 13 decreases as α increases.Table 14 shows that for all 1 < α < 2, convergence order is approximately 2.

Conclusion
In conclusion, the OCFE method using quadratic B-Spline basis functions is computational efficient for solving the fractional space diffusion equation.The method is unconditionally stable and is convergent of order (3 − α), 1 < α < 2 for some of the examples.Various numerical simulations are presented in order to validate the robustness of the method.To the best of knowledge, this is the first time quadratic B-splines have been applied in the context of the OCFE method for solving the space fractional diffusion equation.

where β = 2 −
α, a jk is the coefficient of b k and

∫
Let x be in [x i , x i+1 ] and x ¯i be the midpoint of [x i , x]

Table 1 .
Comparison of L ∞ errors when h = ∆t at t = 1

Table 2 .
Comparison of L ∞ errors when h = ∆t at t = 1

Table 3 .
Table of absolute errors at for different values of α when N = Nt = 50

Table 4 .
Convergence rates for different values of α when N = Nt = 50

Table 5 .
Comparison of L ∞ errors when h = ∆t at t = 1

Table 6 .
Comparison of L ∞ errors when h = ∆t at t = 1

Table 7 .
Comparison of L ∞ errors when h = ∆t at t = 1

Table 8 .
Comparison of L ∞ errors when h = ∆t at t = 1

Table 9 .
Table of absolute errors at for different values of α when N = Nt = 50

Table 10 .
Convergence rates for different values of α when N = Nt = 50

Table 11 .
Table of absolute errors at for different values of α

Table 12 .
Convergence rates for different values of α

Table 13 .
Table of absolute errors at for different values of α when N = Nt = 50

Table 14 .
Convergence rates for different values of α when N = Nt = 50