A FRACTIONAL ANALOG OF CRANK–NICHOLSON METHOD FOR THE TWO SIDED SPACE FRACTIONAL PARTIAL EQUATION WITH FUNCTIONAL DELAY

For two sided space fractional diffusion equation with time functional after-effect, an implicit numerical method is constructed and the order of its convergence is obtained. The method is a fractional analogue of the Crank–Nicholson method, and also uses interpolation and extrapolation of the prehistory of model with respect to time.


Introduction
Diffusion equations with delay of general type, constant or variable, concentrated or distributed, are often applied in the simulation of dynamics processes. There are two effects often combined with these equations: distribution of parameters in space and heredity in time [1,2]. Numerical methods for solving such equations were considered in many papers, for example [3][4][5][6][7]. In paper [8], a technique of study on stability and convergence of numerical algorithms using the general theory of differential schemes was constructed for the heat conduction equation with delay [9] and the theory of numerical methods of the solution of the functional and differential equations were discussed in [10,11]. After that, this technique was applied to research of numerical methods of the solution of the equations of hyperbolic type with delay [12], various equations of parabolic type [13,14] and other types of the equations in partial derivatives with effect of heredity [15]. Here this technique is applied to the equations with partial derivatives of a fractional order effected with functional delay.
Fractional differential equations [16][17][18] gained a great interest in the past decades due to their accuracy in modeling a lot of problems in many fields of science. The equations in partial derivatives of a fractional order were subdivided into two big classes: with a fractional derivative on space and with a fractional derivative on time. Now there are many papers on numerical methods of the solution of such equations [19][20][21]. In this paper, we depend on results of work [22] in which explicit and purely implicit numerical methods for the solution of the equation with two sided space fractional derivative were investigated. In paper [23], an analog of Krank-Nikolson method was investigated for the equation which has a left-sided fractional derivative on space. In our paper, an analog of Krank-Nikolson method for the solution of the equation with two sided space fractional derivative and with functional time-delay is constructed and investigated. We need to extend this implicit scheme to work on the functional delayed form of these equations. The main purpose of this approach is to test the possibility of embedding the constructed implicit difference scheme for the fractional diffusion equation with functional delay into a general difference scheme.
After that, the convergence order is deduced by using the technique of proving similar statements for functional differential equations and methods from the general theory of difference schemes. We note that in work [24] this technique was announced for the numerical method of the solution of the equation from a left-side fractional derivative on space and with delay on time.
We consider the following initial-boundary value problem with functional delay of fractional order diffusion equation 1 < α 2 : is the prehistory function, τ > 0 is the value of delay, with initial conditions: and boundary conditions: The left-sides and the right-sides fractional derivatives are defined in Riemann-Liouville sense [17] We assume that the functions ϕ(x, t), c + (x), c − (x) and the functional f are such that this problem has unique solution u(x, t).
We denote by Q = Q[−τ, 0) the set of functions u(s) that are piecewise continuous on [−τ, 0) with a finite number of points of discontinuity of the first kind and right continuous at the points of discontinuity. We define the norm of functions by the relation We additionally assume that the functional f (x, t, u, v(·)) is given on [0, X] × [t 0 , T ] × R × Q and is Lipschitz in the last two arguments, i. e., there exists a constant L f such that, for all x ∈ [0, X], t ∈ [t 0 , T ], u 1 ∈ R, u 2 ∈ R, v 1 (·) ∈ Q and v 2 (·) ∈ Q the following inequality holds: We will use the approximations based on the shifted Grünwald definitions for the left-sides and right-sides fractional derivatives [22] respectively where N + , N − are positive integers, ∆x + = (x)/N + , ∆x − = (X − x)/N − and the normalized Grünwald weights are given by

Derivation of the fractional Crank-Nicholson method
Let h = X/N, introduce x i = ih, i = 0, . . . , N, and ∆ = (T − t 0 )/M, t j = t 0 + j∆, j = 0, . . . , M, we assume that τ /∆ = m is a positive integer. Denote by u i j approximations of functions u(x i , t j ) at the nodes.
Let us introduce a discrete prehistory for the time points t j , j = 0, . . . , M : will be called the interpolationextrapolation operator of discrete prehistory.
We will use piecewise-linear interpolation with extrapolation by continuation With the shifted Grünwald estimates for left and right fractional derivatives on space, applying usual approximation to a derivative on time and by the use of interpolation which is designed by extrapolation to the prehistory function, the discretized (0.1) at the nodes (x i , t j+1/2 ) will give us the analog of Crank-Nicholson method with the initial and boundary conditions then (1.1) we can be written as follows , we receive the Crank-Nicholson scheme for the variable coefficient heat conductivity equation with functional delay, because

Stability and convergence of the method
Denote by ε i j = u(x i , t j ) − u i j the error of the method at the nodes. We will say that the method converges with order h p + ∆ q if there exists a constant C independent of h and ∆ such that |ε i j | ≤ C(h p + ∆ q ) for all i = 0, 1, · · · , N and j = 0, 1, · · · , M. The residual (without interpolation) of the method (1.2) is, by definition, the grid function

Lemma 1. Suppose that the exact solution u(x, t) is twice continuously differentiable in t and
four times continuously differentiable in x and also the fractional derivatives Following [22], we find that We will use piecewise-linear interpolation with extrapolation by continuation of exact solution The residual (with interpolation) of methods (1.2) can be defined using the following grid function Lemma 2. Assume that the conditions of Lemma 1 are satisfied, then the residual of the method with piecewise-linear interpolation and extrapolation by continuation has order ∆ 2 + h.
P r o o f. Since the piecewise linear interpolation operator with extrapolation by continuation has second order [8], therefore there exists a constant C, such that for all t The residual with interpolation (2.5) and the residual without interpolation (2.1) are connected with the ratio From (0.2), (2.6), (2.7) and from Lemma 1 we get that the proof is completed.
We can copy the method (1.2) as follows (2.8) In order to reduce it to the general scheme [8], introduce y j = (u 1 j , u 2 j , · · · , u N −1 j ) ∈ Y, such that Y is a vector space of dimension N − 1 with norm Then (2.8) can be rewritten as where elements of a matrix A with dimension N − 1 × N − 1 have the form Lemma 3. [18,25] The coefficients g α,s satisfy the following properties

P r o o f. From the values of the elements of matrix A and from Lemma 3, we conclude that
This inequality means that the coefficient matrix E + A is a nonsingular (reversible), strictly diagonally dominant matrix.
We will note that if α is not integer number, then inequality (2.10) is strict and A is strictly diagonally dominant matrix.
As the matrix of E + A is reversible, equation (2.8) can be copied as follows , I is the operator of piecewise-linear interpolation with extrapolation by continuation.
From this lemma, it follows there is such constantŜ, that S n ≤Ŝ (2.12) for any natural degree n. This means that the method (2.11) is (unconditionally) stable. Now, the theorem of a convergence order [8] which is the main result of the general theory of numerical methods of the solution of evolutionary equations with heredity can be applied to the scheme (2.11). As the resulted scheme (2.11) has some different features in comparison with the scheme in [8], in particular, another determination of stability (2.12) is used in comparison with [8] then proof of the following theorem is needed.
Theorem 1. Assume that the condition in Lemma 3 is satisfied, the function Φ satisfies Lipschitz condition in the second argument, interpolation operator I satisfies Lipschitz condition, the residual with interpolation has order ∆ q +h p . Then, the method (2.11) converges with order ∆ q +h p .
P r o o f. Introduce the following vector of exact values z n = (u(x 1 , t n ), u(x 2 , t n ), · · · , u(x N −1 , t n )) and note that γ n = z n − y n , n = −m, ..., M . Then for n = 0, ..., M − 1 we have γ n+1 = Sγ n + ∆γ n + ∆ν n , is the vector residual with interpolation ν n = (ν 1 n , ν 2 n , · · · , ν N −n n ). The assumptions that mapping Φ and I are Lipschitz implies that (2.14) where L = L Φ L I , L Φ is Lipschitz constant of Φ on second argument, L I is Lipschitz constant of operator I. We notice that in the case of the described method L Φ = (E − A) −1 L f , The piecewise linear interpolation operator with extrapolation by continuation has Lipschitz constant L I = 2. Also, we note that γ 0 = 0, as starting values of a method are the exact values, unlike the general scheme [8].
It follows from (2.13) that From (2.15), (2.14) and from the definition of stability, we have (2.16) We use the notation R = max Induction base. If we set n = 0, in (2.18), then Induction step. Let the estimate (2.19) be valid for all indices from 1 to n. We need to show that the estimate is also valid for n + 1. Fix j n. Let i 0 = i 0 (j) be an index for which max j−m i j { δ i } is obtained. By induction assumption, we gain Thus, the following estimate is also valid Using the previous inequality and (2.18), we have

Numerical example
Let us consider the following two sided space fractional equation with varying delay in the time such that 1/2 ≤ x ≤ 3/2, 1 ≤ t ≤ 5, α is a constant, such that 1 < α < 2, with initial and boundary conditions on the form The exact solution of (3.1) is u( Define the maximum error by

Conclusion
A fractional analog of Crank Nicholson method is constructed to deal with the one dimensional space two sided space fractional diffusion equation with functional delay. The method is based on the idea of separating the current state and the prehistory function. To deal with the prehistory function, we introduced a discrete prehistory for the time points and the piecewise-linear interpolation with extrapolation by continuation operator is used. Unconditional stability of the proposed scheme is achieved. A theorem is obtained on the order of convergence of the method, which used the technique of proving similar statements for functional differential equations and methods from the general theory of difference schemes. A numerical example with time varying delay is introduced to test the agreement between theoretical and numerical results.