Singular perturbation problems (SPPs) were first established as a research domain in the early 1990's [1] with the development of the boundary-layer idea in viscous flow [2] and has flourished over the last few years. Despite the large amount of studies that have already been done in this thematic area, more relevant and timely research is still ongoing.

Differential equations whose highest order derivative terms are attached with small positive number(s) are called singularly perturbed problems (SPPs). Singularly perturbation problems appearing with two small parameters are said to be two-parameter singularly perturbed problems. A singularly perturbed delay differential equation (SPDDE) is a differential equation in which the highest derivative is multiplied by a small parameter and containing at least one delay term either at the space variable, time variable, or both.

A lot of real-life physical problems are represented by linear or nonlinear differential models or by SPPs whose solution depends on the magnitude of the perturbation parameter. Singularly perturbed problems (SPPs) occur in the modeling of fluid dynamics, elasticity theory, quantum mechanics, reaction diffusion process, chemical reactor theory, plasma dynamics, meteorology, diffraction theory, aerodynamics, modeling of semi-conductors, hydrodynamics, and in several other applied fields [3–6].

Two-parameter singularly perturbed parabolic differential equations with time delay have many applications in different fields, for example, in engineering such as drift diffusion equation of semi conductor modeling [7] and chemical reactor model [8] in fluid dynamics [9].

Friedrichs and Wasow [10] were the first to use the term singular perturbation problems in their seminar at New York. In such problems, there are often narrow transition regions called boundary layers. In these regions, the solution changes rapidly or jumps abruptly and behaves regularly and slowly away from the layers.

For the solution of singular perturbation problems, one may apply the numerical approach or the asymptotic approach. The asymptotic approach provides the qualitative behavior of the problem and gives only a semi-quantitative information. However, the numerical approach provides quantitative information.

To solve singularly perturbed problems numerically (when analytical solutions are not available or more complicated), one can use finite difference methods, finite elements methods, spline approximation methods, and others, but, unless very fine grids are used, standard finite difference methods can not resolve the layers(s) and may not provide a uniformly convergent solution throughout the given domain.

The two non classical finite difference methods (FDMs) used for solving SPPs are fitted mesh methods (FMFDMs) and fitted operator methods (FOFDMs). In this article, we develop a uniformly convergent and accurate non-standard fitted operator finite difference method (NSFOFDM) based on the methodology of Mickens [11].

As the parameters ε and μ in the problem (1) of section (2) tend to zero, the solution will produce boundary layer(s) at x = 0 and x = 1. When μ = 1, problem (1) is convection- diffusion problem [12–14], and in this case, a boundary layer(s) of width O(ε) will occur around the edge x = 0. Again, when μ = 0, we have a parabolic reaction-diffusion problem [15] and thin boundary layers of width O(ε) appear near x = 0 and x = 1. O'Malley [16] introduced singularly perturbed two-parameter problems and examined asymptotic expansion for their solutions.

O'Malley [16, 17] identified that the nature of these problems is quite affected by the choice of ratio of μ² to ε. O'Malley et al. developed numerical methods to improve the accuracy of the asymptotic methods [16]. The class of time-dependent SPPs of convection-diffusion types with two parameter were studied in Munyakazi [18] using the classical finite difference method. Recently, the numerical solution of second-order two-parametric singularly perturbed ordinary differential equations (ODEs) with smooth data [19–32] and non-smooth data [33, 34] were considered.

Some uniformly convergent numerical methods for singularly perturbed time dependent delay differential equations have been developed in Bashier and Patidar, Kaushik et al., Kumar and Kumar, Erdogan and Cen, Cen, Singh et al., Ansari et al., and Kumar and Kumar [35–42].

In Govindarao et al. [43], a first-order uniformly convergent method was developed for two-parameter time dependent SPPs using an upwind finite difference scheme on Shishkin type meshes.

Solving two-parameter SPPs analytically is either more difficult task or the analytical solution does not exist. This is because of the small parameters attached to the highest order terms of the given problem. These attached small parameters exhibit a layer behavior in the solutions. The classical finite difference methods give unstable solution in the layer region. Moreover, the convergence and stability of the solution in numerical part varies according to the small parameters. From the existing literature we have seen, developing a parameter uniformly convergent numerical method for two-parameter singularly perturbed problems is still a challenging task.

The objective of this study is to analyze the solution when the delay is non-zero and the effect of the delay on the boundary layer solution, as well as investigate problems (1)–(2) with smooth data. We are inspired to develop a parameter uniformly convergent numerical scheme to treat a class of second-order two-parameter singularly perturbed time dependent problem (1)–(2). A non-standard fitted operator finite difference method based on the Crank-Nicolson discretization for time variable comprising a non-standard fitted operator finite difference on uniform mesh for spatial variable. The developed scheme is of second order in time and first order in space but has been improved to second order in both variables by using temporal mesh refinement in Section (5) in Tables 4, 5. Moreover, the comparison of the developed scheme with the existing scheme in Kumar and Kumar [44] is investigated in Section (5) in Table 6. The comparison shows that the maximum point-wise error of our scheme is less than the scheme in Kumar and Kumar [44].

This article is organized as follows. We first discuss the qualitative properties such as the bounds of the analytical solution u(x, t) of problem (1–2) and its derivative bounds in Section (2). The numerical scheme of the continuous problem is presented in Section (3). In this section, we also discuss the time discretization, the space discretization, the continuous problem discretization, and bounds of the discrete solution. The stability and convergence analysis of the scheme is presented in Section (4). In Section (5), we provide numerical example to show uniformly convergence of solution and its accuracy. We present the result and conclusions in Section (6).

We consider the following two families of two-parameter singularly perturbed time-delay problem. Our domain D̄ = D∪∂D, where D = (0, 1) × (0, T] and ∂D = L_l ∪ L_d ∪ L_r with L_d = [0, 1] × [−γ, 0](delay interval), L_l = {0} × (0, T](left side boundary) and L_r = {1} × (0, T](right side boundary). The governing equation is as follows:

where Lu(x, t) = εu_xx(x, t) + μa(x, t)u_x(x, t) − b(x, t)u(x, t), 0 < ε ≤ 1 and 0 ≤ μ ≤ 1 are perturbation parameters and γ is a delay parameter. In problem (1–2), we suppose that a(x, t), b(x, t), c(x, t), f(x, t), Φ_l(t), Φ_r(t), and Φ_d(x, t) for (x,t)∈D¯ are sufficiently smooth functions such that a(x, t) ≥ α > 0, b(x, t) ≥ β > 0, and c(x, t) ≥ ϒ > 0, independent of the perturbation parameters. At the corners, the regularity and compatibility conditions are

for D = (0, 1) × (0, T], and so that Φ_d (x, t) (initial-boundary data) satisfies appropriate compatibility criteria at the two corners, (0, 0) and (1, 0). Based on the above assumptions, the given problem in (1) possesses a unique solution in the considered domain.

Here, we develop the numerical scheme by discretizing the temporal domain, the spatial domain, and the given singularly perturbed problem in (1)–(2).

In this section, we investigate the stability and uniform convergence of the developed scheme.

Lemma 4.1. Discrete minimum principle.

Assume that [LN,MU]ij+1 is the discrete operator given in (10) and φij+1 is any mesh function satisfying φij+1|∂DN,M≥0 and [LN,Mφ]ij+1|DN,M≤0 for 0 ≤ i ≤ N, 0 ≤ j ≤ M and then φij+1|D̄N,M≥0.

Proof. Let s and l be indices such that φsl+1=min(i,j)φij+1 for φij+1∈D̄N,M. Again, assume that φsl+1<0. It is clear to see that (s, l) ∉ {0, N} × {0, M} because φsl+1≥0. It follows that φs+1l+1-φsl+1>0 and φsl+1-φs-1l+1<0.

which is a contradiction to the fact that LN,Mφsl+1≤0. Therefore, φsl+1≥0. The indices s and l being arbitrary, we obtain φij+1≥0 in D̄N,M.     □     □

The immediate consequence of the above lemma is the following lemma which is about a uniform stability estimate.

Lemma 4.2. Uniform stability estimate.

At any time level t_j, if Zij+1 is any mesh function such that Z0j+1=ZNj+1=0, then

Proof. To prove this lemma, we use the concept of barrier functions (φ±)ij and the above discrete minimum principle. Therefore, we define the two barrier functions as

where

⇒LN,K(φ±)ij+1≤0. Therefore, by lemma (4.1) above, we obtain

(φ±)ij+1≥0.     □

This lemma shows the uniform stability of the operator L^{N, M}. In the following two lemmas, we analyze the convergence of scheme in (8) using bounds of the truncation error in both variables.

Theorem 4.1. Error estimate in the spatial discretization.Let Uj+1(xi) and Uij+1 are the solution Equations 5, 8, respectively. If N and C are mesh number and sufficiently large constant, then the following error bound holds.

Proof. To prove this theorem, we use the differential and difference equation, and then define the error as follows:

Now, the Taylor's series expansion of Ui+1j+1, Ui-1j+1 and 1ϕi2 are

Using these substitutions in (12) and applying some simplification gives