# Comparative Study on Numerical Methods for Singularly Perturbed Advanced-Delay Differential Equations

Authors: P. Hammachukiattikul; E. Sekar; A. Tamilselvan; R. Vadivel; N. Gunasekaran; Praveen Agarwal
Journal: Journal of Mathematics
 (2021)

---

## Abstract

In this paper, we consider a class of singularly perturbed advanced-delay differential equations of convection-diffusion type. We use finite and hybrid difference schemes to solve the problem on piecewise Shishkin mesh. We have established almost first- and second-order convergence with respect to finite difference and hybrid difference methods. An error estimate is derived with the discrete norm. In the end, numerical examples are given to show the advantages of the proposed results (Mathematics Subject Classification: 65L11, 65L12, and 65L20).


---

## Body

## 1. Introduction

Differential equations depend both on past and future values (mixed delay) called functional differential equations. It attains many application problems such as optimal control problems [1], nerve conduction theory [2], economic dynamics [3], traveling waves in a spatial lattice [4] and has discussed both linear and nonlinear functional differential equations.The functional differential equation has been multiplied by small parameter0<ε<1 in the highest order derivative term called the singularly perturbed mixed delay differential equations. The main determination for such a problem is the study of biological science, epidemics, and population [5–10].The authors in [11] have considered functional differential equation in singularly perturbed problems, such as(1)σ22y″x+μ−xy′x+λEyx+aE+λIyx−aI−λE+λIyx=−1,and considered the problem of determining the expected time for the generation of action potentials in nerve cells by random synaptic inputs in the dendrites. The general linear second-order functional differential equation with the boundary-value problem arises in the modeling of neuron activation, where σ and μ are the variance and drift parameters and y is the expected first-exit time. The first-order derivative term −xy′x corresponds to exponential decay between synaptic and inputs. The undifferentiated terms correspond to excitatory and inhibitory synaptic inputs modeled as a Poisson process with mean rates λE and λI; they produce jumps in the membrane potential of amounts aE and aI, which are small quantities and could depend on the voltage. The boundary condition is(2)yx=0,∀x∉x1,x2,where the values x=x1 and x=x2 correspond to inhibitory reversal potential and the threshold value of membrane potential for action potential generation. This biological problem motivates the investigation of boundary-value problems for differential-difference equations with mixed shifts. In this biological model, using the Taylor series for the small delay term, provided the delay is of order ε, the small delay problem has oscillatory solution that has been discussed in [12]. The same authors discussed the signal transmission problem in [13].The authors in [14, 15] have considered the singularly perturbed problem with derivative depending on small delay term such as(3)εy″t+aty′t+btyt+ctyt−τ=ft,where0<τ≪1,to solve the boundary-value problem using the following numerical method such as the finite difference scheme [14, 16], fitted mesh B-spline collocation method [17], and hybrid difference scheme [15].The authors in [18, 19] investigated various concepts of singularly perturbed differential equation with derivative depending on both past and future small variables,(4)εy″t+aty′t+btyt+ctyt−τ+dtyt+τ=ft,where0<τ≪1,also proposed a finite difference scheme to solve singular perturbation problems in [18, 20, 21].The authors in [19] have been proposed to solve the singular perturbation problem with mixed small shifts using the fitted operator method. In recent years, the authors in [22–25] considered singular perturbation problem with derivative depending on large delay τ=1 variable, such as(5)εy″t+aty′t+btyt+ctyt−1=ft,has been developed various numerical schemes are finite and hybrid difference method [22], iterative scheme [26], finite element method [27, 28]. The study in [23] proposed solving singularly perturbed delay differential equation with integral boundary condition using finite difference method.Throughout the literature, the researcher concentrates on solving the singular perturbation problem with a small delay or mixed small delay or large delay using finite or hybrid or finite element methods on uniform meshes or nonuniform mesh. To the best of the author’s knowledge, up to now, no theoretical results are given for comparative study on numerical methods for singularly perturbed advanced-delay differential equations. Moreover, we proposed two numerical methods such as the finite and hybrid difference scheme on nonuniform meshes, to solve the singular perturbation problem with mixed large delay using the finite difference scheme and hybrid difference scheme on Shishkin mesh.This paper is structured as follows: Section2 describes the problem statement. Section 3 proves the maximum principle and stability result. Moreover, it introduces the terminology for Shishkin decomposition and proves many inequalities. In Section 4, we introduce the numerical methods to discretize the continuous problem. Error analysis for finite and hybrid difference scheme approximate solution is given in Sections 5 and 6. Finally, Section 7 presents numerical results.Throughout our analysis, we use the following notations:Γ¯=0,3, Γ=0,3, Γ1=0,1, Γ2=1,2, Γ3=2,3, Γ∗=Γ1∪Γ2∪Γ3. Γ¯3N=0,1,2,…,3N, Γ13N=1,2,…,N−1, Γ23N=N+1,…,2N−1, Γ33N=2N+1,…,3N−1. The parameter ε and mesh points 3N are independent of C and C1 are positive constants. The norm is yΓ=supr∈Γyr.


## 2. Statement of the Problem

Consider the following singularly perturbed mixed delay differential equation:(6)Kyr=−εy″r+ary′r+bryr+cryr−1+dryr+1=fr,r∈Γ,yr=ϕr,r∈−1,0,yr=φr,r∈3,4,where ϕr and φr are history function on −1,0 and 3,4. Assume that ar≥α1>α>0, br≥β≥0, γ≤cr≤0, dr≥η≥0, β+γ≥β1>0, and the coefficients are smooth function on r∈Γ¯. The above problem solution satisfies yr∈G=C0Γ¯∩C1Γ∩C2Γ∗. Problem (1) is rewritten as Kyr=gr, where(7)Kyr=K1yr=−εy″r+ary′r+bryr+dryr+1,r∈Γ1,K2yr=−εy″r+ary′r+bryr+cryr−1+dryr+1,r∈Γ2,K3yr=−εy″r+ary′r+bryr+cryr−1,r∈Γ3,(8)gr=fr−crϕr−1,r∈Γ1,fr,r∈Γ2,fr−drφr+1,r∈Γ3,with boundary conditions(9)yr=ϕr,r∈−1,0,y1−=y1+,y′1−=y′1+,y2−=y2+,y′2−=y′2+,yr=φr,r∈3,4.


## 3. Analytical Results

Lemma 1 (maximum principle).
Ifyr∈G such that y0≥0, y3≥0, K1yr≥0, ∀r∈Γ1, K2yr≥0, ∀r∈Γ2, K3yr≥0, ∀r∈Γ3, y′1≤0, and y′2≤0, then yr≥0, ∀r∈Γ¯.
Proof. Let(10)sr=112+r4,r∈0,1,212+r6,r∈1,2,412+r12,r∈2,3.Clearly,sr>0,∀x∈Γ¯, Ksr>0,∀r∈Γ∗, s0>0, s3>0, s′1<0, and s′2<0. Consider that μ=max−ψr/sr:r∈Γ¯; then, there exists r0∈Γ¯ such that ψr0+μsr0=0 and ψr+μsr≥0,∀r∈Γ¯ implies that ψ+μs obtain minimum at t=r0. If μ<0, then ψr≥0.Ifμ>0, then the function ψr nonnegative is not possible. The following cases are easy to prove the contradiction if μ>0.Case (i):r0=0:(11)0<ψ+μs,0=ψ0+μs0=0.Case (ii):r0∈Γ1:(12)0<K1ψ+μsr0≤0.Case (iii):r0=1:(13)0≤ψ+μs′1=ψ′1+μs′1<0.Case (iv):r0∈Γ2:(14)0<K2ψ+μsr0≤0.Case (v):r0=2:(15)0≤ψ+μs′2=ψ′2+μs′2<0.Case (vi):r0∈Γ3:(16)0<K3ψ+μsr0≤0.Case (vii):r0=3:(17)0≤ψ+μs3=ψ3+μs3=0.All the cases are contradiction.Lemma 2 (stability result).
Ifyris a solution of problems (7)–(9), then(18)yr≤Cmaxy0,y3,supr∈Γ∗Kyr,r∈Γ¯.Lemma 3.
Ifyris a solution of problems (7)–(9), then(19)ykrΓ∗≤Cε−k,wherek=1,2,3,4.Proof.
First, to provey′r is bound on Γ1, (20)K1yr=−εy″r+ary′r+bryr+dryr+1=fr−crϕr−1.
Integrating the above equation on both sides, we have(21)−εy′r−y′0=−aryr−a0u0+∫0ra′tytdt−∫0rbtytdt−∫0rdtyt+1dt+∫0rft−ctϕt−1dt.
Therefore,(22)εy′0=εy′r−aryr−a0y0+∫0ra′tytdt−∫0rbtytdt−∫0rdtyt+1dt+∫0rft−ctϕt−1dt.
Using Mean Value Theorem, thenεy′x≤Cyr,fr,ϕ−1,0, for some x∈0,ε and εy′0≤Cyr+fr+ϕr. Then, we have εy′r≤Cmaxyr,fr,ϕ,r∈Γ1.
To provey′r is bound on Γ2, (23)K2yr=−εy″r+ary′r+bryr+cryr−1+dryr+1=fr.
Integrating the above equation on both sides, we have(24)−εy′r−y′0=−aryr−a0u0+∫0ra′tytdt−∫0rbtytdt−∫0rctyt−1dt−∫0rdtyt+1dt+∫0rftdt.
Therefore,(25)εy′0=εy′r−aryr−a0u0+∫0ra′tytdt−∫0rbtytdt−∫0rctyt−1dt−∫0rdtyt+1dt+∫0rftdt.
Using Mean Value Theorem, thenεy′x≤Cyr,fr,ϕ−1,0, for some x∈0,ε and εy′0≤Cyr+fr+ϕr. Then, we have εy′r≤Cmaxyr,fr,ϕ,r∈Γ2.
Next, to provey′r is bound on Γ3,(26)K3yr=−εy″r+ary′r+bryr+cryr−1=fr−drφr+1.
Integrating the above equation on both sides, we have(27)−εy′r−y′0=−aryr−a0u0+∫0ra′tytdt−∫0rbtytdt−∫0rctyt−1dt−∫0rdtyt+1dt+∫0rftdt.
Therefore,(28)εy′0=εy′r−aryr−a0u0+∫0ra′tytdt−∫0rbtytdt−∫0rctyt−1dt−∫0rdtyt+1dt+∫0rftdt.
Using Mean Value Theorem, thenεy′x≤Cyr,fr,ϕ−1,0, for some x∈0,ε and εy′0≤Cyr+fr+ϕr. Then, we have, εy′r≤Cmaxyr,fr,ϕ,r∈Γ3.
Hence,ykrΓ∗≤Cε−k,wherek=2,3,4.


### 3.1. Shishkin Decomposition

The solutionyr is decomposed into vr smooth component and wr-layer component. Furthermore, vr=v0r+εv1r+ε2v2r, where v0r, v1r, and v2r are solutions of the following differential equations.Obtain reduced problem solutionv0r∈X such that(29a)arv0′r+brv0r+crv0r−1=fr,r∈Γ∩Γ1∪Γ2,v0r=ϕr,r∈−1,0.(29b)arv0′r+brv0r+drv0r+1=fr,r∈Γ∩Γ2∪Γ3∪3,v0r=φr,r∈3,4.Ifv1r∈C0Γ¯∩C1Γ∗∪3,(30a)arv1′r+brv1r+crv1r−1=v0″r,r∈Γ∩Γ1∪Γ2,v1r=0,r∈−1,0.(30b)arv1′r+brv1r+drv1r+1=v0″r,r∈Γ∩Γ2∪Γ3∪3,v1r=0,r∈3,4.Ifv2r∈X, (31)−εv2″r+arv2″r+brv2r+crv2r−1+drv2r+1=v1″r,r∈Γ∗v2r=0,r∈−1,0,v2r=0,r∈3,4.Ifvr∈C0Γ¯∩C2Γ∗,(32)Kvr=−εv″r+arv′r+brvr+crvr−1+drvr+1=fr,r∈Γ∗,vr=ϕr,r∈−1,0,v1=v01+εv11+ε2v21,v2=v02+εv12+ε2v22,vr=0,r∈3,4.Also,wr satisfies the following problem: if the singular component wr∈C0Γ¯∩C2Γ∗,(33)Kwr=−εw″r+arw′r+brwr+crwr−1+drwr+1=0,r∈Γ∗,wr=0,r∈−1,0,w′1=−v′1,w′2=−v′2,wr=φr,r∈3,4.Furthermore, we decomposewr as wr=wBr+wI1r+wI2r, where the function wBr is boundary layer component and wI1r, wI2r are interior layer components.If the boundary layerwBr∈X,(34)KwBr=−εwB″r+arwB′r+brwBr+crwBr−1+drwBr+1=0,wBr=0,r∈−1,0,wBr=φr,r∈3,4.If the first interior layerwI1r∈C0Γ¯∩C2Γ∗,(35)KwI1r=−εwI1′′r+arwI1′r+brwI1r+crwI1r−1+drwI1r+1=0,wI1r=0,r∈−1,0,wI1′1=−v′1,wI1r=0,r∈3,4.If the second interior layerwI2r∈C0Γ¯∩C2Γ∗,(36)KwI2r=−εwI2′′r+arwI2′r+brwI2r+crwI2r−1+drwI2r+1=0,wI2r=0,r∈−1,0,wI2′2=−v′2,wI2r=0,r∈3,4.Theorem 1.
Ifyrandv0rare solutions of problems (7)–(9) and (29a)-(29b), then(37)yr−v0r≤C1ε+exp−α3−rε,r∈Γ¯.Proof.
Consider(38)Θ±r=C1εsr+exp−α3−rε±yr−v0r,r∈Γ¯.
Clearly,Θ±r∈C0Γ¯∩C2Γ∗. Note that Θ±0≥0, Θ±3≥0, for a suitable choice of C1>0.
Ifr∈Γ1,(39)K1Θ±r=C1αεar−α+br+drexpαεexp−α3−rε+εars′r+brsr+drsr+1±εv0″r,≥C1αεα1−α+β+ηexpαεexp−α3−rε+εα4+βsr+ηsr+1±Cε≥0.
Ifr∈Γ2,(40)K2Θ±r=C1αεar−α+br+crexp−αε+drexpαεexp−α3−rε+εars′r+brsr+crsr−1+drsr+1±εv0″r,≥C1αεα1−α+β+γ+ηexpαεexp−α3−rε+εα6+βsr+γsr−1+ηsr+1∓Cε≥0.
Following the same process, we haveK3Θ±r≥0.
Using Lemma1, then Θ±r≥0,r∈Γ¯. Therefore, yr−v0r≤C1ε+exp−α3−r/ε.Lemma 4.
Ifvrandwrare the solution of regular and singular component problems (32) and (33), then(41)vkrΓ∗≤C1+ε2−k,fork=0,1,2,3,4,(42)wBkr≤Cε−kexp−α3−rε,r∈Γ∗,(43)wI1kr≤Cε1−kexp−α1−rε,r∈Γ1,ε1−k,r∈Γ2,r∈Γ3,(44)wI2kr≤Cε1−kexp−α2−rε,r∈Γ2,ε1−k,r∈Γ1,r∈Γ3,where k = 0, 1, 2, 3, 4.Proof.
The smooth component derivative bound is easy to prove by using stability result and integrating (30a), (30b), and (31). Next, to prove (42), consider that(45)Φ±r=Cexp−α3−rε±wBr,r∈Γ¯.
Note thatΦ±0≥0, Φ±3≥0, and(46)KΦ±r≥Cαεα1−α+β+γ+ηexpαεexp−α3−rε±KwBr≥0.
By Lemma1,(47)wBr≤Cexp−α3−rε.
Integration of (34) yields the estimates of wB′r. From the differential equations (33), one can derive the rest of the derivative estimates (42).
Inequalities (43) and (44) can be proved, using Theorem 1 and maximum principle for the barrier functions:(48)Φ±r=Cεexp−α1−rε±wI1r,r∈Γ1,rε±wI1r,r∈Γ2,rε±wI1r,r∈Γ3,Φ±r=Cεexp−α2−rε±wI2r,r∈Γ1,rε±wI2r,r∈Γ2,rε±wI2r,r∈Γ3.
Hence, it is proved.Remark.
The following inequalities are easy to prove, using Theorem1 and Lemma 4:(49)yr−vr≤Cε+εexp−α1−rε+exp−α3−rε,r∈Γ1,ε+εexp−α2−rε+exp−α3−rε,r∈Γ2,ε+exp−α3−rε,r∈Γ3.


## 3.1. Shishkin Decomposition

The solutionyr is decomposed into vr smooth component and wr-layer component. Furthermore, vr=v0r+εv1r+ε2v2r, where v0r, v1r, and v2r are solutions of the following differential equations.Obtain reduced problem solutionv0r∈X such that(29a)arv0′r+brv0r+crv0r−1=fr,r∈Γ∩Γ1∪Γ2,v0r=ϕr,r∈−1,0.(29b)arv0′r+brv0r+drv0r+1=fr,r∈Γ∩Γ2∪Γ3∪3,v0r=φr,r∈3,4.Ifv1r∈C0Γ¯∩C1Γ∗∪3,(30a)arv1′r+brv1r+crv1r−1=v0″r,r∈Γ∩Γ1∪Γ2,v1r=0,r∈−1,0.(30b)arv1′r+brv1r+drv1r+1=v0″r,r∈Γ∩Γ2∪Γ3∪3,v1r=0,r∈3,4.Ifv2r∈X, (31)−εv2″r+arv2″r+brv2r+crv2r−1+drv2r+1=v1″r,r∈Γ∗v2r=0,r∈−1,0,v2r=0,r∈3,4.Ifvr∈C0Γ¯∩C2Γ∗,(32)Kvr=−εv″r+arv′r+brvr+crvr−1+drvr+1=fr,r∈Γ∗,vr=ϕr,r∈−1,0,v1=v01+εv11+ε2v21,v2=v02+εv12+ε2v22,vr=0,r∈3,4.Also,wr satisfies the following problem: if the singular component wr∈C0Γ¯∩C2Γ∗,(33)Kwr=−εw″r+arw′r+brwr+crwr−1+drwr+1=0,r∈Γ∗,wr=0,r∈−1,0,w′1=−v′1,w′2=−v′2,wr=φr,r∈3,4.Furthermore, we decomposewr as wr=wBr+wI1r+wI2r, where the function wBr is boundary layer component and wI1r, wI2r are interior layer components.If the boundary layerwBr∈X,(34)KwBr=−εwB″r+arwB′r+brwBr+crwBr−1+drwBr+1=0,wBr=0,r∈−1,0,wBr=φr,r∈3,4.If the first interior layerwI1r∈C0Γ¯∩C2Γ∗,(35)KwI1r=−εwI1′′r+arwI1′r+brwI1r+crwI1r−1+drwI1r+1=0,wI1r=0,r∈−1,0,wI1′1=−v′1,wI1r=0,r∈3,4.If the second interior layerwI2r∈C0Γ¯∩C2Γ∗,(36)KwI2r=−εwI2′′r+arwI2′r+brwI2r+crwI2r−1+drwI2r+1=0,wI2r=0,r∈−1,0,wI2′2=−v′2,wI2r=0,r∈3,4.Theorem 1.
Ifyrandv0rare solutions of problems (7)–(9) and (29a)-(29b), then(37)yr−v0r≤C1ε+exp−α3−rε,r∈Γ¯.Proof.
Consider(38)Θ±r=C1εsr+exp−α3−rε±yr−v0r,r∈Γ¯.
Clearly,Θ±r∈C0Γ¯∩C2Γ∗. Note that Θ±0≥0, Θ±3≥0, for a suitable choice of C1>0.
Ifr∈Γ1,(39)K1Θ±r=C1αεar−α+br+drexpαεexp−α3−rε+εars′r+brsr+drsr+1±εv0″r,≥C1αεα1−α+β+ηexpαεexp−α3−rε+εα4+βsr+ηsr+1±Cε≥0.
Ifr∈Γ2,(40)K2Θ±r=C1αεar−α+br+crexp−αε+drexpαεexp−α3−rε+εars′r+brsr+crsr−1+drsr+1±εv0″r,≥C1αεα1−α+β+γ+ηexpαεexp−α3−rε+εα6+βsr+γsr−1+ηsr+1∓Cε≥0.
Following the same process, we haveK3Θ±r≥0.
Using Lemma1, then Θ±r≥0,r∈Γ¯. Therefore, yr−v0r≤C1ε+exp−α3−r/ε.Lemma 4.
Ifvrandwrare the solution of regular and singular component problems (32) and (33), then(41)vkrΓ∗≤C1+ε2−k,fork=0,1,2,3,4,(42)wBkr≤Cε−kexp−α3−rε,r∈Γ∗,(43)wI1kr≤Cε1−kexp−α1−rε,r∈Γ1,ε1−k,r∈Γ2,r∈Γ3,(44)wI2kr≤Cε1−kexp−α2−rε,r∈Γ2,ε1−k,r∈Γ1,r∈Γ3,where k = 0, 1, 2, 3, 4.Proof.
The smooth component derivative bound is easy to prove by using stability result and integrating (30a), (30b), and (31). Next, to prove (42), consider that(45)Φ±r=Cexp−α3−rε±wBr,r∈Γ¯.
Note thatΦ±0≥0, Φ±3≥0, and(46)KΦ±r≥Cαεα1−α+β+γ+ηexpαεexp−α3−rε±KwBr≥0.
By Lemma1,(47)wBr≤Cexp−α3−rε.
Integration of (34) yields the estimates of wB′r. From the differential equations (33), one can derive the rest of the derivative estimates (42).
Inequalities (43) and (44) can be proved, using Theorem 1 and maximum principle for the barrier functions:(48)Φ±r=Cεexp−α1−rε±wI1r,r∈Γ1,rε±wI1r,r∈Γ2,rε±wI1r,r∈Γ3,Φ±r=Cεexp−α2−rε±wI2r,r∈Γ1,rε±wI2r,r∈Γ2,rε±wI2r,r∈Γ3.
Hence, it is proved.Remark.
The following inequalities are easy to prove, using Theorem1 and Lemma 4:(49)yr−vr≤Cε+εexp−α1−rε+exp−α3−rε,r∈Γ1,ε+εexp−α2−rε+exp−α3−rε,r∈Γ2,ε+exp−α3−rε,r∈Γ3.


## 4. The Discrete Problem

### 4.1. Shishkin Mesh

Problems (7)–(9) are convection-diffusion type containing delay term. Then, the layers occur in boundary at t=3 and interior at t=1 and t=2.The intervals0,1, 1,2, and 2,3 are partitioned into 0,1−σ, 1−σ,1, 1,2−σ, 2−σ,2, 2,3−σ, and 3−σ,3 for each interval N/2 mesh points and σ=min1/2,2ε/αlnN is transition parameter.The interior of points is denoted byΓ¯3N=r0,r1,…,r3N. Then, the mesh widths are(50)hri=H=21−σN,fori=1toN2,i=N+1to3N2,i=2N+1to5N2,h=2σN,fori=N2+1toN,i=3N2+1to2N,i=5N2+1to3N.


### 4.2. Finite Difference Method

The discrete scheme corresponding to the original problems (7)–(9) is as follows:(51)KNYri=K1NYri=−εδ2Yri+ariD−Yri+briYri+driYri+N=fi−ciϕi−N,ri∈Γ13N,K2NYri=−εδ2Yri+ariD−Yri+briYri+criYri−N+driYri+N=fi,ri∈Γ23N,K3NYri=−εδ2Yri+ariD−Yri+briYri+criYri−N=fi−diφi+N,ri∈Γ33N,with(52)Yr0=ϕ0,D−UN=D+UN,D−U2N=D+U2N,Yr3=φ3N.


### 4.3. Hybrid Difference Scheme

The hybrid scheme corresponding to the original problems (7)–(9) is as follows:(53)K1NYri=−εδ2Yri+ai−1/2D−Yri+briU^ri+driU^ri+N=fi−1/2−ciϕ^i−N,i=1toN2,−εδ2Yri+ariD0Yri+briYri+driYri+N=fi−ciϕi−N,i=N2+1toN−1,(54)K2NYri=−εδ2Yri+ai−1/2D−Yri+briU^ri+criU^ri−N+driU^ri+N=fi−1/2,i=N+1to3N2,−εδ2Yri+ariD0Yri+briYri+criYri−N+driYri+N=fi,i=3N2+1to2N−1,(55)K3NYri=−εδ2Yri+ai−1/2D−Yri+briU^ri+criU^ri−N=fi−1/2−diφ^i+N,i=2N+1to5N2,−εδ2Yri+ariD0Yri+briYri+criYri−N=fi−diφi+N,i=5N2+1to3N−1,(56)KNYri=Yi−2−4Yi−1+3Yi2h−−Yi+2+4Yi+1−3Yi2H,i=N,2N,where(57)δ2Yri=2hi+hi+1Yri+1−Yrihi+1−Yri−Yri−1hi,U^ri=Yri+Yri−12,D0Yri=Yri+1−Yri−1hi+hi+1,D−Yri=Yri−Yri−1hi,ai−1/2=ari−1+ri2.


## 4.1. Shishkin Mesh

Problems (7)–(9) are convection-diffusion type containing delay term. Then, the layers occur in boundary at t=3 and interior at t=1 and t=2.The intervals0,1, 1,2, and 2,3 are partitioned into 0,1−σ, 1−σ,1, 1,2−σ, 2−σ,2, 2,3−σ, and 3−σ,3 for each interval N/2 mesh points and σ=min1/2,2ε/αlnN is transition parameter.The interior of points is denoted byΓ¯3N=r0,r1,…,r3N. Then, the mesh widths are(50)hri=H=21−σN,fori=1toN2,i=N+1to3N2,i=2N+1to5N2,h=2σN,fori=N2+1toN,i=3N2+1to2N,i=5N2+1to3N.


## 4.2. Finite Difference Method

The discrete scheme corresponding to the original problems (7)–(9) is as follows:(51)KNYri=K1NYri=−εδ2Yri+ariD−Yri+briYri+driYri+N=fi−ciϕi−N,ri∈Γ13N,K2NYri=−εδ2Yri+ariD−Yri+briYri+criYri−N+driYri+N=fi,ri∈Γ23N,K3NYri=−εδ2Yri+ariD−Yri+briYri+criYri−N=fi−diφi+N,ri∈Γ33N,with(52)Yr0=ϕ0,D−UN=D+UN,D−U2N=D+U2N,Yr3=φ3N.


## 4.3. Hybrid Difference Scheme

The hybrid scheme corresponding to the original problems (7)–(9) is as follows:(53)K1NYri=−εδ2Yri+ai−1/2D−Yri+briU^ri+driU^ri+N=fi−1/2−ciϕ^i−N,i=1toN2,−εδ2Yri+ariD0Yri+briYri+driYri+N=fi−ciϕi−N,i=N2+1toN−1,(54)K2NYri=−εδ2Yri+ai−1/2D−Yri+briU^ri+criU^ri−N+driU^ri+N=fi−1/2,i=N+1to3N2,−εδ2Yri+ariD0Yri+briYri+criYri−N+driYri+N=fi,i=3N2+1to2N−1,(55)K3NYri=−εδ2Yri+ai−1/2D−Yri+briU^ri+criU^ri−N=fi−1/2−diφ^i+N,i=2N+1to5N2,−εδ2Yri+ariD0Yri+briYri+criYri−N=fi−diφi+N,i=5N2+1to3N−1,(56)KNYri=Yi−2−4Yi−1+3Yi2h−−Yi+2+4Yi+1−3Yi2H,i=N,2N,where(57)δ2Yri=2hi+hi+1Yri+1−Yrihi+1−Yri−Yri−1hi,U^ri=Yri+Yri−12,D0Yri=Yri+1−Yri−1hi+hi+1,D−Yri=Yri−Yri−1hi,ai−1/2=ari−1+ri2.


## 5. Numerical Estimates for the Finite Difference Method

Lemma 5 (discrete maximum principle).
IfUrisatisfiesUr0≥0,Ur3N≥0,K1NUri≥0,K2NUri≥0,K3NUri≥0,D+UrN−D−UrN≤0, andD+Ur2N−D−Ur2N≤0, thenUri≥0,∀ri∈Γ¯3N.Proof.
(58)DefineSri=112+ri4,ri∈0,1∩Γ¯3N,212+ri6,ri∈1,2∩Γ¯3N,412+ri12,ri∈2,3∩Γ¯3N.
It is easy to see thatsri>0,∀ri∈Γ¯3N, Ksri>0,∀ri∈Γ13N∪Γ23N∪Γ33N, sr0>0, sr3N>0, D+srN−D−srN<0, and D+sr2N−D−sr2N<0. Let μ=min−ψri/sri:ri∈Γ¯3N.
Then, there existsrk∈Γ¯3N such that ψrk+μsrk=0 and ψri+μsri≥0,∀ri∈Γ¯3N. Then, ψ+μs attains its maximum at ri=rk. If μ<0, then ψ≥0. Suppose μ>0.Case (i):rk=r0:(59)0<ψ+μsr0=0.Case (ii):rk∈Γ13N:(60)0<K1ψ+μsrk=−εδ2ψ+μsrk+arkD−ψ+μsrk+brkψ+μsrk+drkψ+μsrk+N≤0.Case (iii):rk=rN:(61)0≤Dψ+μsrN<0.Case (iv):rk∈Γ23N:(62)0<K2ψ+μsrk=−εδ2ψ+μsrk+arkD−ψ+μsrk+brkψ+μsrk+crkψ+μsrk−N+drkψ+μsrk+N≤0.Case (v):rk=r2N:(63)0≤Dψ+μsr2N<0.Case (vi):rk∈Γ33N:(64)0<K3ψ+μsrk=−εδ2ψ+μsrk+arkD−ψ+μsrk+brkψ+μsrk+crkψ+μsrk−N≤0.Case (vii):rk=r3N:(65)0<ψ+μsr3N=0.
All the cases are a contradiction.Lemma 6.
The discrete solution of (51) and (52) is bounded:(66)Uri≤CmaxUr0,Ur3N,maxi∈Γ13N∪Γ23N∪Γ33NKNUri.Proof.
Considerψ±ri=CMsri±Uri,0≤i≤3N, where M=maxUr0,Ur3N,maxi∈Γ13N∪Γ23N∪Γ33NKNUri.
Observeψ±r0≥0 and ψ±r3N≥0:(67)K1Nψ±ri≥0,∀ri∈Γ1N,K2Nψ±ri≥0,∀ri∈Γ2N,K3Nψ±ri≥0,∀ri∈Γ3N,D+ψ±rN−D−ψ±rN≤0.D+ψ±r2N−D−ψ±r2N≤0.
Using Lemma5, ψ±ri≥0, ∀ri∈Γ¯3N.
To decompose numerical solutionYri into Vri and Wri satisfy the following equations, respectively:(68)KNVri=−εδ2Vri+ariD−Vri+briVri+criVri−N+driVri+N=fi,i∈Γ¯3N\0,N,2N,3N,Vr0=v0,DVrN=v′1,DVr2N=v′2,Vr3N=v3,(69)KNWri=−εδ2Wri+ariD−Wri+briWri+criWri−N+driWri+N=fi,i∈Γ¯3N\0,N,2N,3N,Wr0=w0,DWrN=−DVrN,DWr2N=−DVr2N,Wr3N=w3.Theorem 2.
IfYriandVriare a solution of discretization problem (51), (52), and (68), thenYri−Vri≤CN−1.
Proof. Consider(70)θ±ri=CN−1sri±Yri−Vri,∀ri∈Γ¯3N.Note thatθ±r0≥0 and θ±r3N≥0:(71)K1Nθ±ri≥0,for alli∈1,2,…,N−1,K2Nθ±ri≥0,for alli∈N+1,…,2N−1,K3Nθ±ri≥0,for alli∈2N+1,…,3N−1,Dθ±rN=−CN−1/6±v′1<0, Dθ±r2N<0, using Lemma 5; then, the theorem has been proved.Theorem 3.
The error estimates for smooth components bounded byCN−1:(72)vri−Vri≤CN−1,ri∈Γ¯3N.Proof.
The proof of Theorem3 has the same idea in [29]:(73)KNvri−Vri≤CN−1,i∈Γ13N∪Γ23N∪Γ33N.
Using Lemma6, then(74)vri−Vri≤CN−1,i∈Γ13N∪Γ23N∪Γ33N.
Therefore, we getvri−Vri≤CN−1,i∈Γ¯3N.Theorem 4.
Derive the error estimates for singular components bounded byCN−1log2N:(75)wri−Wri≤CN−1log2N,ri∈Γ¯3N.Proof.
Note that(76)wri−Wri≤yri−Yri+vri−Vri.
Then, by (49) and from Theorems 1, 3, we have(77)yri−Yri=Yri−Vri+vri−Vri+yri−vri,yri−Yri≤Yri−Vri+vri−Vri+yri−vri.
Now,(78)wri−Wri≤Yri−Vri+2vri−Vri+yri−vri≤C1N−1+C1exp−α3−rε+ε≤C1N−1+C1exp−α3−rε+C1N−1≤C1exp−ασε+C1N−1≤CN−1,i=0to5N2.
Consider the mesh functions:(79)Φ±ri=C1N−1sri+C1N−1σε2ri−3−σ±wri−Wri,ri∈3−σ,3∩Γ¯3N.
Observe thatΦ±r5N/2≥0 and Φ±r3N≥0, and KNΦ±ri≥0.
Then, by the Lemma5, we have Φ±ri≥0,ri∈Γ¯3N. Therefore,(80)wri−Wri≤CN−1log2N,ri∈Γ¯3N.Theorem 5.
IfyriandYriare a solution of (7)–(9) and (51), (52),(81)yri−Yri≤CN−1log2N,ri∈Γ¯3N.That is, the order of convergence is almost one.Proof.
The proof of Theorem5 follows from yk=vk+wk, Yk=Vk+Wk, and Theorems 3 and 4.


## 6. Numerical Estimates for the Hybrid Difference Method

Assume the following inequality:(82)NlnN≥2α1α.Lemma 7.
Assume (78) holds true. LetΨrisatisfyΨr0≥0,Ψr3N≥0; the operatorKNdefined by (53)–(55) satisfiesK1NΨri≥0,K2NΨri≥0,K3NΨri≥0; and thenΨri≥0, ∀ri∈Γ¯3N.Lemma 8.
IfΨriis discrete solution of problems (53)–(55), then(83)Ψri≤CmaxΨr0,Ψr3N,maxi∈Γ13N∪Γ23N∪Γ33NKNΨri.


### 6.1. Error Estimate

To decompose the numerical solutionYri into Vri and Wri, satisfy the following equations, respectively:(84)KNVri=−εδ2Vri+ai−1/2D−Vri+briV^ri+criV^ri−N+driU^ri+N=fi−1/2,−εδ2Vri+ariD0Vri+briVri+criVri−N+driVri+N=fi.(85)KNWri=−εδ2Wri+ai−1/2D−Wri+briW^ri+criW^ri−N+driW^ri+N=0,−εδ2Wri+ariD0Wri+briWri+criWri−N+driWri+N=0.Lemma 9.
Derive the error estimation of discretization original problems (53)–(56) and regular problem (84) solutions:(86)Yri−Vri≤CN−2.Proof.
The proof of Lemma9 has the same idea in Lemma 7:(87)θ±ri=CN−2sri±Yri−Vri,∀ri∈Γ¯3N.Lemma 10.
The error estimates for smooth components are bounded byCN−2:(88)Vri−vri≤CN−2,ri∈Γ¯3N.Proof.
Utilizing the method adopted in [30],(89)KY−yri≤CεHy3+H2y3+y2,i=1toN2,εh2y4+h2aiy3,i=N2+1toN−1.
Usingε≤CN−1 and the above equation, the bounds on the derivatives of v can be written as(90)K1V−vri≤CεHv3+H2v3+v2,i=1toN2,εh2v4+h2aiv3+h2civ3,i=N2+1toN−1,K1V−vri≤CN−1ε+N−1,i=1toN2N−2,i=N2+1toN−1.
Then, we haveK1V−vri≤CN−2. Similarly, KjV−vri≤CN−2,j=2,3, KjV−vri≤CN−2,j=1,2,3,i∈Γ¯3N\0,N,2N,3N, and, by Lemma 8, we have(91)Vri−vri≤CN−2,ri∈Γ¯3N.Lemma 11.
Derive the error estimates for singular components bounded byCN−2log2N:(92)wri−Wri≤CN−2log2N,ri∈Γ¯3N.Proof.
Note that(93)wri−Wri≤yri−Yri+vri−Vri,yri−Yri≤Yri−Vri+vri−Vri+yri−vri.
Now,(94)wri−Wri≤Yri−Vri+2vri−Vri+yri−vri,≤C1exp−ασε+C1N−2≤CN−2,i=0to5N2.
Consider the mesh functions(95)Φ±ri=C1N−2sri+C1N−2σε2ri−3−σ±wri−Wri,ri∈3−σ,3∩Γ¯3N.
Clearly,Φ±r5N/2≥0 and Φ±r3N≥0, for a suitable choice of C1>0.(96)KNΦ±ri≥0.
Then, by Lemma7, we have Φ±ri≥0,ri∈Γ¯3N. Therefore,(97)wri−Wri≤CN−2log2N,ri∈Γ¯3N.Theorem 6.
IfyriandYriare the solution of (7)–(9) and (53)–(56), then(98)yri−Yri≤CN−2log2N,ri∈Γ¯3N.Proof.
The proof of Theorem6 follows from yk=vk+wk and Yk=Vk+Wk and using Theorems 3 and 4


## 6.1. Error Estimate

To decompose the numerical solutionYri into Vri and Wri, satisfy the following equations, respectively:(84)KNVri=−εδ2Vri+ai−1/2D−Vri+briV^ri+criV^ri−N+driU^ri+N=fi−1/2,−εδ2Vri+ariD0Vri+briVri+criVri−N+driVri+N=fi.(85)KNWri=−εδ2Wri+ai−1/2D−Wri+briW^ri+criW^ri−N+driW^ri+N=0,−εδ2Wri+ariD0Wri+briWri+criWri−N+driWri+N=0.Lemma 9.
Derive the error estimation of discretization original problems (53)–(56) and regular problem (84) solutions:(86)Yri−Vri≤CN−2.Proof.
The proof of Lemma9 has the same idea in Lemma 7:(87)θ±ri=CN−2sri±Yri−Vri,∀ri∈Γ¯3N.Lemma 10.
The error estimates for smooth components are bounded byCN−2:(88)Vri−vri≤CN−2,ri∈Γ¯3N.Proof.
Utilizing the method adopted in [30],(89)KY−yri≤CεHy3+H2y3+y2,i=1toN2,εh2y4+h2aiy3,i=N2+1toN−1.
Usingε≤CN−1 and the above equation, the bounds on the derivatives of v can be written as(90)K1V−vri≤CεHv3+H2v3+v2,i=1toN2,εh2v4+h2aiv3+h2civ3,i=N2+1toN−1,K1V−vri≤CN−1ε+N−1,i=1toN2N−2,i=N2+1toN−1.
Then, we haveK1V−vri≤CN−2. Similarly, KjV−vri≤CN−2,j=2,3, KjV−vri≤CN−2,j=1,2,3,i∈Γ¯3N\0,N,2N,3N, and, by Lemma 8, we have(91)Vri−vri≤CN−2,ri∈Γ¯3N.Lemma 11.
Derive the error estimates for singular components bounded byCN−2log2N:(92)wri−Wri≤CN−2log2N,ri∈Γ¯3N.Proof.
Note that(93)wri−Wri≤yri−Yri+vri−Vri,yri−Yri≤Yri−Vri+vri−Vri+yri−vri.
Now,(94)wri−Wri≤Yri−Vri+2vri−Vri+yri−vri,≤C1exp−ασε+C1N−2≤CN−2,i=0to5N2.
Consider the mesh functions(95)Φ±ri=C1N−2sri+C1N−2σε2ri−3−σ±wri−Wri,ri∈3−σ,3∩Γ¯3N.
Clearly,Φ±r5N/2≥0 and Φ±r3N≥0, for a suitable choice of C1>0.(96)KNΦ±ri≥0.
Then, by Lemma7, we have Φ±ri≥0,ri∈Γ¯3N. Therefore,(97)wri−Wri≤CN−2log2N,ri∈Γ¯3N.Theorem 6.
IfyriandYriare the solution of (7)–(9) and (53)–(56), then(98)yri−Yri≤CN−2log2N,ri∈Γ¯3N.Proof.
The proof of Theorem6 follows from yk=vk+wk and Yk=Vk+Wk and using Theorems 3 and 4


## 7. Numerical Experiments

In this section, consider two examples for constant and variable coefficient problems and apply both of the numerical methods to find error and rate of convergence. The exact solution is not easy to find in these problems. Therefore, we use the double mesh principle:(99)DεN=max0≤i≤3NUiN−U2i2N.We compute the uniform error and the rate of convergence as(100)DN=maxεDεN,PN=log2DND2N.To solve the following numerical examples, we use two computational methods such as finite and hybrid difference scheme on the nonuniform mesh.Example 1.
(101)−εy″r+5y′r+2yr−yr−1+yr+1=1,forr∈Γ∗,yr=1,forr∈−1,0,yr=2,forr∈3,4.Example 2.
(102)−εy″r+r+5y′r+2yr−yr−1+x2yr+1=et,forr∈Γ∗,yr=1,forr∈−1,0,yr=2,forr∈3,4.
We proved that the error is of orderON−1lnN and ON−2ln2N. The theory has been validated with two examples; referring to these numerical results, it can be observed that the proposed method has been effective and applicable.


## 8. Discussion

In the literature, many authors have considered singular perturbation problem mixed delayτ≪1 differential equation. In this paper, we consider a singular perturbation problem with mixed delay τ=1 differential equation. We suggested two computational methods such as finite and hybrid difference scheme. We proved that the error is of order ON−1lnN and ON−2ln2N. Finally, two numerical examples are also presented to validate the theoretical results of this study. Maximum pointwise errors and order of convergence of Examples 1 and 2 are given in Tables 1 and 2, respectively.Table 1
ComputedPN rate of convergence and DN maximum errors for Example 1.
Number of mesh points 3Nε1632641282565121024Finite difference method10−32.8530e − 031.6423e − 038.6257e − 044.3424e − 042.1521e − 041.0694e − 045.4208e − 0510−43.0817e − 031.8315e − 039.9232e − 045.1418e − 042.6100e − 041.3150e − 046.6324e − 0510−53.1526e − 031.8900e − 031.0321e − 035.3860e − 042.7489e − 041.3887e − 046.9901e − 0510−63.1750e − 031.9083e − 031.0446e − 035.4624e − 042.7922e − 041.4116e − 047.1008e − 0510−73.1820e − 031.9141e − 031.0486e − 035.4864e − 042.8058e − 041.4188e − 047.1355e − 0510−83.1842e − 031.9159e − 031.0498e − 035.4940e − 042.8101e − 041.4211e − 047.1465e − 0510−93.1849e−031.9165e−031.0502e − 035.4965e − 042.8115e− 041.4218e − 047.1499e − 0510−103.1851e − 031.9167e − 031.0503e − 035.4972e − 042.8119e − 041.4220e − 047.1510e − 05DN3.1851e − 031.9167e − 031.0503e − 035.4972e − 042.8119e − 041.4220e − 047.1510e − 05PN7.327520e − 018.6771e − 019.3412e − 019.6712e − 019.8357e − 019.9178e − 01—Hybrid difference method10−31.6679e − 026.7615e − 032.4201e − 038.5271e − 042.6786e − 047.8677e − 052.6259e − 0510−41.6684e − 026.7758e − 032.4331e − 038.6109e − 042.7258e − 048.3394e − 052.5535e − 0510−51.6684e − 026.7772e − 032.4344e−038.6194e − 042.7307e − 048.4126e − 052.6113e − 0510−61.6684e − 026.7774e − 032.4345e−038.6202e − 042.7312e − 048.4200e − 052.6192e − 0510−71.6684e − 026.7774e − 032.4345e − 038.6203e − 042.7312e − 048.4207e − 052.6200e − 0510−81.6684e − 026.7774e − 032.4345e − 038.6203e − 042.7312e − 048.4208e − 052.6201e − 0510−91.6684e − 026.7774e − 032.4345e − 038.6204e − 042.7312e − 048.4209e − 052.6200e − 0510−101.6685e − 026.7774e − 032.4346e − 038.6201e − 042.7308e − 048.4224e − 052.6206e − 05DN1.6685e − 026.7774e − 032.4346e − 038.6204e − 042.7312e − 048.4224e − 052.6259e − 05PN1.2997e + 001.4770e + 001.4979e + 001.6581e + 001.6972e + 001.6814e + 00—Table 2
ComputedPN rate of convergence and DN maximum errors for Example 2.
Number of mesh points 3Nε1632641282565121024Finite difference method10−35.5384e − 032.5272e − 031.1887e − 035.6672e − 042.7241e − 041.3205e − 046.4844e − 0510−45.8841e − 032.7585e − 031.3350e − 036.5461e − 043.2299e − 041.6003e − 047.9627e − 0510−55.9911e − 032.8296e − 031.3798e − 036.8149e − 043.3841e − 041.6851e − 048.4085e − 0510−66.0248e − 032.8519e − 031.3939e − 036.8989e − 043.4321e − 041.7116e − 048.5469e − 0510−76.0354e − 032.8590e − 031.3983e − 036.9253e − 043.4473e − 041.7199e − 048.5904e − 0510−86.0387e − 032.8612e − 031.3997e − 036.9337e − 043.4521e − 041.7225e − 048.6041e − 0510−96.0398e − 032.8619e − 031.4001e − 036.9363e − 043.4536e − 041.7233e − 048.6085e − 0510−106.0401e − 032.8621e − 031.4003e − 036.9371e − 043.4541e − 041.7236e − 048.6098e − 05DN6.0401e − 032.8621e − 031.4003e − 036.9371e − 043.4541e − 041.7236e − 048.6098e − 05PN1.0774e + 001.0313e + 001.0133e + 001.0060e + 001.0028e + 001.0013e + 00—Hybrid difference method10−32.8953e − 021.1470e − 023.9972e − 031.4195e − 034.4816e − 041.3262e − 043.963906e − 0510−42.8954e − 021.1483e − 024.0105e − 031.4282e − 034.5308e − 041.3799e − 044.255447e − 0510−52.8954e − 021.1484e − 024.0119e − 031.4291e − 034.5358e − 041.3871e − 044.303554e − 0510−62.8954e − 021.1485e − 024.0120e − 031.4292e − 034.5363e − 041.3879e − 044.308411e − 0510−72.8954e − 021.1485e − 024.0120e − 031.4292e − 034.5364e − 041.3879e − 044.308899e − 0510−82.8954e − 021.1485e − 024.0120e − 031.4292e − 034.5364e − 041.3879e − 044.308959e − 0510−92.8954e − 021.1485e − 024.0120e − 031.4292e − 034.5363e − 041.3880e − 044.309005e − 0510−102.8954e − 021.1485e − 024.0121e − 031.4292e − 034.5358e − 041.3882e − 044.309514e − 05DN2.8954e − 021.1485e − 024.0121e − 031.4292e − 034.5364e − 041.3882e − 044.309514e − 05PN1.3340e + 001.5173e + 001.4891e + 001.6556e + 001.7082e + 001.6876e + 00—



---

*Source: 6636607-2021-06-04.xml*