# Higher-Order Uniformly Convergent Numerical Scheme for Singularly Perturbed Differential Difference Equations with Mixed Small Shifts

Authors: Mesfin Mekuria Woldaregay; Gemechis File Duressa
Journal: International Journal of Differential Equations
 (2020)

---

## Abstract

This paper deals with numerical treatment of singularly perturbed differential difference equations involving mixed small shifts on the reaction terms. The highest-order derivative term in the equation is multiplied by a small perturbation parameterε taking arbitrary values in the interval 0,1. For small values of ε, the solution of the problem exhibits exponential boundary layer on the left or right side of the domain and the derivatives of the solution behave boundlessly large. The terms having the shifts are treated using Taylor’s series approximation. The resulting singularly perturbed boundary value problem is solved using exponentially fitted operator FDM. Uniform stability of the scheme is investigated and analysed using comparison principle and solution bound. The formulated scheme converges uniformly with linear order before Richardson extrapolation and quadratic order after Richardson extrapolation. The theoretical analysis of the scheme is validated using numerical test examples for different values of ε and mesh number N.


---

## Body

## 1. Introduction

Differential equations play a prominent role in many disciplines including engineering, physics, economics, and biology. Currently different authors are working on analytical and numerical solutions of differential equations using different techniques [1, 2]. Differential difference equations (DDEs) are differential equations where the evolution of the system not only depends on the present state of the system but also depends on the past history. Singularly perturbed differential difference equations are differential equations in which the highest-order derivative term is multiplied by a small perturbation parameter ε and involves at least one term with delay. In general, when the perturbation parameter tends to zero, the smoothness of the solution of the singularly perturbed differential difference equations (SPDDEs) deteriorates and it forms boundary layer [3]. Such type of equations has applications in the study of variational problems in control theory [4] and in modelling of neuronal variability [5].The presence of singular perturbation parameterε in the equation leads to oscillation in the computed solution, while using standard numerical methods like FDM, FEM, and spline method [6]. To avoid this oscillation, an unacceptably large number of mesh points are required when ε is very small. This is not practical and leads to round-off error. So, to overcome the drawbacks associated with standard numerical methods, different authors have developed schemes that converge uniformly.Numerical treatments of a class of SPDDEs have received a great deal of attention recently because of their wide applications. It is of theoretical and practical interest to consider numerical methods for such problems. Owing to this, here we present some prior studies on numerical solution of the considered problem. Lange and Miura in [7–10] studied a class of second-order DDEs in which the second derivative term is multiplied by a small parameter. The authors extend the method of matched asymptotic expansions initially developed for solving boundary value problems to obtain approximate solution for SPDDEs. In a series of papers [11–14], Kadalbajoo and Sharma developed uniformly convergent numerical methods using fitted mesh FDMs techniques. Swamy et al. [15–17] considered the problem and developed a numerical scheme using fitted operator finite difference techniques. Melesse et al. [18] applied initial value technique for treating the considered SPDDEs. Ranjan and Prasad [19] used modified fitted operator FDM for solving the problem. Sirisha et al. [20] developed fitted operator finite difference scheme using the procedure of domain decomposition. A number of authors have developed numerical scheme using exponentially fitted method for solving SPDDEs. To the authors’ knowledge, none of them show the uniform convergence of their schemes. This motivates to treat the considered SPDDEs and formulate the uniform convergence analysis of the scheme. Our contribution in this paper is to develop higher-order uniformly convergent numerical scheme using exponentially fitted FDM and to analyse the uniform convergence of the proposed scheme.Notation 1.
The symbolC is used to denote positive constant independent of ε and N. The norm . denotes the maximum norm.


## 2. Statement of the Problem

A class of second-order singularly perturbed differential difference equations having mixed shift on reaction terms is given by(1)−εu″x+axu′x+αxux−δ+ωxux+βxux+η=fx,x∈Ω=0,1,with the interval conditions(2)ux=ϕx,x∈−δ,0,ux=ψx,x∈1,1+η,where ε∈0,1 is singular perturbation parameter and δ and η are delay and advance parameters satisfying δ,η<ε. The functions ax,αx,ωx,βx,fx,ϕx, and ψx are assumed to be sufficiently smooth, bounded, and independent of ε for guaranteeing the existence of unique solution. The coefficient functions αx,ωx, and βx are assumed to satisfy(3)αx+ωx+βx≥q∗>0,∀x∈Ω¯,for some constant q∗.For the case whereδ,η<ε, using Taylor’s series approximation for the terms with shift is appropriate [21]. So, we approximate ux−δ and ux+η as(4)ux−δ≈ux−δu′x+δ22u″x+Oδ3,ux+η≈ux+ηu′x+η22u″x+Oη3.Using the approximation in (4) into (1)-(2) gives(5)−cεu″x+pxu′x+qxux=fx,x∈0,1,with the boundary conditions(6)u0=ϕ0,u1=ψ1,where cε=ε2−δ2/2α−η2/2β, px=ax−δαx+ηβx, and qx=αx+βx+ωx for α and β are lower bounds of αx and βx, respectively. For small values of δ,η, equations (1)-(2) and (5)-(6) are asymptotically equivalent, since the truncated term in (5)-(6) is order of Oδ3,η3.Let us denote differential operatorL for the differential equation in (5) as(7)Lux=−cεu″x+pxu′x+qxux.The solution of the problem in (5)-(6) exhibits regular boundary layer of thickness Ocε and the position of the boundary layer depends on the sign of px. If px<0, left boundary layer exists and, for px>0, right boundary layer exists. If px,x∈0,1 changes sign interior layer will occur [22].Settingcε=0 in equations (5)-(6) gives the reduced problem. For the case where px>0, it is given by(8)pxu0′x+qxu0x=fx,∀x∈Ω,u00=ϕ0,and, for the case where px<0, it is given by(9)pxu0′x+qxu0x=fx,∀x∈Ω,u01=ψ1.It is a first-order initial value problem; for small values ofcε, the solution of (5)-(6) is very close to the solution of (8) or (9).


### 2.1. Properties of the Analytical Solution

Lemma 1 (The maximum principle).
Letz be a sufficiently smooth function defined on Ω, which satisfies z0≥0 and z1≥0. Then Lzx>0,∀x∈Ω, implies that zx≥0,∀x∈Ω¯.Proof.
Letx∗∈Ω¯ be such that zx∗=minx∈Ω¯zx and suppose that zx∗<0. It is clear that x∗∉0,1. Since zx∗=minx∈Ω¯zx form extrema values in calculus we have z′x∗=0 and z″x∗≥0 giving that Lzx∗<0, which is contradiction to the assumption made above: Lzx∗>0,∀x∈Ω. Therefore, zx≥0,∀x∈Ω.Lemma 2 (Stability estimate).
Letu be the solution of (5)-(6); then, it satisfies the bound(10)ux≤fq∗+maxϕ0,ψ1,where q∗ is lower bound of qx.Proof.
Let us define barrier functionsϑ±x as ϑ±x=f/q∗+maxϕ0,ψ1±ux and apply the maximum principle to obtain the required bound.
On the boundary points, we have(11)ϑ±0=fq∗+maxϕ0,ψ1±u0≥0,ϑ±1=fq∗+maxϕ0,ψ1±u1≥0.
On the differential operator,(12)Lϑ±x=−cεϑ±″x+pxϑ±′x+qxϑ±x=−cε0±u″x+px0±u′x+qxfq∗+maxϕ0,ψ1±ux=qxfq∗+maxϕ0,ψ1±fx≥0,sinceqx≥q∗>0,which implies that Lϑ±x≥0. Hence, using the maximum principle in Lemma 1, we obtain ϑ±x≥0,∀x∈Ω¯.Lemma 3.
The derivatives of the solutions of (5)-(6) satisfy the bound(13)ukx≤C1+cε−kexp−p∗xcε,x∈Ω¯,0≤k≤4,for left boundary layer problems and(14)ukx≤C1+cε−kexp−p∗1−xcε,x∈Ω¯,0≤k≤4,for right boundary layer problems.Proof.
See [23].


## 2.1. Properties of the Analytical Solution

Lemma 1 (The maximum principle).
Letz be a sufficiently smooth function defined on Ω, which satisfies z0≥0 and z1≥0. Then Lzx>0,∀x∈Ω, implies that zx≥0,∀x∈Ω¯.Proof.
Letx∗∈Ω¯ be such that zx∗=minx∈Ω¯zx and suppose that zx∗<0. It is clear that x∗∉0,1. Since zx∗=minx∈Ω¯zx form extrema values in calculus we have z′x∗=0 and z″x∗≥0 giving that Lzx∗<0, which is contradiction to the assumption made above: Lzx∗>0,∀x∈Ω. Therefore, zx≥0,∀x∈Ω.Lemma 2 (Stability estimate).
Letu be the solution of (5)-(6); then, it satisfies the bound(10)ux≤fq∗+maxϕ0,ψ1,where q∗ is lower bound of qx.Proof.
Let us define barrier functionsϑ±x as ϑ±x=f/q∗+maxϕ0,ψ1±ux and apply the maximum principle to obtain the required bound.
On the boundary points, we have(11)ϑ±0=fq∗+maxϕ0,ψ1±u0≥0,ϑ±1=fq∗+maxϕ0,ψ1±u1≥0.
On the differential operator,(12)Lϑ±x=−cεϑ±″x+pxϑ±′x+qxϑ±x=−cε0±u″x+px0±u′x+qxfq∗+maxϕ0,ψ1±ux=qxfq∗+maxϕ0,ψ1±fx≥0,sinceqx≥q∗>0,which implies that Lϑ±x≥0. Hence, using the maximum principle in Lemma 1, we obtain ϑ±x≥0,∀x∈Ω¯.Lemma 3.
The derivatives of the solutions of (5)-(6) satisfy the bound(13)ukx≤C1+cε−kexp−p∗xcε,x∈Ω¯,0≤k≤4,for left boundary layer problems and(14)ukx≤C1+cε−kexp−p∗1−xcε,x∈Ω¯,0≤k≤4,for right boundary layer problems.Proof.
See [23].


## 3. Numerical Scheme

First, let us discretize the domainΩ¯=0,1 into N equal number of subintervals with mesh length h=1/N as Ω¯N=xi=ih0N. Let ux be smooth function on the domain Ω¯=0,1; then, using Taylor series approximation, we have(15)uxi+1=ui+1=ui+hui′+h22!ui″+h33!ui3+h44!ui4+h55!ui5+h66!ui6+Oh7,uxi−1=ui−1=ui−hui′+h22!ui″−h33!ui3+h44!ui4−h55!ui5+h66!ui6+Oh7.Taking the difference in (15), we obtain(16)ui−1−2ui+ui+1=2h22!ui″+2h44!ui4+2h66!ui6+Oh8.Differentiating (16) two times gives(17)ui−1″−2ui″+ui+1″=2h22!ui4+2h44!ui6+2h66!ui8+Oh8.Now, multiplying (17) by −h2/12 and adding with (16) to eliminate the term with ui4 gives(18)ui−1−2ui+ui+1=h212ui−1″+10ui′+ui+1″+τ,where τ=Oh6. Evaluating (5) at xi−1,xi, and xi+1, respectively, we obtain(19)−cεui−1″=−pxi−1ui−1′−qxi−1ui−1+fxi−1,−cεui″=−pxiui′−qxiui+fxi,−cεui+1″=−pxi+1ui+1′−qxi+1ui+1+fxi+1.Next, approximate the first derivative termsui−1′, ui′, and ui−1′ in (19), using the right shifted, central, and left shifted finite difference approximations as(20)ui−1′=ui−1−4ui+3ui+12h+hui″+Oh2,ui′=ui+1−ui−12h+Oh2,ui+1′=−3ui−1+4ui−ui+12h−hui″+Oh2.Substituting (20) into (19) and then (19) into (18) gives(21)−cε−hpi−112+hpi+112Ui−1−2Ui+Ui+1h2+pi−124h−3Ui−1+4Ui−Ui+1+10pi24hUi+1−Ui−1+pi+124hUi−1−4Ui+3Ui+1+qi−112Ui−1+10qi12Ui+qi+112Ui+1=112fi−1+10fi+fi+1+Oh2,i=1,2,…,N−1,where pi−1,pi, and pi+1 are denoted for pxi−1,pxi, and pxi+1, respectively. We denote Ui for the approximate solution of uxi in the above discretization.To get small truncation error in boundary layer region, we apply exponentially fitted operator finite difference method (FOFDM). For developing the FOFDM, we use the theory developed in asymptotic method for treating singularly perturbed BVPs. Let us consider and treat the left and the right boundary layer cases separately.


### 3.1. Left Boundary Layer Problems

In this case, the boundary layer occurs nearx=0. From the theory of singular perturbation given in [24], the zeroth-order asymptotic solution of (5)-(6) is given by(22)ux=u0x+p0pxϕ0−u00exp−∫0xpxcε−qxpxdx+Ocε,where u0 is the solution of the reduced problem. Using Taylor’s series approximation for u0x,px, and qx centred at xi=ih up to first order and considering cε⟶0, the discretized form of (22) becomes(23)uih=u0ih+ϕ0−u00exp−ρpi,where ρ=h/cε,h=1/N. Similarly, we write(24)ui+1h=u0ih+ϕ0−u00exp−ρpii+1,ui−1h=u0ih+ϕ0−u00exp−ρpii−1.To handle the effect of the perturbation parameter, exponential fitting factorσ1ρ is multiplied on the term containing the perturbation parameter as(25)−σ1ρcε−hpi−112+hpi+112Ui−1−2Ui+Ui+1h2+pi−124h−3Ui−1+4Ui−Ui+1+10pi24hUi+1−Ui−1+pi+124hUi−1−4Ui+3Ui+1+qi−112Ui−1+10qi12Ui+qi+112Ui+1=112fi−1+10fi+fi+1,i=1,2,…,N−1.Multiplying both sides of (25) by h and substituting ρ for h/cε and taking the limit as h⟶0 give(26)−limh⟶0σ1ρρUi−1−2Ui+Ui+1+p024limh⟶0−3Ui−1+4Ui−Ui+1+10p024limh⟶0Ui+1−Ui−1+p024limh⟶0Ui−1−4Ui+3Ui+1=0.From (23) and (24), we obtain(27)limh⟶0Ui−1−2Ui+Ui+1=ϕ0−u00e−p0iρep0ρ−2+e−p0ρ,limh⟶0−3Ui−1+4Ui−Ui+1=ϕ0−u00e−p0iρ−3ep0ρ+4−e−p0ρ,limh⟶0Ui−1−4Ui+3Ui+1=ϕ0−u00e−p0iρep0ρ−4+3e−p0ρ,limh⟶0Ui+1−Ui−1=ϕ0−u00e−p0iρe−p0ρ−ep0ρ.Substituting (27) into (26) and simplifying give(28)σ1ρρep0ρ−2+e−p0ρ=p024−3ep0ρ−e−p0ρ+10e−p0ρ−ep0ρ+ep0ρ+3e−p0ρ.The exponential fitting factor is obtained as(29)σ1ρ=ρp02cothρp02.Hence, the required finite difference scheme becomes(30)LLhUi=112fi−1+10fi+fi+1,i=1,2,…,N−1,where(31)LLhUi=−σ1ρcε−hpi−112+hpi+112Ui−1−2Ui+Ui+1h2+pi−124h−3Ui−1+4Ui−Ui+1+10pi24hUi+1−Ui−1+pi+124hUi−1−4Ui+3Ui+1+qi−112Ui−1+10qi12Ui+qi+112Ui+1,with the boundary values U0=ϕ0 and UN=ψ1.


### 3.2. Right Boundary Layer Problems

In this case, the boundary layer occurs nearx=1. From [24], the zeroth-order asymptotic solution of (5)-(6) is given by(32)ux=u0x+p1pxψ1−u01exp−∫x1pxcε−qxpxdx+Ocε,where u0 is the solution of the reduced problem.Using Taylor’s series approximation foru0x,px, and qx centred at xi=ih up to first order and considering cε⟶0, the discretized form of (24) becomes(33)uih=u0ih+ϕ0−u00exp−pi1cε−ρ.Using similar procedures as the left boundary layer case, the exponential fitting factor is obtained as(34)σ2ρ=ρp12cothρp12.Hence, the required finite difference scheme becomes(35)LRhUi=112fi−1+10fi+fi+1,i=1,2,…,N−1,where(36)LRhUi=−σ2ρcε−hpi−112+hpi+112Ui−1−2Ui+Ui+1h2+pi−124h−3Ui−1+4Ui−Ui+1+10pi24hUi+1−Ui−1+pi+124hUi−1−4Ui+3Ui+1+qi−112Ui−1+10qi12Ui+qi+112Ui+1.


### 3.3. Convergence Analysis

In this section, we show the stability and convergence analysis for the right boundary layer problems. In similar manner, it is proved for the left boundary layer case. First, we need to prove the discrete comparison principle for the scheme in (35) for guaranteeing existence of unique discrete solution.Lemma 4 (Discrete comparison principle).
Assume that, for mesh functionUi there exists a comparison function Vi such that LRhUi≤LRhVi,i=1,2,…,N−1 and if U0≤V0 and UN≤VN, then Ui≤Vi,∀i,i=0,1,2,…,N.Proof.
The matrix associated with operatorLRh is of size N+1×N+1 and satisfies the property of M-matrix. See the detailed proof in [23].
This lemma gives guarantee for the existence of unique discrete solution. In the next lemma, we discuss the uniform stability of the discrete solution.Lemma 5 (Discrete stability estimate).
The solutionUi of the discrete scheme in (35) satisfies the following bound:(37)Ui≤LRhUiq∗+maxU0,UN.Proof.
LetQ=LRhUi/q∗+maxU0,UN and define barrier functions ϑi± by ϑi±=Q±Ui.
On the boundary points, we obtain(38)ϑ0±=Q±U0=LRhUiq∗+maxU0,UN±ϕ0≥0,ϑN±=Q±UN=LRhUiq∗+maxU0,UN±ψ1≥0.
On the discretized spatial domainxi,1<i<N−1, we obtain(39)LRhϑi±=−cεσρ−hpi−112+hpi+112Q±Ui+1−2Q±Ui+Q±Ui−1h2+pi−124h−3Q+Ui−1+4Q+Ui−Q+Ui+1+10pi24hQ+Ui+1−Q+Ui−1+pi+124hQ+Ui−1−4Q+Ui+3Q+Ui+1+qi−112Q+Ui−1+10qi12Q+Ui+qi+112Q+Ui+1=∓cεσρ−hpi−112+hpi+112Ui+1−2Ui+Ui−1h2±pi−124h−3Ui−1+4Ui−Ui+1±10pi24hUi+1−Ui−1±pi+124hUi−1−4Ui+3Ui+1±qi−112Ui−1±10qi12Ui±qi+112Ui+1+qi−112+10qi12+qi+112Q=qi−112+10qi12+qi+112LRhUiq∗+maxU0,UN±112fi−1+10fi+fi+1≥0,sinceqxi≥q∗>0.
By the discrete comparison principle in Lemma4, we obtain ϑi±≥0,∀xi∈Ω¯N. Hence, the required bound is satisfied.
Now, let us denote the right shifted, centred, and left shifted finite differences, respectively, as(40)ΔRuxi=ui−1−4ui+3ui+12h+hui″,ΔCuxi=ui+1−ui−12h,ΔLuxi=−3ui−1+4ui−ui+12h−hui″.
Using Taylor’s series approximation, we obtain the bound(41)ui−1′−ΔLuxi≤CN−2u4xi,ui′−ΔCuxi≤CN−2u″xi,ui+1′−ΔRuxi≤CN−2u4xi,where ukxi=maxx0≤xi≤xNukxi,k=2,4. Similarly, we have(42)u″xi−D+D−uxi≤CN−2u4xi,D+D−uxi≤Cu″xi.
Now, forρ>0, C1 and C2 are constants, and we have ρcothρ−1≤C1ρ2, for ρ≤1. For ρ⟶∞, since limρ⟶∞cothρ=1, ρcothρ−1≤C1ρ is given. In general, for all ρ>0, we write(43)C1ρ2ρ+1≤ρcothρ−1≤C2ρ2ρ+1,implying that(44)cεpxiρ2cothp1ρ2−1≤cεN−1/cε2N−1/cε+1=N−2N−1+cε.
The following theorem gives truncation error bound of the proposed scheme.Theorem 1.
Letuxi and Ui be solutions of (5)-(6) and (35), respectively. Then the following error estimate holds:(45)LRhuxi−Ui≤CN−2N−1+cε1+cε−4maxx1≤xi≤xN−1exp−p∗1−xicε.Proof.
Consider the truncation error that is given by(46)Lhuxi−LRhUi≤−cεu″xi−σρD+D−uxi+pi−112ui−1′−ΔLuxi+10pi12ui′−ΔCuxi+pi+112ui+1′−ΔRuxi≤−cεp1ρ2cothp1ρ2−1D+D−uxi+cεu″xi−D+D−uxi+pi−112ui−1′−ΔLuxi+10pi12ui′−ΔCuxi+pi+112ui+1′−ΔRuxi.
Using the bounds in (44), (41), and (42) gives(47)Lhuxi−Ui≤CN−2N−1+cεu″xi+cεCN−2u4xi+CN−2u4xi≤CN−2N−1+cεu″xi+CN−2u4xi.
Substituting the bounds for the derivatives of the solution in Lemma3, we obtain(48)Lhuxi−Ui≤CN−2N−1+cε1+cε−2exp−p∗1−xicε+CN−21+cε−4exp−p∗1−xicε.
Sincecε−2≤cε−4, we obtain(49)Lhuxi−Ui≤CN−2N−1+cε1+cε−4maxx1≤xi≤xN−1exp−p∗1−xicε.Lemma 6.
Forcε⟶0, and for given fixed N, we obtain(50)limcε⟶0maxxjexp−p∗xj/cεcεm=0,limcε⟶0maxxjexp−p∗xj/cεcεm=0,m=1,2,3,…,where xj=jh,h=1/N,∀j=1,2,…,N−1.Proof.
See [25–27].Theorem 2.
Under the hypothesis of boundedness of discrete solution, the solution of the discrete schemes in (30) satisfies the following uniform error bound:(51)sup0<cε≪1uxi−Ui≤CN−1.Proof.
Substituting the results in Lemma6 into Theorem 1, applying Lemma 5 gives the required bound.


### 3.4. Richardson Extrapolation

We apply the Richardson extrapolation technique to accelerate the rate of convergence of the proposed scheme. Richardson extrapolation is a convergence acceleration technique that involves combination of two computed approximations of solution. Interested reader can see the details of Richardson extrapolation in [28]. From (34) and Lemma 6, we obtain(52)LRhuxi−Ui≤CN−1+N−2,where uxi and Ui are the exact and approximate solutions of (5)-(6), respectively. Applying Lemma 5 in (52) gives(53)uxi−Ui≤CN−1+ON−2.Let us denoteUi2N for an approximate solution on 2N number of mesh points by including the mid-points. Applying the same procedure as in (52) and (53) on 2N number of mesh points, we obtain(54)uxi−Ui2N≤C2N−1+ON−2.Combining (53) and (54) for removing the term CN−1 results in(55)uxi−2Ui2N−Ui≤CN−2,giving (56)Uiext=2Ui2N−Ui,as the extrapolated solution. The error bound for the extrapolated solution in (56) becomes(57)sup0<cε≪1uxi−Uiext≤CN−2.


## 3.1. Left Boundary Layer Problems

In this case, the boundary layer occurs nearx=0. From the theory of singular perturbation given in [24], the zeroth-order asymptotic solution of (5)-(6) is given by(22)ux=u0x+p0pxϕ0−u00exp−∫0xpxcε−qxpxdx+Ocε,where u0 is the solution of the reduced problem. Using Taylor’s series approximation for u0x,px, and qx centred at xi=ih up to first order and considering cε⟶0, the discretized form of (22) becomes(23)uih=u0ih+ϕ0−u00exp−ρpi,where ρ=h/cε,h=1/N. Similarly, we write(24)ui+1h=u0ih+ϕ0−u00exp−ρpii+1,ui−1h=u0ih+ϕ0−u00exp−ρpii−1.To handle the effect of the perturbation parameter, exponential fitting factorσ1ρ is multiplied on the term containing the perturbation parameter as(25)−σ1ρcε−hpi−112+hpi+112Ui−1−2Ui+Ui+1h2+pi−124h−3Ui−1+4Ui−Ui+1+10pi24hUi+1−Ui−1+pi+124hUi−1−4Ui+3Ui+1+qi−112Ui−1+10qi12Ui+qi+112Ui+1=112fi−1+10fi+fi+1,i=1,2,…,N−1.Multiplying both sides of (25) by h and substituting ρ for h/cε and taking the limit as h⟶0 give(26)−limh⟶0σ1ρρUi−1−2Ui+Ui+1+p024limh⟶0−3Ui−1+4Ui−Ui+1+10p024limh⟶0Ui+1−Ui−1+p024limh⟶0Ui−1−4Ui+3Ui+1=0.From (23) and (24), we obtain(27)limh⟶0Ui−1−2Ui+Ui+1=ϕ0−u00e−p0iρep0ρ−2+e−p0ρ,limh⟶0−3Ui−1+4Ui−Ui+1=ϕ0−u00e−p0iρ−3ep0ρ+4−e−p0ρ,limh⟶0Ui−1−4Ui+3Ui+1=ϕ0−u00e−p0iρep0ρ−4+3e−p0ρ,limh⟶0Ui+1−Ui−1=ϕ0−u00e−p0iρe−p0ρ−ep0ρ.Substituting (27) into (26) and simplifying give(28)σ1ρρep0ρ−2+e−p0ρ=p024−3ep0ρ−e−p0ρ+10e−p0ρ−ep0ρ+ep0ρ+3e−p0ρ.The exponential fitting factor is obtained as(29)σ1ρ=ρp02cothρp02.Hence, the required finite difference scheme becomes(30)LLhUi=112fi−1+10fi+fi+1,i=1,2,…,N−1,where(31)LLhUi=−σ1ρcε−hpi−112+hpi+112Ui−1−2Ui+Ui+1h2+pi−124h−3Ui−1+4Ui−Ui+1+10pi24hUi+1−Ui−1+pi+124hUi−1−4Ui+3Ui+1+qi−112Ui−1+10qi12Ui+qi+112Ui+1,with the boundary values U0=ϕ0 and UN=ψ1.


## 3.2. Right Boundary Layer Problems

In this case, the boundary layer occurs nearx=1. From [24], the zeroth-order asymptotic solution of (5)-(6) is given by(32)ux=u0x+p1pxψ1−u01exp−∫x1pxcε−qxpxdx+Ocε,where u0 is the solution of the reduced problem.Using Taylor’s series approximation foru0x,px, and qx centred at xi=ih up to first order and considering cε⟶0, the discretized form of (24) becomes(33)uih=u0ih+ϕ0−u00exp−pi1cε−ρ.Using similar procedures as the left boundary layer case, the exponential fitting factor is obtained as(34)σ2ρ=ρp12cothρp12.Hence, the required finite difference scheme becomes(35)LRhUi=112fi−1+10fi+fi+1,i=1,2,…,N−1,where(36)LRhUi=−σ2ρcε−hpi−112+hpi+112Ui−1−2Ui+Ui+1h2+pi−124h−3Ui−1+4Ui−Ui+1+10pi24hUi+1−Ui−1+pi+124hUi−1−4Ui+3Ui+1+qi−112Ui−1+10qi12Ui+qi+112Ui+1.


## 3.3. Convergence Analysis

In this section, we show the stability and convergence analysis for the right boundary layer problems. In similar manner, it is proved for the left boundary layer case. First, we need to prove the discrete comparison principle for the scheme in (35) for guaranteeing existence of unique discrete solution.Lemma 4 (Discrete comparison principle).
Assume that, for mesh functionUi there exists a comparison function Vi such that LRhUi≤LRhVi,i=1,2,…,N−1 and if U0≤V0 and UN≤VN, then Ui≤Vi,∀i,i=0,1,2,…,N.Proof.
The matrix associated with operatorLRh is of size N+1×N+1 and satisfies the property of M-matrix. See the detailed proof in [23].
This lemma gives guarantee for the existence of unique discrete solution. In the next lemma, we discuss the uniform stability of the discrete solution.Lemma 5 (Discrete stability estimate).
The solutionUi of the discrete scheme in (35) satisfies the following bound:(37)Ui≤LRhUiq∗+maxU0,UN.Proof.
LetQ=LRhUi/q∗+maxU0,UN and define barrier functions ϑi± by ϑi±=Q±Ui.
On the boundary points, we obtain(38)ϑ0±=Q±U0=LRhUiq∗+maxU0,UN±ϕ0≥0,ϑN±=Q±UN=LRhUiq∗+maxU0,UN±ψ1≥0.
On the discretized spatial domainxi,1<i<N−1, we obtain(39)LRhϑi±=−cεσρ−hpi−112+hpi+112Q±Ui+1−2Q±Ui+Q±Ui−1h2+pi−124h−3Q+Ui−1+4Q+Ui−Q+Ui+1+10pi24hQ+Ui+1−Q+Ui−1+pi+124hQ+Ui−1−4Q+Ui+3Q+Ui+1+qi−112Q+Ui−1+10qi12Q+Ui+qi+112Q+Ui+1=∓cεσρ−hpi−112+hpi+112Ui+1−2Ui+Ui−1h2±pi−124h−3Ui−1+4Ui−Ui+1±10pi24hUi+1−Ui−1±pi+124hUi−1−4Ui+3Ui+1±qi−112Ui−1±10qi12Ui±qi+112Ui+1+qi−112+10qi12+qi+112Q=qi−112+10qi12+qi+112LRhUiq∗+maxU0,UN±112fi−1+10fi+fi+1≥0,sinceqxi≥q∗>0.
By the discrete comparison principle in Lemma4, we obtain ϑi±≥0,∀xi∈Ω¯N. Hence, the required bound is satisfied.
Now, let us denote the right shifted, centred, and left shifted finite differences, respectively, as(40)ΔRuxi=ui−1−4ui+3ui+12h+hui″,ΔCuxi=ui+1−ui−12h,ΔLuxi=−3ui−1+4ui−ui+12h−hui″.
Using Taylor’s series approximation, we obtain the bound(41)ui−1′−ΔLuxi≤CN−2u4xi,ui′−ΔCuxi≤CN−2u″xi,ui+1′−ΔRuxi≤CN−2u4xi,where ukxi=maxx0≤xi≤xNukxi,k=2,4. Similarly, we have(42)u″xi−D+D−uxi≤CN−2u4xi,D+D−uxi≤Cu″xi.
Now, forρ>0, C1 and C2 are constants, and we have ρcothρ−1≤C1ρ2, for ρ≤1. For ρ⟶∞, since limρ⟶∞cothρ=1, ρcothρ−1≤C1ρ is given. In general, for all ρ>0, we write(43)C1ρ2ρ+1≤ρcothρ−1≤C2ρ2ρ+1,implying that(44)cεpxiρ2cothp1ρ2−1≤cεN−1/cε2N−1/cε+1=N−2N−1+cε.
The following theorem gives truncation error bound of the proposed scheme.Theorem 1.
Letuxi and Ui be solutions of (5)-(6) and (35), respectively. Then the following error estimate holds:(45)LRhuxi−Ui≤CN−2N−1+cε1+cε−4maxx1≤xi≤xN−1exp−p∗1−xicε.Proof.
Consider the truncation error that is given by(46)Lhuxi−LRhUi≤−cεu″xi−σρD+D−uxi+pi−112ui−1′−ΔLuxi+10pi12ui′−ΔCuxi+pi+112ui+1′−ΔRuxi≤−cεp1ρ2cothp1ρ2−1D+D−uxi+cεu″xi−D+D−uxi+pi−112ui−1′−ΔLuxi+10pi12ui′−ΔCuxi+pi+112ui+1′−ΔRuxi.
Using the bounds in (44), (41), and (42) gives(47)Lhuxi−Ui≤CN−2N−1+cεu″xi+cεCN−2u4xi+CN−2u4xi≤CN−2N−1+cεu″xi+CN−2u4xi.
Substituting the bounds for the derivatives of the solution in Lemma3, we obtain(48)Lhuxi−Ui≤CN−2N−1+cε1+cε−2exp−p∗1−xicε+CN−21+cε−4exp−p∗1−xicε.
Sincecε−2≤cε−4, we obtain(49)Lhuxi−Ui≤CN−2N−1+cε1+cε−4maxx1≤xi≤xN−1exp−p∗1−xicε.Lemma 6.
Forcε⟶0, and for given fixed N, we obtain(50)limcε⟶0maxxjexp−p∗xj/cεcεm=0,limcε⟶0maxxjexp−p∗xj/cεcεm=0,m=1,2,3,…,where xj=jh,h=1/N,∀j=1,2,…,N−1.Proof.
See [25–27].Theorem 2.
Under the hypothesis of boundedness of discrete solution, the solution of the discrete schemes in (30) satisfies the following uniform error bound:(51)sup0<cε≪1uxi−Ui≤CN−1.Proof.
Substituting the results in Lemma6 into Theorem 1, applying Lemma 5 gives the required bound.


## 3.4. Richardson Extrapolation

We apply the Richardson extrapolation technique to accelerate the rate of convergence of the proposed scheme. Richardson extrapolation is a convergence acceleration technique that involves combination of two computed approximations of solution. Interested reader can see the details of Richardson extrapolation in [28]. From (34) and Lemma 6, we obtain(52)LRhuxi−Ui≤CN−1+N−2,where uxi and Ui are the exact and approximate solutions of (5)-(6), respectively. Applying Lemma 5 in (52) gives(53)uxi−Ui≤CN−1+ON−2.Let us denoteUi2N for an approximate solution on 2N number of mesh points by including the mid-points. Applying the same procedure as in (52) and (53) on 2N number of mesh points, we obtain(54)uxi−Ui2N≤C2N−1+ON−2.Combining (53) and (54) for removing the term CN−1 results in(55)uxi−2Ui2N−Ui≤CN−2,giving (56)Uiext=2Ui2N−Ui,as the extrapolated solution. The error bound for the extrapolated solution in (56) becomes(57)sup0<cε≪1uxi−Uiext≤CN−2.


## 4. Examples and Numerical Results

In this section, we consider numerical examples to illustrate the theoretical findings of the developed schemes.Example 1.
Consider the problem(58)−εu″x+u′x+2ux−δ+5ux−ux+η=0,with interval conditions ux=1,−δ≤x≤0,ux=−1,1≤x≤1+η.Example 2.
Consider the problem(59)−εu″x+1+expx2u′x+xux−δ−x2ux+1−exp−xux+η=−1,with interval conditions ux=1,−δ≤x≤0,ux=−1,1≤x≤1+η.Example 3.
Consider the problem(60)−εu″x−0.5u′x+3ux−δ+2ux−2ux+η=−1,with interval conditions ux=1,−δ≤x≤0,ux=0,1≤x≤1+η.Example 4.
Consider the problem(61)εu″x+2.5u′x−2expxux−δ−ux−xux+η=1,with interval conditions ux=1,−δ≤x≤0,ux=1,1≤x≤1+η.
The exact solution of the constant coefficient interval value problem(62)εu″x+au′x+αux−δ+ωux+βux+η=fx,on Ω=0,1, with the interval conditions ux=ϕx,x∈−δ,0,ux=ψx,x∈1,1+η, is given by(63)ux=c1em1x+c2em2x+fc,where c=α+ω+β, c1=−f+ψ1c+em2f−ϕ0c/cem1−em2, c2=−f+ψ1c+em1f−ϕ0c/cem1−em2, m1=−a−δα+ηβ+a−δα+ηβ2−4cε/2ε, and m2=−a−δα+ηβ−a−δα+ηβ2−4cε/2ε.
The exact solutions of the variable coefficient problems are not known. So, we use the procedure of the double mesh technique to calculate maximum absolute error. The maximum absolute error is defined as(64)EεN=max0≤i≤NUiN−Ui2N,where UiN denotes the solution of the problem on N number of mesh points and Ui2N denotes the numerical solution on 2N number of mesh points by including the mid-points xi+1/2 into the mesh numbers. The uniform error estimate is defined as(65)EN=maxεEεN.
The rate of convergence of the scheme is given by(66)RεN=log2EεNEε2N,and the uniform rate of convergence is given as(67)RN=log2ENE2N.
In Tables1–8, the maximum absolute error of Examples 1–4 using the proposed scheme is given. In Tables 1, 3, 5, and 7, the maximum absolute error before the Richardson extrapolation is given, and in Tables 2, 4, 6, and 8, the maximum absolute error after the Richardson extrapolation is given. As one observes in the tables, for each number of mesh interval N as ε⟶0, the maximum absolute error becomes stable and uniform. This indicates that the proposed scheme convergence is independent of the perturbation parameter ε. In the last two rows of each table, we observe the ε-uniform error and the ε-uniform rate of convergence of the scheme. The scheme before the extrapolation gives first-order uniform convergence and the extrapolated scheme gives second-order uniform convergence. In Tables 9 and 10, we compare the maximum absolute error of the proposed scheme with recently published papers in [14, 15, 19]. As one observes, the proposed scheme gives more accurate result.
In Figure1, the influence of the delay parameter on the behaviour of the solution of Examples 3 and 4 is shown for ε=2−3 and δ=0.1ε,0.5ε, and 0.9ε. From Figure 2, we observe the numerical solution of Examples 1–4 for different values of the perturbation parameter ε=2−5,2−6, and 2−7. As observed in the figures, for ε going small, strong boundary layer is created.Table 1
Maximum absolute error of Example1 before Richardson extrapolation for δ=0.6ε​andη=0.5ε.
ε↓N⟶8163264128256201.3808e − 033.4045e − 048.4818e − 052.1193e − 055.2970e − 061.3242e − 062−41.7200e − 025.0151e − 031.4739e − 033.5993e − 049.0018e − 052.2464e − 052−84.1828e − 022.6301e − 021.4357e − 026.0620e − 031.8873e − 035.0407e − 042−124.1711e − 022.6194e − 021.4838e − 027.9779e − 034.1434e − 032.1118e − 032−164.1704e − 022.6186e − 021.4834e − 027.9750e − 034.1418e − 032.1124e − 032−204.1703e − 022.6186e − 021.4833e − 027.9748e − 034.1417e − 032.1123e − 032−244.1703e − 022.6186e − 021.4833e − 027.9748e − 034.1417e − 032.1123e − 032−284.1703e − 022.6186e − 021.4833e − 027.9748e − 034.1417e − 032.1123e − 032−324.1703e − 022.6186e − 021.4833e − 027.9748e − 034.1417e − 032.1123e − 032−364.1703e − 022.6186e − 021.4833e − 027.9748e − 034.1417e − 032.1123e − 032−404.1703e − 022.6186e − 021.4833e − 027.9748e − 034.1417e − 032.1123e − 03EN4.1703e − 022.6186e − 021.4833e − 027.9748e − 034.1417e − 032.1123e − 03RN0.67140.82000.89530.94520.9714–Table 2
Maximum absolute error of Example1 after Richardson extrapolation for δ=0.6ε​η=0.5ε.
ε↓N⟶8163264128256203.4045e − 042.1193e − 051.3242e − 068.2744e − 085.3296e − 091.3016e − 092−45.0151e − 033.5993e − 042.2464e − 051.4036e − 068.7720e − 082.2310e − 082−82.6301e − 026.0620e − 035.0407e − 043.3479e − 052.0797e − 064.2147e − 072−122.6194e − 027.9779e − 032.1118e − 034.0826e − 043.2897e − 055.2346e − 062−162.6186e − 027.9750e − 032.1124e − 035.3615e − 041.3446e − 043.2541e − 052−202.6186e − 027.9748e − 032.1123e − 035.3614e − 041.3446e − 043.2541e − 052−242.6186e − 027.9748e − 032.1123e − 035.3614e − 041.3446e − 043.2541e − 052−282.6186e − 027.9748e − 032.1123e − 035.3614e − 041.3446e − 043.2541e − 052−322.6186e − 027.9748e − 032.1123e − 035.3614e − 041.3446e − 043.2541e − 052−362.6186e − 027.9748e − 032.1123e − 035.3614e − 041.3446e − 043.2541e − 052−402.6186e − 027.9748e − 032.1123e − 035.3614e − 041.3446e − 043.2541e − 05EN2.6186e − 027.9748e − 032.1123e − 035.3614e − 041.3446e − 043.2541e − 05RN1.71531.91661.97811.99542.0468–Table 3
Maximum absolute error of Example2 before Richardson extrapolation for δ=0.6εandη=0.5ε.
ε↓N⟶8163264128256203.4492e − 031.6925e − 038.3747e − 044.1601e − 042.0724e − 041.0343e − 042−42.2763e − 027.6376e − 031.4739e − 031.5941e − 037.8163e − 043.8784e − 042−82.3670e − 021.1054e − 025.2737e − 032.5716e − 031.2314e − 034.9403e − 042−122.3670e − 021.1054e − 025.2735e − 032.5654e − 031.2638e − 036.2704e − 042−162.3670e − 021.1054e − 025.2735e − 032.5654e − 031.2638e − 036.2704e − 042−202.3670e − 021.1054e − 025.2735e − 032.5654e − 031.2638e − 036.2704e − 042−242.3670e − 021.1054e − 025.2735e − 032.5654e − 031.2638e − 036.2704e − 042−282.3670e − 021.1054e − 025.2735e − 032.5654e − 031.2638e − 036.2704e − 042−322.3670e − 021.1054e − 025.2735e − 032.5654e − 031.2638e − 036.2704e − 042−362.3670e − 021.1054e − 025.2735e − 032.5654e − 031.2638e − 036.2704e − 042−402.3670e − 021.1054e − 025.2735e − 032.5654e − 031.2638e − 036.2704e − 04EN2.3670e − 021.1054e − 025.2735e − 032.5654e − 031.2638e − 036.2704e − 04RN1.09851.06771.03961.02141.0111–Table 4
Maximum absolute error of Example2 after Richardson extrapolation for δ=0.6εandη=0.5ε.
ε↓N⟶8163264128256201.6925e − 034.1601e − 041.0343e − 042.5819e − 056.4524e − 061.1321e − 062−47.6376e − 031.5941e − 033.8784e − 049.6460e − 052.4085e − 057.6001e − 062−81.1054e − 022.5716e − 034.9403e − 041.0586e − 042.6498e − 058.0231e − 062−121.1054e − 022.5654e − 036.2704e − 041.5625e − 043.0982e − 058.2543e − 062−161.1054e − 022.5654e − 036.2704e − 041.5584e − 043.8901e − 059.2541e − 062−201.1054e − 022.5654e − 036.2704e − 041.5584e − 043.8901e − 059.2541e − 062−241.1054e − 022.5654e − 036.2704e − 041.5584e − 043.8901e − 059.2541e − 062−281.1054e − 022.5654e − 036.2704e − 041.5584e − 043.8901e − 059.2541e − 062−321.1054e − 022.5654e − 036.2704e − 041.5584e − 043.8901e − 059.2541e − 062−361.1054e − 022.5654e − 036.2704e − 041.5584e − 043.8901e − 059.2541e − 062−401.1054e − 022.5654e − 036.2704e − 041.5584e − 043.8901e − 059.2541e − 06EN1.1054e − 022.5654e − 036.2704e − 041.5584e − 043.8901e − 059.2541e − 06RN2.10732.03262.00852.00222.0716–Table 5
Maximum absolute error of Example3 before Richardson extrapolation for δ=0.6εandη=0.5ε.
ε↓N⟶8163264128256201.1318e − 032.7968e − 046.9551e − 051.7366e − 054.3405e − 061.0850e − 062−41.4786e − 023.6872e − 038.5968e − 042.1118e − 045.2578e − 051.3152e − 052−84.7613e − 022.9531e − 021.4879e − 025.5516e − 031.4704e − 033.3684e − 042−124.8375e − 023.0553e − 021.7528e − 029.4597e − 034.9243e − 032.4698e − 032−164.8423e − 023.0586e − 021.7549e − 029.4718e − 034.9323e − 032.5185e − 032−204.8426e − 023.0589e − 021.7550e − 029.4725e − 034.9327e − 032.5187e − 032−244.8426e − 023.0589e − 021.7550e − 029.4725e − 034.9328e − 032.5187e − 032−284.8426e − 023.0589e − 021.7550e − 029.4725e − 034.9328e − 032.5187e − 032−324.8426e − 023.0589e − 021.7550e − 029.4725e − 034.9328e − 032.5187e − 032−364.8426e − 023.0589e − 021.7550e − 029.4725e − 034.9328e − 032.5187e − 032−404.8426e − 023.0589e − 021.7550e − 029.4725e − 034.9328e − 032.5187e − 03EN4.8426e − 023.0589e − 021.7550e − 029.4725e − 034.9328e − 032.5187e − 03RN0.66280.80150.88970.94130.9697–Table 6
Maximum absolute error of Example3 after Richardson extrapolation for δ=0.6εandη=0.5ε.
ε↓N⟶8163264128256202.7968e − 041.7366e − 051.0850e − 066.7812e − 084.2386e − 091.0101e − 092−43.6872e − 032.1118e − 041.3152e − 058.2173e − 075.0203e − 081.2341e − 082−82.9531e − 025.5516e − 033.3684e − 042.0489e − 051.2809e − 063.7125e − 072−123.0553e − 029.4597e − 032.4698e − 033.8494e − 042.2756e − 054.8133e − 062−163.0586e − 029.4718e − 032.5185e − 036.3984e − 041.5770e − 043.8042e − 052−203.0589e − 029.4725e − 032.5187e − 036.3989e − 041.6062e − 043.9001e − 052−243.0589e − 029.4725e − 032.5187e − 036.3989e − 041.6062e − 043.9001e − 052−283.0589e − 029.4725e − 032.5187e − 036.3989e − 041.6062e − 043.9001e − 052−323.0589e − 029.4725e − 032.5187e − 036.3989e − 041.6062e − 043.9001e − 052−363.0589e − 029.4725e − 032.5187e − 036.3989e − 041.6062e − 043.9001e − 052−403.0589e − 029.4725e − 032.5187e − 036.3989e − 041.6062e − 043.9001e − 05EN3.0589e − 029.4725e − 032.5187e − 036.3989e − 041.6062e − 043.9001e − 05RN1.69121.91111.97681.99422.0421–Table 7
Maximum absolute error of Example4 before Richardson extrapolation for δ=0.6εandη=0.5ε.
ε↓N⟶8163264128256209.0320e − 034.5788e − 032.2918e − 031.1469e − 035.7338e − 042.8666e − 042−43.7979e − 021.6150e − 026.1639e − 032.4525e − 031.0562e − 034.8456e − 042−84.3631e − 022.4943e − 021.3421e − 026.9160e − 033.1649e − 031.2055e − 032−124.3723e − 022.4996e − 021.3451e − 027.0011e − 033.5749e − 031.8068e − 032−164.3728e − 022.4999e − 021.3452e − 027.0020e − 033.5755e − 031.8071e − 032−204.3729e − 022.4999e − 021.3452e − 027.0021e − 033.5755e − 031.8071e − 032−244.3729e − 022.4999e − 021.3452e − 027.0021e − 033.5755e − 031.8071e − 032−284.3729e − 022.4999e − 021.3452e − 027.0021e − 033.5755e − 031.8071e − 032−324.3729e − 022.4999e − 021.3452e − 027.0021e − 033.5755e − 031.8071e − 032−364.3729e − 022.4999e − 021.3452e − 027.0021e − 033.5755e − 031.8071e − 032−404.3729e − 022.4999e − 021.3452e − 027.0021e − 033.5755e − 031.8071e − 03EN4.3729e − 022.4999e − 021.3452e − 027.0021e − 033.5755e − 031.8071e − 03RN0.80670.89400.94200.96960.9845–Table 8
Maximum absolute error of Example4 after Richardson extrapolation for δ=0.6εandη=0.5ε.
ε↓N⟶8163264128256204.5788e − 031.1469e − 032.8666e − 047.1658e − 051.7914e − 054.4784e − 062−41.6150e − 022.4525e − 034.8456e − 041.1290e − 042.7710e − 056.5434e − 062−82.4943e − 026.9160e − 031.2055e − 031.6886e − 043.2527e − 058.9632e − 062−122.4996e − 027.0011e − 031.8068e − 034.5088e − 047.6356e − 051.2456e − 052−162.4999e − 027.0020e − 031.8071e − 034.5547e − 041.1410e − 042.8541e − 052−202.4999e − 027.0021e − 031.8071e − 034.5548e − 041.1410e − 042.8541e − 052−242.4999e − 027.0021e − 031.8071e − 034.5548e − 041.1410e − 042.8541e − 052−282.4999e − 027.0021e − 031.8071e − 034.5548e − 041.1410e − 042.8541e − 052−322.4999e − 027.0021e − 031.8071e − 034.5548e − 041.1410e − 042.8541e − 052−362.4999e − 027.0021e − 031.8071e − 034.5548e − 041.1410e − 042.8541e − 052−402.4999e − 027.0021e − 031.8071e − 034.5548e − 041.1410e − 042.8541e − 05EN2.4999e − 027.0021e − 031.8071e − 034.5548e − 041.1410e − 042.8541e − 05RN1.83601.95411.98821.99711.9992–Table 9
Example3: comparison of maximum absolute error of proposed scheme and result in [15] for δ=0.6εandη=0.5ε.
ε↓N⟶8163264128256Proposed scheme10−43.0589e − 029.4725e − 032.5187e − 036.3989e − 041.6062e − 043.9001e − 0510−53.0589e − 029.4725e − 032.5187e − 036.3989e − 041.6062e − 043.9001e − 0510−63.0589e − 029.4725e − 032.5187e − 036.3989e − 041.6062e − 043.9001e − 05Result in [15]10−40.102076120.062816380.035561780.019067400.009898950.0050491410−50.102102360.062833820.035570700.019070160.009897940.0050460410−60.102104990.062835570.035571590.019070440.009897980.00504673Table 10
Example4: comparison of maximum absolute error of proposed scheme and results in [14, 19] for δ=0.6εandη=0.5ε.
ε↓N⟶101001010010100Result in [19]Result in [14]Proposed scheme10−11.5339e − 021.9170e − 040.015967000.000172549.2052e − 035.8912e − 0410−22.8175e − 021.8656e − 030.080165820.001995052.0378e − 021.0991e − 0310−32.8534e − 023.3897e − 030.101683690.009999872.0543e − 022.2978e − 0310−42.8570e − 023.3954e − 030.104068040.012304942.0554e − 022.3059e − 03Figure 1
Effect of delay on the solution behaviour forε=2−3. (a) Example 3. (b) Example 4.
(a)(b)Figure 2
Boundary layer formation asε goes small. (a) Example 1. (b) Example 2. (c) Example 3. (d) Example 4.
(a)(b)(c)(d)


## 5. Conclusion

In this paper, singularly perturbed differential difference equations having mixed small shifts on reaction terms of the equation are considered. The considered problem exhibits boundary layer for small values of the perturbation parameter. The bounds and the behaviour of the continuous solution are discussed. Numerical scheme is developed using the technique of exponentially fitted finite difference method. Stability of the scheme is investigated using comparison principle and solution bound. The proposed scheme converges uniformly with rate of convergence of one before Richardson extrapolation and of two after Richardson extrapolation is applied. Test examples exhibiting boundary layers are considered to validate the theoretical finding. The finding in the computation agrees well with the theoretical findings. The proposed scheme gives more accurate results than existing research findings in the literature. In future works, we extend the proposed scheme for singularly perturbed parabolic differential difference equations and singularly perturbed problems with degenerate coefficients.



---

*Source: 6661592-2020-12-28.xml*