Abstract

A numerical solution for neutral delay fractional order partial differential equations involving the Caputo fractional derivative is constructed. In line with this goal, the drift term and the time Caputo fractional derivative are discretized by a finite difference approximation. The energy method is used to investigate the rate of convergence and unconditional stability of the temporal discretization. The interpolation of moving Kriging technique is then used to approximate the space derivative, yielding a meshless numerical formulation. We conclude with some numerical experiments that validate the theoretical findings.

1. Introduction

Partial differential equations (PDEs) with time delay play an important role in the mathematical modeling of complex phenomena and processes whose states depend not only on a given moment in time but also on one or more previous moments. We can mention a simple scenario involving the hemodynamic behavior of a person suffering from low or high glucose decompensation. This person can then be given intravenous insulin to compensate for the low level. Because the drug must be introduced into the bloodstream for it to take effect, the preceding scenario can be interpreted as a delay problem. As a result, there is a growing interest in studying biological and physical models with delay. The solutions of delay PDEs may represent voltage, concentrations, temperature, or various particle densities such as bacteria, cells, animals, and chemicals [13].

Delay PDEs with fractional derivatives have recently been studied using various numerical and analytical techniques such as [48]. It was pointed out in [9] that the derivatives of the dependent variable in the neutral type delay differential equations are both with and without time delays. Delay differential equations of neutral type appear in a variety of new phenomena, and its theory is even more complicated than the theory of nonneutral delay differential equations. From both a theoretical and practical standpoint, the oscillatory behavior of neutral system solutions is important. For some applications, such as the population growth, motion of radiating electrons, and the spread of epidemics in networks with lossless transmission lines, we refer the interested reader to [914].

A consideration of the following fractional PDE with a constant delay is the goal of this paper. For that end, we introduce with initial and boundary conditions where is the Caputo fractional derivative which is defined by

A novel interpolating element-free Galerkin approach to approximate the solution of the two-dimensional elastoplasticity problems was constructed in [15] using the interpolating moving least squares scheme for obtaining the shape function. Moreover, an improved element-free Galerkin scheme to solve nonlinear elastic large deformation problems was considered in [16]. The interpolating moving least squares approach using a nonsingular weight function is employed in [17] to approximate the solution of the problem of inhomogeneous swelling of polymer gels, and also the penalty scheme is used to enforce the displacement boundary condition; thus, an improved element-free Galerkin approach was constructed.

The interpolating element free Galerkin method has been developed to solve a variety of problems, including two-dimensional elastoplasticity problems [15, 18], two-dimensional potential problems [19], two- and three-dimensional Stokes flow problems [20], two-dimensional large deformation problems [21], incompressible Navier-Stokes equation [22], steady heat conduction problems [23], two-dimensional transient heat conduction problems [24], three-dimensional wave equations [25], two-dimensional Schrödinger equation [26], two-dimensional large deformation of inhomogeneous swelling of gels [27], biological populations [28], two-dimensional elastodynamics problems [29], and two-dimensional unsteady state heat conduction problems [30]. The theoretical analysis for the complex moving least squares approximation, the properties of its shape function, and its stability was analyzed in [31]. In [32], a variational multiscale interpolating element-free Galerkin scheme was established for solving the Darcy flow. For the numerical solution of generalised Oseen problems, a novel variational multiscale interpolating element-free Galerkin scheme was developed in [33] based on moving Kriging interpolation for obtaining shape functions using the Kronecker delta function. Zaky and Hendy [34] constructed a finite difference/Galerkin spectral approach for solving the Higgs boson equation in the de Sitter spacetime universe, which can inherit the discrete energy dissipation property. A high-order efficient difference/Galerkin spectral approach was proposed in [35] for solving the time-space fractional Ginzburg-Landau equation. Hendy and Zaky [36] proposed a finite difference/spectral method based on the formula on nonuniform meshes for time stepping and the Legendre-Galerkin spectral approach for solving a coupled system of nonlinear multiterm time-space fractional diffusion equations.

This paper is built up as follows. In Section 2, the temporal discretization is discussed. The analysis of the temporal discretization scheme is constructed in Section 3. The moving Kriging technique and its implementation are demonstrated in Section 4. Finally, numerical experiments are presented in Section 5 to illustrate the analysis of the obtained scheme.

2. Temporal Discretization

Assume that such that is a positive integer. Take and . Also, to make being grid points, the time-variable step size should be surrounded by instead of for . Thus, for . Here, we present a time-discrete scheme for Equation (1). For any function , we set

Lemma 1 (see [37]). Assume and . Then in which

Let be the exact solution of (1) and where . Thus, Equation (1) at can be rewritten as

Making use of Taylor expansion yields

Employing Lemma 1 and putting give

Furthermore, there is a constant that

Substituting the above result into (10) arrives at in which there exists such that

Removing yields

In the current paper, is an approximation of exact solution .

3. Analysis of the Temporal Discretization

In the current section, we check the stability of the numerical procedure.

Lemma 2 (see [38]). Let be a nonnegative sequence, and the sequence satisfies Then, for and , we have

Lemma 3 (see [37]). For any and , we obtain

Theorem 4. Let ; then scheme (15) is unconditionally stable.

Proof. We define . Now, we have Multiplying relation (19) by , integrating over and then summing from to give Recalling the left hand side of the above relation, invoking Schwartz inequality and Lemma 3 yields Moreover, for the first term in the right hand side of Equation (20), we have On the other hand, according to some simple mathematical actions, we have Also, for the delay term, we arrive at Replacing the above relations in Equation (20) yields or Now, Equation (26) can be simplified as Changing index from to arrives at The use of Equation (29) and Lemma 2 yields Thus, there exists that

4. Moving Kriging Interpolation and Its Implementation

Following [39, 40], we will invoke the technique of moving Kriging. Up to our knowledge and armed by the fact of the advantage of less CPU time consuming needed for the element free Galerkin approach based on the shape functions of moving Kriging than what needed for the element free Galerkin approach based on the shape functions of moving least squares approximation. In the meantime, the shape functions of moving Kriging interpolation can be deduced, which is analogous to moving least squares approximation over subdomain . Let is the approximate solution of on . The local approximation is formulated for any subdomain as such that and are monomial basis functions and monomial coefficients, respectively. Also, be the realization of a stochastic process. The covariance matrix of is given as in which (i) is the correlation matrix(ii) is the correlation function between any pair of nodes located at and

The correlation function is defined as [39, 40]

such that is a value of the correlation parameter [39, 40]. Using the best linear unbiased (BLUP) [39], we can write Equation (31) as follows [39, 40] in which

We will introduce some notations. The vector of known functions can be written as follows [39, 40]

and the matrix of defined function values at the nodes has the following representation [39, 40]

The correlation matrix is given as [39, 40]

The correlation vector at the nodes has the following form

The matrices and are given as

where is the identity matrix. Accordingly, Equation (34) can be written as follows [39, 40] or where the moving Kriging approach’s shape functions are as follows [39, 40]:

Now, we are ready to implement that kind of interpolation to the problem under consideration. Let the approximation solution of this equation be

in which are shape functions of moving Kriging approximation. Substituting Equation (44) in relation (15) gives

By collocating a set of arbitrary distributed nodes in the computational domain concludes

After doing some simplifications, we have

where . Now, the above formulation yields the following system of equations in which

5. Numerical Verification

In the current section, we investigate the convergence, capability, and stability of the developed numerical procedure. Also, the computational rate is calculated by

We consider the following problem in which the initial conditions are and also with no-flux boundary condition. The exact solution is .

Table 1 shows the results obtained based on the 500 collocation points, , , , and different values of . Table 1 confirms that the theoretical order (TO) in temporal direction is near to the computational order, i.e., . Figure 1 demonstrates the approximate solutions (a) and its absolute errors (b) on square domains (top figures), (middle figures), and (bottom figures) with , , and also 1000 collocation points. Figure 2 illustrates the approximate solutions (a) and its absolute errors (b) on irregular domains where (top figures), (middle figures), and (bottom figures) with , , and also 1000 collocation points. Figure 3 presents the approximate solutions (a) and its absolute errors (b) on irregular domains with , , and also 1000 collocation points. Figure 4 presents the approximate solutions (a) and its absolute errors (b) on irregular domains with , , and also 1000 collocation points.

6. Conclusion

The current paper presented a new numerical procedure for solving fractional damped diffusion-wave equations with delay. In this process, the time derivative is discretized by a finite difference scheme, and we constructed a time-discrete scheme. The stability and convergence of the proposed numerical formulation are studied, analytically and numerically. Then, the moving Kriging interpolation technique, as a meshless method, is used to get a fully discrete scheme. The proposed numerical method is flexible to simulate a wide range of PDEs including delay PDEs on irregular computational domains. Finally, an example is provided to demonstrate the stability and convergence of the new technique.

Data Availability

The data used to support the findings of this study are included within the article.

Conflicts of Interest

The authors declare that they have no competing interests.

Authors’ Contributions

M. Abbaszadeh contributed to the conceptualization, writing—original draft, data curation, figure preparation, and methodology. M. Dehghan contributed to the supervision and formal analysis and provided expert advice and validation. M. A. Zaky contributed to the writing—review and editing and carried out the experiments and formal analysis. A.S. Hendy contributed to the writing—review and editing and carried out the experiments and formal analysis. All coauthors contributed to editing the manuscript. All authors contributed significantly in writing this article. All authors read and approved the final manuscript.