Abstract

Numerical computation for the class of singularly perturbed delay parabolic reaction diffusion equations with integral boundary condition has been considered. A parameter-uniform numerical method is constructed via the nonstandard finite difference method for the spatial direction, and the backward Euler method for the resulting system of initial value problems in temporal direction is used. The integral boundary condition is treated using numerical integration techniques. Maximum absolute errors and the rate of convergence for different values of perturbation parameter and mesh sizes are tabulated for two model examples. The proposed method is shown to be parameter-uniformly convergent.

1. Introduction

Singularly perturbed delay differential equations model physical problems for which the evaluation depends on the present state of the system and also on the past history. In general, when the perturbation parameter approaches zero, the smoothness of the solution of the singularly perturbed delay differential equations deteriorates and it forms a boundary layer(s) [1]. The time-dependent problem, i.e., singularly perturbed parabolic delay differential equations, has applications in the study of modeling of physiological kinetics, production of blood cell, and neuronal variability (see [24] and the references therein). But, for more applications, the theory and numerical treatment of delay differential equations can be found in [5].

It is known that standard numerical methods for solving singular perturbation problems are unstable and fail to give accurate results when the perturbation parameter is small. The accuracy and convergence of the methods need attention, because the treatment of singular perturbation problems is not trivial, and the solution depends on perturbation parameter and mesh size [6]. This suggests that the numerical treatment should be improved. The presence of the singular perturbation parameter leads to occurrences of oscillations or divergence in the computed solutions while using standard numerical methods [7]. To overcome these oscillations or divergence, a large number of mesh points are required as goes to zero. This is difficult and sometimes impossible to handle such cases. Therefore, it is necessary to develop suitable numerical methods which are uniformly convergent to solve this type of differential equations.

Many researchers have proposed different numerical methods to solve singularly perturbed time delay parabolic differential equations; for instance, in [8], Ansari et al. proposed a parameter robust finite difference method for singularly perturbed delay parabolic partial differential equations of reaction diffusion type with Dirichlet boundary conditions on Shishkin mesh. They made use of the results of Miller et al. [9]. In [10], Bashier and Patidar improved the order of accuracy for the method suggested in [8]. Avudai Selvi and Ramanujam [11] studied singularly perturbed time delay parabolic reaction diffusion equation subject to mixed-type boundary conditions. A second-order fitted operator finite difference method for reaction diffusion type equation was suggested by Bashier and Patidar [12]. Singh et al. [13] developed an overlapping Schwarz domain decomposition method for parabolic reaction diffusion problems with time delay. An extended cubic B-spline collocation method on Shishkin mesh has been developed for singularly perturbed delay parabolic reaction diffusion problem with mixed-type boundary conditions in [14]. A numerical method is developed for solving a singularly perturbed parabolic delay differential equation whose solution exhibits a boundary layer behaviour in [15] and singularly perturbed parabolic differential difference equations arising in computational neuroscience in [16]. The authors in [17] gave a parameter robust numerical scheme to the solution of second-order singularly perturbed parabolic partial differential equations with large delay by using the finite difference scheme on the piecewise-uniform Shishkin mesh. In [18], Kumar and Kumari constructed a parameter-uniform implicit scheme for a class of parabolic singularly perturbed reaction diffusion initial-boundary value problems with large delay in the spatial direction, where the solution exhibits twin boundary layers and an interior layer on a piecewise-uniform mesh.

But, in recent years, there has been a growing interest in numerical methods on singularly perturbed delay differential equations with integral boundary conditions. In [19, 20], the authors consider a first-order singularly perturbed quasilinear boundary value problem with integral boundary conditions on Shishkin mesh. The authors in [21, 22] studied singularly perturbed boundary value problem for second-order ordinary differential equations with nonlocal boundary condition on uniform mesh. In [23], Sekar and Tamilselvan have discussed singularly perturbed delay differential equations of convection-diffusion type with integral boundary condition, using the finite difference method on Shishkin mesh. The method has almost first-order convergence with respect to . Debela and Duressa [24] proposed first-order numerical methods using the exponentially fitted finite difference method on a uniform mesh. They also proposed the accelerated fitted finite difference method for singularly perturbed delay differential equations with nonlocal boundary condition in [25]. In [26], a uniformly convergent numerical method is constructed and a numerical integration method is used to treat the nonlocal boundary condition. However, up to the best of our knowledge, except the work in [27], not much work has been done to solve the problem under consideration.

In the present paper, motivated by the studies of [27], we construct and analyze an accurate and uniformly convergent numerical scheme for solving singularly perturbed delay parabolic reaction diffusion equations with integral boundary condition. The proposed scheme uses the procedures of method of line (i.e., first discretizing in spatial direction followed by discretization in temporal direction) using the backward Euler method in temporal direction and nonstandard fitted operator FDM on spatial direction.

The content of this paper is as follows: In Section 2, definition of the problem and the properties of continuous solution are given. In Section 3, discretizing the spatial domain and techniques of nonstandard finite difference are discussed, and the uniform convergence of the discrete problem is proved. Next, the backward Euler method used for the system of initial value problems resulted from spatial discretization and the convergence of the discrete scheme are discussed. In Section 4, numerical results and discussion are given to validate the theoretical analysis, and finally, in Section 5, the conclusion of the work done is presented.

Notation 1. The notation for jump in any function is at a point ‘’ with and C is a positive constant independent of and . The norm for a given function defined on the domain is calculated by .

2. Statement of the Problem

Let us consider the following singularly perturbed delay parabolic differential equation with integral boundary condition on , where is some fixed positive time and , where and are the left and the right sides of the rectangular domain corresponding to and , respectively, and is the base of the domain and given by .

Subject to initial conditionsand boundary conditionswhere and is a given constant. and are sufficiently smooth, bounded functions that satisfy

Furthermore, is a nonnegative, monotonic function and satisfies .

Problems (1)–(3) can be rewritten aswherewith boundary conditionswhere .

2.1. Properties of Continuous Problems

In the study of the numerical aspects of singularly perturbed problems, their analytical aspects play an important role. Setting the value , the reduced problem corresponding to (6)–(8) is

As need not satisfy and , the solution exhibits boundary layers at . Furthermore, as need not be equal to , the solution exhibits interior layers at . Problems (1)–(3) satisfy the following lemma.

Lemma 1 (continuous maximum principle). Let such that , , and . Then, .

Proof. Define a test functionNote that , and . LetThen, there exists such that and , . Therefore, the function attains its minimum at . Suppose the theorem does not hold true, then :Case (i): . It is a contradiction.Case (ii): . It is a contradiction.Case (iii): . It is a contradiction.Case (iv): . It is a contradiction.Case (v): . It is a contradiction.Hence, the proof of the theorem is verified.

Lemma 2 (stability result). The solution of problems (6)–(8) satisfies the bound

Proof. It can be easily proved using the maximum principle (Lemma 1) and the barrier functionswhere and are the test functions as in Lemma 1.
The existence and uniqueness of a solution for problems (6)–(8) can be established by assuming that the data are continuous and imposing appropriate compatibility conditions at the corner points (0, 0), (2, 0), (−1, 0), and (1, 0) (see [28]). Then, the required compatibility conditions areNote that, are assumed to be smooth for (15) to make sense that .
The following theorem gives sufficient conditions for the existence of a unique solution of problems (6)–(8) (see [28]).

Theorem 1. Let Assume that compatibility conditions (14) and (15) are fulfilled. Then, problems (6)–(8) have a unique solution and .

Proof. For the proof, refer to [28].

Theorem 2. Let . Assume that compatibility conditions (14) and (15) are fulfilled. Then, problems (6)–(8) have a unique solution and . Furthermore, the derivatives of the solution satisfywhere the constant C is independent of .

Proof. For the proof, refer to [27].
The nonclassical bounds in singular and regular components and their derivatives are established in the following theorem since the bounds in Theorem 2 are not adequate for the proof of parameter-uniform error estimate.

Theorem 3. Let the data , and . Assume that compatibility conditions (14) and (15) are satisfied. Then, we havewhere C is the constant independent of .

Proof. One may refer [27] for the details.

Theorem 4. The partial derivative of satisfies such that .

Proof. The proof of the theorem is completed by using the estimates of (18) and (19) and the decomposition .

3. Description of the Numerical Method

In this section, the spatial and temporal semidiscretization was described and also their convergence was analyzed separately.

3.1. Spatial Semidiscretization

The theoretical basis of the nonstandard discrete modeling method is based on the development of the exact finite difference method. In [29], Mickens presented techniques and rules for developing nonstandard FDMs for different problem types. In Mickens rules, to develop a discrete scheme, a denominator function for the discrete derivatives must be expressed in terms of more complicated functions of step sizes than those used in the standard procedures. This complicated function constitutes a general property of the schemes, which is useful while designing reliable schemes for such problems. In line to this, we need to derive a nonstandard finite difference method which captures the layer behaviour of the problem on a uniform mesh. Mickens in [29, 30] gave a novel approach of nonstandard finite difference method, and the basic idea of this method is to replace the denominator function of the second-order derivative with a suitable function in the differential equation, which comes under the category of the fitted operator finite difference method. On the spatial domain [0, 2], we introduce the equidistant meshes with uniform mesh length such thatwhere is the step size and is the number of mesh points in the spatial direction.

Moreover, we consider the domain which is discretized into equal number of subintervals, each of length . For the problem of the form in (6)–(8), in order to construct the exact finite difference scheme, we follow the procedures used in [31]. Consider the constant coefficient subequations given in (22) by ignoring the time variable aswhere . Thus, problem (22) has two independent solutions, namely, , with

We denote the approximate solution of at by Now, our objective is to calculate a difference equation which has the same general solution as problem (22) has at the grid point given by . Using the procedures used in [31], we have

Simplifying (24), we obtain thatis an exact difference scheme for (22).

After doing the arithmetic manipulation and rearrangement on (25), we obtain

The denominator function becomes . Adopting this denominator function for the variable coefficient problem, we write it asCase 1: consider (6)–(8) on the domain .Let denote the approximation of . By using the nonstandard finite difference approximation, at this stage, the problem in (6)–(8) reduces to the semidiscrete form aswhereCase 2: since the problem involves nonlocal boundary conditions, we considered the following to obtain an equation at the end boundary condition.

For , the composite Simpson’s integration rule approximates the nonlocal condition in (1) as follows:where ,

Now, using (31) and left boundary condition, we can obtain the following formula for :

Therefore, on the whole domain , the basic schemes to solve (6)–(8) are the schemes given in (29), (30), and (32).

The time domain is continuous at this stage, and the system of IVPs in (29), (30), and (32) can be written in the compact form aswhere is a matrix of size and are vectors of size . The entries of the coefficient matrix are, respectively, given byand the corresponding right-hand-side vector in the equation system has the entries

Before we proceed with the convergence analysis, we highlight some properties of the discrete problem to original (1) in the form of lemma which plays a major role in the said analysis.

For ,

For ,

Subject to the boundary conditions,where

Problems (29), (30), and (32) satisfy the following well-known semidiscrete maximum principle on .

Lemma 3 (semidiscrete maximum principle). Assume thatand mesh function satisfies Then, , , and imply that .

Proof. Define a test function asNote thatLet . Then, there exists such that and . Therefore, the function attains its minimum at . Suppose the theorem does not hold true, then :Case (i): . It is a contradiction.Case (ii): . It is a contradiction.Case (iii): . It is a contradiction.Case (iv): . It is a contradiction.Case (v): . It is a contradiction.Hence, the theorem is proved.

Lemma 4 (semidiscrete stability result). Let be any mesh function; then,

Proof. Consider the barrier functionswhereand is the test function as in Lemma 3.
From (44), it is clear that , and andUsing Lemma 3, .

3.2. Convergence Analysis

Now, let us analyze the convergence of the spatial semidiscretization. We proved above the discrete problem satisfying the maximum principle and the uniform stability estimate. Note that is denoted for the spatial discretization approximate solution to the exact solution at .

Lemma 5. For a fixed mesh and for , we havewhere .

Proof. Consider the partition of the interval [0, 1], with points . By using , then we haveThen, applying L’Hospital’s rule repeatedly results inHence, the proof is completed.
Now, the truncation error of scheme (29) is given byTaylor series expansions of the terms are given as follows:Using the truncated Taylor series expansions of the terms yieldsNext, we use a truncated Taylor series expansion of the denominator function of order five:Now, substituting (53) into (52), we obtainUsing the bounds on the derivatives and Lemma 5 giveswhere we have used the relation The discrete problem satisfies the following bound:Using Lemma 5, we get the result:

Remark 1. A similar analysis for convergence may be carried out for the finite difference scheme on (30).
Next, the truncation error of scheme (31) at the point is given byNow, we take the norm of :where .
Applying Lemma 4, we have the bound .
Hence,

Theorem 4. The error associated with the spatial semidiscretization of problems (29), (30), and (32) satisfieswhere C is a constant independent of and .

3.3. Temporal Semidiscretization

IVPs (29), (30), and (32) are discretized by the backward Euler method on a uniform mesh. Now, we denote the approximation of by . We perform the time discretization as follows: let be a positive integer and , and , where denotes the number of mesh points in time direction; then,with initial condition . Rearranging equation (61) gives

Lemma 6. The local truncation error associated with the time integration satisfieswhere is a constant independent of and .

Proof. The local truncation error is defined asA Taylor series expansion of takes the formThe local truncation error by subtracting yieldsSince the matrix is invertible, using the relation for small and , we obtainwhich completes the proof.

Lemma 7. In this temporal direction, the global error estimates are given bywhere is the global error in the temporal direction at time level.

Proof. One may refer [15] for the details.
Since and are independent of the perturbation parameter , taking the supremum for all , we obtainThis shows that the discretization in time direction is consistent and global error is bounded, with the error bound .

3.4. Fully Discretized Problem

Now, using (60) and (69), the uniform convergence of the fully discrete scheme was proved as

Hence, we obtain the required bound as follows:

Thus, the inequality in (71) shows the parameter-uniform convergence of the proposed scheme with order of convergence, second order in spatial direction and first order in temporal direction.

4. Numerical Results and Discussion

To validate the established theoretical results, we perform numerical experiments using the proposed numerical scheme on the problem given in (6)–(8). We consider two numerical examples to verify the parameter-uniform convergence of the proposed scheme.

Exact solution is not available for these model problems; therefore, maximum nodal errors are calculated by using the double mesh technique aswhere N is the number of mesh points in the x direction. For any value of the mesh points N and , the uniform (parameter-uniform) error estimate is calculated by

The rate of convergence of the method is given by

The numerical results are presented for the value of the perturbation parameter:

Example 1. Consider the following singularly perturbed problem [27]:Subject to the initial and boundary conditions,

Example 2. Consider the following singularly perturbed problem:Subject to the initial and boundary conditions,In Tables 1 and 2 the numerical results of the proposed method are tabulated in terms of maximum absolute errors, numerical rate of convergence, and uniform errors. As one observes from tables, for each number of mesh intervals as , the maximum absolute error becomes stable and uniform. This shows that convergence of the proposed method is independent of . The uniform error and uniform rate of convergence of the method are observed from the last row of each table. In Tables 3 and 4, the performance of the proposed scheme is investigated by comparing it with recently published paper in [27]. As one can see, the proposed method gives more accurate results.
In Figures 14, we observe the behaviour of the numerical solution, and for very small values of , the problem under consideration exhibits strong parabolic boundary layer near ; furthermore, an interior layer at is created.

5. Conclusions

In this paper, a parameter-uniform numerical scheme has been developed for solving singularly perturbed delay parabolic differential equations with integral boundary condition exhibiting parabolic boundary layers and an interior layer for small values of the perturbation parameter. The behaviour of the analytical solution and its boundedness are discussed. The numerical scheme is developed based on the method of a line that constitutes the nonstandard finite difference for the spatial discretization, and a backward Euler method is used in the temporal direction for the system of initial value problems resulting from the spatial discretization. To treat the integral boundary condition, the numerical integration technique is applied. The stability and convergence of the proposed scheme are analyzed. Two model examples have been considered to validate the applicability of the method. The computational results are presented in terms of tables and figures. The proposed numerical scheme could be the product for new equations (see e.g. [3234]), is shown to be first-order convergent in time and second-order convergent in space, and gives more accurate results than the previously developed methods existing in the literature.

Data Availability

No data were used to support the study.

Conflicts of Interest

The authors declare no conflicts of interest.

Acknowledgments

The authors would like to acknowledge Jimma University, College of Natural Sciences, for partial financial support.