## ORIGINAL RESEARCH article

Front. Appl. Math. Stat., 31 August 2023
Sec. Numerical Analysis and Scientific Computation
Volume 9 - 2023 | https://doi.org/10.3389/fams.2023.1244490

# Fitted computational method for singularly perturbed convection-diffusion equation with time delay

• 1Department of Applied Mathematics, Adama Science and Technology University, Adama, Ethiopia
• 2Department of Mathematics, Jimma University, Jimma, Ethiopia

A uniformly convergent numerical scheme is proposed to solve a singularly perturbed convection-diffusion problem with a large time delay. The diffusion term of the problem is multiplied by a perturbation parameter, ε. For a small ε, the problem exhibits a boundary layer, which makes it challenging to solve it analytically or using standard numerical methods. As a result, the backward Euler scheme is applied in the temporal direction. Non-symmetric finite difference schemes are applied for approximating the first-order derivative terms, and a higher-order finite difference method is applied for approximating the second-order derivative term. Furthermore, an exponential fitting factor is computed and induced in the difference scheme to handle the effect of the small parameter. Using the discrete maximum principle, the stability of the scheme is examined and analyzed. The developed scheme is parameter-uniform with a linear order of convergence in both space and time. To examine the accuracy of the method, two model examples are considered. Further, the boundary layer behavior of the solutions is given graphically.

## 1. Introduction

Delay differential equations (DDEs) are differential equations in which the evolution of the system is influenced by its past history. DDEs are called retarded types if the delay argument does not appear in the highest-order derivative term; otherwise, they are neutral types. DDEs play an important role in a variety of fields, including robotics, biosciences [1], economics, epidemiology and mechanics [2], fluid dynamics, reaction-diffusion equations [3], and population dynamics [4].

A singularly perturbed delay differential equation (SPDDE) is a delay differential equation in which its higher-order derivative term is multiplied by a small perturbation parameter (0 < ε ≪ 1) and contains at least one delay parameter on the term different from the highest derivative. In contrast to the magnitude of the delay parameter with the perturbation parameter, the delay is classified as a large delay or a small delay. If the magnitude of the delay parameter of the SPDDE is smaller than the perturbation parameter, then the equation is said to be a singularly perturbed delay differential equation with a small delay, whereas when the magnitude of the delay parameter is higher than the perturbation parameter, it is said to be a singularly perturbed delay differential equation with a large delay [5]. A singularly perturbed problem, which arises as a time delay, occurs in many application areas of science and engineering, for instance, in the simulation of oil extraction from underground reservoirs, chemical processes, fluid flows, water quality problems in river networks, and mechanical systems [6].

The presence of ε as a multiple of the higher-order derivative term causes a boundary layer. The boundary layer is an asymptotically narrow region located in the neighborhood of the endpoints of the domain, where the solution has a steep gradient as ε tends to zero [6]. With the rapidly changing behavior of the solution in the boundary layer, one encounters computational difficulties in treating a singularly perturbed problem using analytically or classical numerical schemes. On the contrary, classical numerical schemes lead to spurious non-physical oscillations in the numerical solution, unless an unacceptably large number of mesh points are considered, which leads to a massive computational cost [7]. In response to this, different authors have to look for sounding numerical schemes which converge uniformly regardless of ε.

Recently, the authors in [8], proposed the implicit-Euler scheme in the time direction and the central difference scheme in the space direction. The authors in [7, 9], proposed the implicit-Euler scheme in the time direction and the hybrid scheme by a proper combination of the midpoint upwind in the outer region and the central difference scheme in the inner region in the spatial direction on the Shishkin mesh. Moreover, this method is addressed in [10], for the two-parameter problem. In [11], the authors proposed the implicit-Euler scheme in the time direction and the hybrid scheme on a generalized Shishkin mesh in the spatial direction. Gowrisankar and Natisan in [12] developed the backward Euler scheme in time direction and the upwind finite difference scheme in the spatial direction using a piecewise uniform mesh. The implicit Euler scheme in the time direction and the upwind scheme in the spatial direction are considered in [13]. In [14], the implicit trapezoidal scheme in the time direction and the hybrid scheme by proper combination of the midpoint upwind in the outer region and the central difference in the inner region in the spatial direction are used.

The implicit Euler scheme in the time direction and the central difference scheme in the space direction are used in [4]. The extended cubic B-spline is considered in [15]. A domain decomposition method is considered in [16, 17]. The authors in [18, 19] proposed hybrid scheme on both Shishkin and Bakhvalov meshes. Podila and Kumar [20] proposed a new stable finite difference scheme on a uniform mesh and also on an adaptive mesh. The backward Euler scheme in the time direction and exponentially fitted difference method is considered in [21]. The Crank-Nicolson method in the time direction and a novel fitted finite difference scheme in spatial direction are proposed in [22]. The Crank-Nicolson method in the time direction and an exponentially fitted spline in the spatial direction are discussed in [23]. The implicit Euler scheme in the time direction and the non-standard finite difference method in the space direction are considered in [24]. In [25], the authors proposed Crank-Nicolson method in the time direction and the operator compact implicit (OCI) method on the Shishkin mesh in the space direction. The backward Euler in the time direction and method of line following Micken's type discretization for the space derivatives are used in [26]. Sahoo and Gupta [27] used higher-order difference with an identity expansion (HODIE) on a piecewise uniform mesh. A similar technique was also used in [28] for a coupled system of singularly perturbed problems. The authors [29, 30] proposed the numerical schemes that work for both cases when the delay term is large or small.

The main aim of this work was to develop a ε-uniform numerical scheme for the class of singularly perturbed convection-diffusion problem with a large time delay. The method comprises the backward Euler scheme in the time direction and an exponentially fitted higher-order finite difference scheme in the spatial direction. Error bound and uniform convergence of the developed scheme is investigated and proved. The proposed scheme gives more accurate, stable, and uniformly convergent results.

In this study, C has been considered as a generic positive constant, which does not depend on Δs, Δt, and ε. The maximum norm is denoted by ‖.‖, which is defined by ‖γ‖ = maxs,t∈Ω|γ(s, t)|.

## 2. Continuous problem

Let Ω = Ωs × Ωt = (0, 1) × (0, 𝕋] for 𝕋 > 0, we consider SPDDE of the form

where 𝔏εz(s, t) = −εzss(s, t) + β(s, t)zs(s, t) + α(s, t)z(s, t).

Here, ε ∈ (0, 1] and δ > 0 are the perturbation parameter and the delay parameter, respectively. We pretended that the functions β(s, t), α(s, t), κ(s, t), γ(s, t) on $\stackrel{̄}{\Omega }=\left[0,1\right]×\left[0,𝕋\right]$ and ψb(s, t), ψl(t), ψr(t) on η = ηl ∪ ηr ∪ ηb are smooth enough and bounded which meet α(s, t) ≥ ϖ > 0, κ(s, t) ≥ φ > 0, β(s, t) ≥ μ > 0 on $\stackrel{̄}{\Omega }$. These conditions assure that problem (1) has a boundary layer near s = 1 [7].

### 2.1. A priori bounds

Under the premises that the data are Hölder continuous and satisfy the following compatibility conditions at the corner points and the delay terms [31], we confirm the existence and uniqueness of the solution of (1)

These assumptions and conditions are fulfilled. Then, the problem (1) admits a unique solution [31].

Setting ε = 0, the reduced problem of (1) is given as

where z0(s, t) is the solution of the reduced problem.

Lemma 2.1. Let z(s, t) be the solution of (1). Then, we have

where C does not depends on ε.

Proof: The proof is considered in [7].

The operator $𝔏=\left(\frac{\partial }{\partial t}+{𝔏}_{\epsilon }\right)$ in (1) satisfies the next lemma.

Lemma 2.2. (Maximum principle). Let ν(s, t) ∈ C2(Ω) ∪ C0(η) satisfies ν(s, t) ≥ 0 (s, t) ∈ η. If 𝔏ν(s, t) ≥ 0, (s, t) ∈ Ω, then $\nu \left(s,t\right)\ge 0,\left(s,t\right)\in \stackrel{̄}{\Omega }$.

Proof: The proof is considered in [14].

Lemma 2.3. (Stability result). Let z(s, t) be the solution of (1). Then, we have

where ϖ ≤ α(s, t).

Proof: The proof is considered in [22].

Lemma 2.4. The derivative of the solution z(s, t) of (1) with respect to s and t satisfy

where μ ≤ β(s, t).

Proof: The proof is considered in [7].

## 3. Numerical scheme

### 3.1. Temporal semi-discretization

The time domain [0, 𝕋] is discretized uniformly with step size Δt as ${\Omega }_{t}^{M}=\left\{{t}_{m}=m\Delta t,m=0,1,2,...,M,{t}_{M}=𝕋,\Delta t=𝕋/M,\right\}$ and ${\Omega }_{t}^{j}=\left\{{t}_{m}=m\Delta t,m=0,1,2,...,j,{t}_{j}=\delta ,\Delta t=\delta /j,\right\}$ with M + 1 mesh points in [0, 𝕋] and j + 1 mesh points in [−δ, 0]. We have 𝕋 = for some positive integer r.

Applying the backward Euler scheme for time derivative, we get

Simplifying (8), we have

where ${𝔏}^{\Delta t}{Z}^{m}\left(s\right)=-\epsilon \frac{{d}^{2}{Z}^{m}\left(s\right)}{d{s}^{2}}+\beta \left(s,{t}_{m}\right)\frac{d{Z}^{m}\left(s\right)}{ds}+P\left(s,{t}_{m}\right){Z}^{m}\left(s\right)$ and $P\left(s,{t}_{m}\right)=\left(\frac{1}{\Delta t}+\alpha \left(s,{t}_{m}\right)\right)$ with the boundaries

Now, (9) rewrite as

where ${𝔏}^{\Delta t}Z\left(s\right)=-\epsilon \frac{{d}^{2}Z\left(s\right)}{d{s}^{2}}+\beta \left(s,{t}_{m}\right)\frac{dZ\left(s\right)}{ds}+P\left(s,{t}_{m}\right)Z\left(s\right)$ and and $Q\left(s\right)={Q}^{m}\left(s\right)=Q\left(s,{t}_{m}\right)$.

The local truncation error in the time direction is given as ${e}_{m}\left(s\right):=z\left(s,{t}_{m}\right)-{Z}^{m}\left(s\right),m=0\left(1\right)M.$

Lemma 3.1. The local error em at tm satisfies the bound

Lemma 3.2. The global error Em at tm satisfies the bound

Proof: Using Lemma 3.1, the global error Em bound at mth time step is given as

Lemma 3.3. For every m = 0(1)M − 1, the solution Zm(s) of (9)-(10) satisfies the estimate

Proof: From (11), $-\epsilon \frac{{d}^{2}Z\left(s\right)}{d{s}^{2}}+\beta \left(s,{t}_{m}\right)\frac{dZ\left(s\right)}{ds}=g$, where g = Q(s) − P(s, tm)Z(s).

Now, we integrate twice and we obtain

$Z(s)=ZP(s)+C1+C2∫s1exp(-ε-1(B(1)-B(y)))dy,$

where ${Z}_{P}\left(s\right)=-{\int }_{s}^{1}u\left(y\right)dy,u\left(s\right)={\int }_{s}^{1}{\epsilon }^{-1}g\left(y\right)exp\left(-{\epsilon }^{-1}\left(B\left(y\right)-B\left(s\right)\right)\right)dy,B\left(s\right)={\int }_{0}^{s}\beta \left(y\right)dy.$

Using inequality

$exp(-ε-1(B(y)-B(s)))≤exp(-ε-1μ(y-s)),s≤y,$

and the bound

Hence, |Zp(s)| ≤ C. Here, ${Z}^{\prime }\left(1\right)=-{C}_{2}$. The boundary condition Z(1) = 0 yields C1 = 0. Now, the constants C1 and C2 must satisfy

$C2∫01exp(-ε-1(B(1)-B(y)))dy=-Zp(0).$

Since B(s) is bounded on (0, 1), B(1) − B(y) ≤ C(1 − y). Then,

$∫01exp(-ε-1(B(1)-B(y)))dy≥∫01exp(-ε-1-μ(1-y))dy≥Cε.$

It follows that $|{C}_{2}|\le C{\epsilon }^{-1}$. Hence, $|{Z}^{\prime }\left(1\right)|=|{C}_{2}|\le C{\epsilon }^{-1}$. Finally,

$Z′(s)=u(s)-C2exp(-ε-1(B(1)-B(s)))$

implies that

The proof is done for i = 1. For i > 1 follows by induction and repeated differentiation. For the details, refer [32].

### 3.2. Spatial discretization

We discretize the spatial domain [0, 1] into N equal number of sub-intervals with the length of h as 0 = s0, s1, ..., sN = 1, and sn = nh, n = 0(1)N. Consider a smooth function Z(s) in the interval [0, 1]. From Taylor's series approximation, we get

Following a similar relation of (15), it holds

From (16), we have

Substituting (17) into (16) and simplifying, we obtain

where $𝔗=\frac{{h}^{4}}{20}{Z}_{n}^{\left(4\right)}-\frac{13{h}^{8}}{302,400}{Z}_{n}^{\left(8\right)}+O\left({h}^{10}\right).$

From (11), we draw

where

Using (19), we have

From the Taylor series approximations of ${Z}_{n-1}^{{}^{\prime }},{Z}_{n}^{{}^{\prime }}$ and ${Z}_{n+1}^{{}^{\prime }}$, we get

Substituting (21) into (20), we have

From (18), we draw

Substituting (22) into (23) and rearranging, we obtain

where

#### 3.2.1. Computing the exponential fitting factor

We introduce the exponential fitting factor σ to handle the effect of ε in the layer. From the singular perturbation theory stated in [33], the zero order asymptotic solution of the problem of the form

is given as

From Taylor's series, approximation for β(s) and α(s) restricting to their first terms about s = 1 is given as

where Z0(s) is the solution of reduced problem. Taking h → 0 and solving (28) at sn−1, sn, and sn+1, we get

where $\rho =\frac{h}{\epsilon }$. Multiplying (24) by h and the term containing ε by σ and evaluating the limit of the resulting equation as h → 0, we get

From (29), we have

Substituting (31) into (30) and simplifying yields

$σρ(eβ(1)ρ+e-β(1)ρ-2)=β(sn,tm)60(30eβ(1)ρ-30e-β(1)ρ).$

Then, we get the fitting factor σ

Therefore, the required scheme is taken as

where

In the explicit form, it becomes

where

### 3.3. Stability and uniform convergence analysis

Lemma 3.4. (Discrete maximum principle). Assume that ${Z}_{0}^{m}\ge 0,{Z}_{N}^{m}\ge 0$ and ${𝔏}^{\Delta t,h}{Z}_{n}^{m}\ge 0,\forall n=1\left(1\right)N-1$, then ${Z}_{n}^{m}\ge 0,\forall n=0\left(1\right)N$.

Proof: Assume that there is k ∈ {0, 1, 2, …, N}, such that ${Z}_{k}^{m}={min}_{0\le n\le N}{Z}_{n}^{m}<0$. Assume that ${Z}_{k}^{m}<0$ and from the assumption, it is shown that k ∉ {0, 1}. So, we have ${Z}_{k+1}^{m}-{Z}_{k}^{m}>0$ and ${Z}_{k}^{m}-{Z}_{k-1}^{m}<0$. Then, we get ${𝔏}^{\Delta t,h}{Z}_{k}^{m}<0$ for k = 1(1)N − 1. So, the assumption is wrong. Therefore, ${Z}_{n}^{m}\ge 0$ and ∀n = 0(1)N.

Lemma 3.5. (Uniform stability result). Let ${Z}_{n}^{m}$ be the solution of (33), then we have

where P(sn, tm) ≥ ζ > 0.

Proof: Let and define the barrier functions ${\vartheta }_{n,m}^{±}$ by ${\vartheta }_{n,m}^{±}=R±{Z}_{n}^{m}$. On the boundaries, we get ψl(tm) ≥ 0 and . For sn, n = 1(1)N − 1, we obtain

By Lemma 3.4, we get . Therefore, the needed bound is obtained.

From Taylor's series expansion, we get

where

The next theorem provides the truncation error estimate for the developed scheme.

Theorem 3.6. Let the coefficients α(s, tm), β(s, tm), and κ(s, tm) of (9)-(10) be sufficiently smooth such that Zm(s) ∈ C4[0, 1]. Then, the solution ${Z}_{n}^{m}$ of (33) satisfies the next bound

Proof: The error estimate in the spatial direction is given as

where $\sigma =\beta \left({s}_{n},{t}_{m}\right)\frac{\rho }{2}coth\left(\beta \left(1\right)\frac{\rho }{2}\right)$ and $\rho =\frac{h}{\epsilon }$.

For the constants C1 and C2, we have $|\rho coth\left(\rho \right)-1|\le {C}_{1}{\rho }^{2}$ for ρ ≤ 1. For ρ → ∞, since ${lim}_{\rho \to \infty }coth\left(\rho \right)=1$ which gives |ρ coth(ρ) − 1| ≤ C1ρ.

Generally, ∀ρ > 0, we express as

and we have

From the bounds in (38), (40), and (42), we have

By Lemma 3.3, we have

Obviously, ε−3 ≥ ε−2, then we draw

thus, it gives the wanted bound.

Lemma 3.7. For a fixed mesh and as ε → 0, it holds

where sn = nh, 1 ≤ nN − 1.

Proof: The proof is in [22].

Theorem 3.8. Let ${Z}_{n}^{m}$ be the solution of (33), then we have the following uniform error bound

Proof: Substituting Lemma 3.7 into (39), we arrive at

Hence, the result leads

Using the sup over all ε ∈ (0, 1], we get

From (46), when ε > h, the obtained method uniformly converges uniformly with order two in the space direction. When ε ≪ h, the method converges uniformly with order one in the space direction.

Theorem 3.9. Let z and Z are the solutions of (1) and (33), respectively, then we have the following uniform error bound

Proof: The proof is considered by combining of Lemma 3.1 and Theorem 3.8.

## 4. Numerical results

Considering two test examples we carry out some numerical inquiries to confirm the developed scheme is ε-uniform convergent. Since the exact solution of the examples are not known, we used a variant of double mesh principle is applied for the numerical inquiries. So, we calculate the maximum pointwise error by ${E}_{\epsilon }^{N,M}={max}_{n,m}|{Z}_{n,m}^{N,M}-{Z}_{n,m}^{2N,2M}|$, the ε-uniform error by ${E}^{N,M}={max}_{n,m}\left({E}_{\epsilon }^{N,M}\right)$, the rate of convergence by ${r}_{\epsilon }^{N,M}=log2\left({E}_{\epsilon }^{N,M}/{E}_{\epsilon }^{2N,2M}\right),$ and the ε-uniform rate of convergence by rN,M = log2(EN,M/E2N, 2M).

### 4.1. Example

Consider the problem [7]

(s, t) ∈ (0, 1) × (0, 2] with interval condition z(s, t) = 0, on (s, t) ∈ [0, 1] × [−1, 0] and the boundary conditions z(0, t) = 0 and z(1, t) = 0, t ∈ [0, 2].

### 4.2. Example

Consider the problem [13]

$∂z∂t-ε∂2z∂s2+(2-s2)∂z∂s+(s+1)(t+1)z(s,t)+z(s,t-δ)=10t2exp(-t)s(1-s),$

(s, t) ∈ (0, 1) × (0, 2] with interval condition z(s, t) = 0, on (s, t) ∈ [0, 1] × [−1, 0] and the boundary conditions z(0, t) = 0 and z(1, t) = 0, t ∈ [0, 2].

For distinguishable values of ε and N, the obtained results for the model Examples 4.1 and 4.2, respectively, ${E}_{\epsilon }^{N,M},{E}^{N,M}$, ${r}_{\epsilon }^{N,M}$, and rN,M of the developed scheme are delineated in Tables 1, 2. From these tables, one can observe that the maximum absolute error decreases as the step sizes decrease for every value of ε, and as ε approaches to zero, the maximum absolute error after getting large becomes constant, which displays ε-uniform convergence of the proposed scheme regardless of ε. On the other hand, the calculated EN,M and the corresponding rN,M using the proposed scheme are given in the last two rows, which confirms that the theoretical finding of the developed scheme is order one in both space and time direction.

TABLE 1

Table 1. ${E}_{\epsilon }^{N,M},{E}^{N,M},{r}_{\epsilon }^{N,M}$, and rN,M for Example 4.1.

TABLE 2

Table 2. ${E}_{\epsilon }^{N,M},{E}^{N,M},{r}_{\epsilon }^{N,M}$, and rN,M for Example 4.2.

In Figures 1, 2, the numerical solutions of the method for Examples 4.1 and 4.2 for different values of ε are given, respectively, for N = 80 and M = 40. Figure 3 displays the effect ε on the solutions profile of the developed scheme for Examples 4.1 and 4.2. From the figures, we see that a strong boundary layer is created on the right side of the spatial domain as ε close to zero. Furthermore, in Figure 4, the maximum point wise errors of the scheme is shown by the log-log scale plot. From these figures, one can observe that maximum absolute error decreases as the step sizes decrease for every values of ε, which confirm ε-uniform convergence of the proposed scheme.

FIGURE 1

Figure 1. Numerical solution of Example 4.1 for ε = 2−6 and ε = 2−20, respectively.

FIGURE 2

Figure 2. Numerical solution of Example 4.2 for ε = 2−6 and ε = 2−20, respectively.

FIGURE 3

Figure 3. Effect of ε on solution profiles for Examples 4.1 and 4.2, respectively.

FIGURE 4

Figure 4. Maximum point wise error in log-log scale plot for Examples 4.1 and 4.2, respectively.

In Table 3, the comparison with results of the developed method with the existing recently published studies of [23, 29] are given for Example 4.1. In Table 4, the comparison with results of the developed method with the existing number of recently published studies of [15, 24, 29, 30] are given for Example 4.2. As one follows, the developed scheme holds more accurate.

TABLE 3

Table 3. EN,M and rN,M for Example 4.1.

TABLE 4

Table 4. EN,M and rN,M for Example 4.2.

## 5. Conclusion

We have developed a numerical method for solving singularly perturbed parabolic convection-diffusion equation with a large time delay. The solution of the problem exhibits a boundary layer on the right side of the domain. The solution has a steep gradient in the layer region due to the presence of ε. In the rapidly changing behavior of the solution in the layer region, one encounters computational difficulties to find the solution using analytically or using classical numerical methods. To handle this effect, we developed method comprises of the backward Euler scheme in the time direction and an exponentially fitted higher order finite difference scheme in the spatial direction. Using comparison principle, the stability of the discrete scheme is analyzed. The stability and uniformly convergent of the method are discussed theoretically. Numerical results are delineated by applying maximum point wise error, ε-uniform error and ε-uniform rate of convergence in tables which are in acceptable agreement with the theoretical analysis. The developed method contributes more accurate, stable, and ε-uniform with a linear order of convergence in the spatial and in the time direction. The proposed scheme can be extended for singularly perturbed turning point problems.

## Data availability statement

The original contributions presented in the study are included in the article/supplementary material, further inquiries can be directed to the corresponding author.

## Author contributions

ST and MW carried out the scheme development, algorithms writing, MATLAB code writing, the numerical simulations, and write final version of the manuscript. GD and TD planned the problem, design, wrote draft of the manuscript, and revised the manuscript. All authors read, commented, and approved the submitted version of the manuscript.

## Acknowledgments

The authors would like to thank the referees for their constructive comments that improved the quality of this article.

## Conflict of interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

## Publisher's note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

## References

1. Baker CT, Bocharov GA, Rihan FA. A Report on the Use of Delay Differential Equations in Numerical Modelling in the Biosciences. Citeseer (1999).

2. Gopalsamy K. Stability and Oscillations in Delay Differential Equations of Population Dynamics. vol. 74. Springer Science and Business Media (1992).

3. Bestehorn M, Grigorieva EV. Formation and propagation of localized states in extended systems. Ann Phys. (2004) 516:423–31. doi: 10.1002/andp.200451607-806

4. Govindarao L, Mohapatra J, Das A. A fourth-order numerical scheme for singularly perturbed delay parabolic problem arising in population dynamics. J Appl Math Comput. (2020) 63:171–95. doi: 10.1007/s12190-019-01313-7

5. Tian H. Numerical Treatment of Singularly Perturbed Delay Differential Equations. Manchester, UK: The University of Manchester (2000).

6. Ansari AR, Bakr SA, Shishkin GI. A parameter-robust finite difference method for singularly perturbed delay parabolic partial differential equations. J Comput and Appl Math. (2007) 205:552–66. doi: 10.1016/j.cam.2006.05.032

7. Das A, Natesan S. Uniformly convergent hybrid numerical scheme for singularly perturbed delay parabolic convection–diffusion problems on Shishkin mesh. Appl Math Comput. (2015) 271:168–86. doi: 10.1016/j.amc.2015.08.137

8. Gowrisankar S, Natesan S. A robust numerical scheme for singularly perturbed delay parabolic initial-boundary-value problems on equidistributed grids. Electron Trans Numer Anal. (2014) 41:376–95. doi: 10.1016/j.cpc.2014.04.004

9. Gupta V, Kumar M, Kumar S. Higher order numerical approximation for time dependent singularly perturbed differential-difference convection-diffusion equations. Num Methods Partial Differ Equ. (2018) 34:357–80. doi: 10.1002/num.22203

10. Sumit, Kumar S, Kuldeep, Kumar M. A robust numerical method for a two-parameter singularly perturbed time delay parabolic problem. Comput Appl Math. (2020) 39:1–25. doi: 10.1007/s40314-020-01236-1

11. Kumar S, Kumar M. High order parameter-uniform discretization for singularly perturbed parabolic partial differential equations with time delay. Comput Math Appl. (2014) 68:1355–67. doi: 10.1016/j.camwa.2014.09.004

12. Gowrisankar S, Natesan S. ε-Uniformly convergent numerical scheme for singularly perturbed delay parabolic partial differential equations. Int J Comput Math. (2017) 94:902–21. doi: 10.1080/00207160.2016.1154948

13. Das A, Natesan S. Second-order uniformly convergent numerical method for singularly perturbed delay parabolic partial differential equations. Int J Comput Math. (2018) 95:490–510. doi: 10.1080/00207160.2017.1290439

14. Govindarao L, Mohapatra J. A second order numerical method for singularly perturbed delay parabolic partial differential equation. Eng Comput. (2018) 36:420–44. doi: 10.1108/EC-08-2018-0337

15. Kumar D, Kumari P. A parameter-uniform scheme for singularly perturbed partial differential equations with a time lag. Num Methods Partial Differ Equ. (2020) 36:868–86. doi: 10.1002/num.22455

16. Singh J, Kumar S, Kumar M. A domain decomposition method for solving singularly perturbed parabolic reaction-diffusion problems with time delay. Num Methods Partial Differ Equ. (2018) 34:1849–66. doi: 10.1002/num.22256

17. Rao SCS, Kumar S. An almost fourth order uniformly convergent domain decomposition method for a coupled system of singularly perturbed reaction–diffusion equations. J Comput Appl Math. (2011) 235:3342–54. doi: 10.1016/j.cam.2011.01.047

18. Kumar S, Kumar M. A second order uniformly convergent numerical scheme for parameterized singularly perturbed delay differential problems. Num Algor. (2017) 76:349–60. doi: 10.1007/s11075-016-0258-9

19. Kumar S, Kumar M. Analysis of some numerical methods on layer adapted meshes for singularly perturbed quasilinear systems. Num Algor. (2016) 71:139–50. doi: 10.1007/s11075-015-9989-2

20. Podila PC, Kumar K. A new stable finite difference scheme and its convergence for time-delayed singularly perturbed parabolic PDEs. Comput Appl Math. (2020) 39:1–16.

21. Tesfaye SK, Woldaregay MM, Dinka TG, Duressa GF. Fitted computational method for solving singularly perturbed small time lag problem. BMC Res Notes. (2022) 15:318. doi: 10.1186/s13104-022-06202-0

22. Woldaregay MM, Aniley WT, Duressa GF. Novel numerical scheme for singularly perturbed time delay convection-diffusion equation. Adv Math Phys. (2021) 2021:6641236. doi: 10.1155/2021/6641236

23. Negero NT, Duressa GF. An efficient numerical approach for singularly perturbed parabolic convection-diffusion problems with large time-lag. J Math Model. (2022) 10:173–90. doi: 10.22124/jmm.2021.19608.1682

24. Babu G, Bansal K. A high order robust numerical scheme for singularly perturbed delay parabolic convection diffusion problems. J Appl Math Comput. (2021) 68:363–89. doi: 10.1007/s12190-021-01512-1

25. Salama AA, Al-Amery DG. A higher order uniformly convergent method for singularly perturbed delay parabolic partial differential equations. Int J Comput Math. (2017) 94:2520–46. doi: 10.1080/00207160.2017.1284317

26. Negero NT, Duressa GF. A method of line with improved accuracy for singularly perturbed parabolic convection–diffusion problems with large temporal lag. Results Appl Math. (2021) 11:100174. doi: 10.1016/j.rinam.2021.100174

27. Sahoo SK, Gupta V. Parameter robust higher-order finite difference method for convection-diffusion problem with time delay. Num Methods Partial Differ Equ. (in press). doi: 10.1002/num.23039

28. Rao SCS, Kumar S. Second order global uniformly convergent numerical method for a coupled system of singularly perturbed initial value problems. Appl Math Comput. (2012) 219:3740–53. doi: 10.1016/j.amc.2012.09.075

29. Kumar D, Kumari P. A parameter-uniform numerical scheme for the parabolic singularly perturbed initial boundary value problems with large time delay. J Appl Math Comput. (2019) 59:179–206. doi: 10.1007/s12190-018-1174-z

30. Negero NT, Duressa GF. Uniform convergent solution of singularly perturbed parabolic differential equations with general temporal-lag. Iran J Sci Technol Trans A Sci. (2022) 46:507–24. doi: 10.1007/s40995-021-01258-2

31. Protter M, Weinberger H. Maximum Principles in Differential Equations. Englewood Cliffs, NJ: Prentice-Hall Inc. (1967).

32. Kellogg RB, Tsan A. Analysis of some difference approximations for a singular perturbation problem without turning points. Math Comp. (1978) 32:1025–39.

33. O'Malley RE. Singular Perturbation Methods for Ordinary Differential Equations. vol. 89. New York, NY: Springer (1991).

Keywords: singularly perturbed, delay differential equation, exponentially fitted finite difference, non-symmetric finite difference, uniform convergence

Citation: Tesfaye SK, Duressa GF, Woldaregay MM and Dinka TG (2023) Fitted computational method for singularly perturbed convection-diffusion equation with time delay. Front. Appl. Math. Stat. 9:1244490. doi: 10.3389/fams.2023.1244490

Received: 22 June 2023; Accepted: 07 August 2023;
Published: 31 August 2023.

Edited by:

Vikas Gupta, LNM Institute of Information Technology, India

Reviewed by:

Sunil Kumar, Indian Institute of Technology (BHU), India
Sanjay Ku Sahoo, LNM Institute of Information Technology, India

Copyright © 2023 Tesfaye, Duressa, Woldaregay and Dinka. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.

*Correspondence: Sisay Ketema Tesfaye, sisayk12@gmail.com

ORCID: Sisay Ketema Tesfaye orcid.org/0009-0005-0602-3013
Gemechis File Duressa orcid.org/0000-0003-1889-4690
Mesfin Mekuria Woldaregay orcid.org/0000-0002-6555-7534