Uniform Hybrid Difference Scheme for Singularly Perturbed Differential-Difference Turning Point Problems Exhibiting Boundary Layers

In this paper, a class of linear second-order singularly perturbed differential-difference turning point problems with mixed shifts exhibiting two exponential boundary layers is considered. For the numerical treatment of these problems, first we employ a second-order Taylor’s series approximation on the terms containing shift parameters and obtain a modified singularly perturbed problem which approximates the original problem. *en a hybrid finite difference scheme on an appropriate piecewise-uniform Shishkin mesh is constructed to discretize the modified problem. Further, we proved that the method is almost second-order ε-uniformly convergent in the maximum norm. Numerical experiments are considered to illustrate the theoretical results. In addition, the effect of the shift parameters on the layer behavior of the solution is also examined.


Introduction
Many real-life phenomena in different fields of science are modeled mathematically by delay differential or differentialdifference equations (DDEs). Equations of this type arise widely in scientific fields such as physics, biosciences, ecology, control theory, economics, material science, medicine, and robotics, in which time evolution depends not only on present states but also on the states at or near a given time in the past. DDEs are also prominent in describing several aspects of infectious disease dynamics such as primary infection, drug therapy, and immune response. In addition, statistical analysis of ecological data has shown that there is evidence of delay effects in the population dynamics of many species; for the detail theory of DDEs, one can refer the books [1,2].
If we restrict the class of DDEs in which the highest derivative is multiplied by a small parameter, then it is said to be a singularly perturbed differential-difference equations (SPDDEs) [3]. In the past, less attention had been given for the solutions of SPDDEs. However, in recent years, there has been a growing interest in the treatment of such problems.
is is due to their importance in the modeling of processes in various fields such as optical bistable devices [4], variational problems in control theory [5,6], the hydrodynamics of liquid helium [7], the first exit-time problem [8], describing the human pupil-light reflex [9], microscale heat transfer [10], and a variety of models for physiological processes or diseases [11,12]. e study of different classes of SPDDEs was initiated by Lange and Miura [8,13,14], where they used extension of the method of matched asymptotic expansions for approximating the solution. But in all the cases, they excluded the occurrence of turning points and left it for future study. On the other hand, Kadalbajoo and Sharma [3,15,16] initiated the numerical study of SPDDEs with mixed shifts by constructing a variety of numerical schemes. In recent years, different scholars further developed numerical schemes for SPDDEs with mixed shifts, to mention few [17][18][19][20][21]. Most of the works developed so far focuses only on SPDDEs without turning points. In contrast, there are few works on singularly perturbed delay turning point problems. e papers by Ria and Sharma [22][23][24][25] are the first and also the only notable works in the treatment of such problems when the solutions exhibit both interior and boundary layers, where the authors used a fitted mesh and fitted operator methods and obtained an almost first-order uniform convergence. erefore, it is natural to develop a robust numerical method for turning point problems having a better accuracy and efficiency.
In this paper, we consider the following second-order linear singularly perturbed differential-difference problem containing mixed shifts and with a turning point at x � 0: where a(x), b(x), c(x), d(x), f(x), ϕ(x), and ψ(x) are sufficiently smooth functions on Ω � (− 1, 1), 0 < ε ≪ 1 is the singular perturbation parameter, and 0 < δ ≪ 1 and 0 < η ≪ 1 are the delay and advance parameters, respectively, together with the following assumptions: under these assumptions, problems (1) and (2) possesses a unique solution having boundary layers of exponential type at x � ± 1, i.e., at both end points [25].
Here for the numerical treatment of (1) and (2), we propose an hybrid finite difference scheme on an appropriate piecewise-uniform Shishkin mesh and analyze the stability and uniform convergence the proposed method. Further, we investigate the effect of the shift parameters on the behavior of the solution. roughout this paper, M (sometimes subscripted) denotes a generic positive constant independent of the singular perturbation parameter ε and in the case of discrete problems, also independent of the mesh parameter N. e maximum norm (i.e., ‖f‖ � max − 1≤x≤1 |f(x)|) is used for studying the convergence of the approximate solution to the exact solution of the problem.

The Continuous Problem
Using Taylor's series expansion to approximate the terms containing shift arguments gives us Substituting (7) into (1) and (2) and simplifying gives the following asymptotically equivalent two-point boundary value problem Since (8) is an approximation version of (1) and (2), it is better to use different notation (say u(x)) for the solution of this approximate equation. us, (8) can be written as where C ε (x) � (ε + (δ 2 /2)c(x) . Moreover, the terms a(x), c(x), d(x), δ, and η are such that |A(x) |≥ 2α > 0, for τ < |x| ≤ 1, for some τ > 0, and later on, we will use the term C ε to denote the constant part of C ε (x) (since c(x), d(x) are bounded and δ, η are small parameters, we have C ε � O(ε)). e solution of problem (9) is an approximation to the solution of the original problem (1) and (2).
We establish some a priori results about the solutions and their derivatives for the modified problem (9). Hereinafter, we divide the interval Ω in to three subintervals as First, we consider the following property of the operator L of (9).
and for A(x) < 0 on Ω 3 , we have Proof. For the proof of this theorem, the reader can refer [25].

□
If λ � B(0)/a′(0) < 0, then the solution u(x) is smooth near the turning point x � 0 [26]. Using this, the following theorem gives the bound for the derivatives of the solution in the interval Ω 2 which contains the turning point x � 0.
□ Finally, to prove the uniform convergence of the proposed numerical method, we need to consider the following theorem which provides bounds for the smooth and singular components of the exact solution u of problem (9). Theorem 3. Let A, B and f ∈ C 4 (Ω), and assume that the solution u(x) of problem (9) is decomposed in to the smooth and singular components as en for all i, 0 ≤ i ≤ 4 the smooth component satisfies and the singular component satisfies Proof. For the proof of this theorem, the reader can refer [25,27].

Discrete Problem
In this section, we describe the piecewise-uniform Shishkin mesh for the discretization of the domain and study the behavior of the hybrid difference scheme used for the modified problem (9).
Since the TPP (9) has two boundary layers at x � ±1, we construct a piecewise-uniform Shishkin mesh by subdividing the domain Ω into three subintervals en the discrete mesh Ω N is obtained by putting a uniform mesh with N/4 mesh elements in both Ω L and Ω R , and a uniform mesh with N/2 mesh elements in . ., N, denotes the variable step size. Since the mesh is piecewise-uniform, then the mesh elements are given by Abstract and Applied Analysis where h � 4τ/N and H � 4(1 − τ)/N are uniform mesh lengths for the piecewise-uniform meshes.
If C ε > MN − 1 , then the mesh becomes equally spaced and the transition parameter becomes τ � 1/2 such that On the other hand, for C ε ≤ MN − 1 , the mesh is piecewise-

Hybrid Difference Scheme.
Before describing the scheme, for a given mesh function y(x i ) � y i , we define the forward, backward, and central difference operators D + , D − , and D 0 by respectively, and the second-order central difference operator δ 2 by where Further, we define the midpoint upwind schemes L N M ± and the classical central difference scheme L N C used to approximate the continuous operator L as where A i±1/2 � (A i + A i±1 )/2 and similarly for C ε,i±1/2 , B i±1/2 , and f i±1/2 . Now, we propose the hybrid difference scheme to discretize (9), which consists of the classical central difference scheme when C ε > MN − 1 , and a proper combination of the midpoint upwind schemes in the outer region Ω C and the central difference scheme in the layer regions Ω L and Ω R , whenever C ε ≤ MN − 1 . Hence, the proposed hybrid scheme on Ω N takes the following form: and the right hand side vector f i as After rearranging the terms in (27), we obtain the following system of equations: where the coefficients are given by In general, central difference schemes can be unstable on coarser meshes, but we use this scheme only on the fine part of the Shishikn mesh and thus attain stability under the mild assumption on the minimum number of mesh points N, considered in the following lemma.
holds true. en the discrete operator defined by (27) (27) has a unique solution.
Proof. To prove the results, it is enough to show that the operator given by (30) is an M-matrix. For this, we need to show that (30) satisfies for all the operators defined in (31)-(33). Here, we separately consider the following two cases based on the relation between C ε and N.
Case 1: when C ε > MN − 1 , the mesh is uniform, and we used central difference scheme on the entire domain. us, for M ≥ ‖A‖ and using (22) into (31), we get Abstract and Applied Analysis and some simple calculations gives Case 2: when C ε ≤ MN − 1 , different operators are used in the layer regions and the outer regions.
In the layer regions, it is apparent that r + i < 0 and r − i < 0, for 1 ≤ i < N/4 and 3N/4 < i ≤ N − 1, respectively. Further, using the first assumption of (34) and (23) in (31) we get In both the layer regions, we simply obtain . In addition, using H ≤ 4N − 1 and the second assumption of (34) in (32) and (33) gives us for N/4 ≤ i ≤ N/2 and N/2 + 1 ≤ i ≤ 3N/4, respectively. Moreover, for all N/4 ≤ i ≤ 3N/4, it is easy to verify that r − i + r c i + r + i > 0. For all the cases, it is verified that the operator (30) satisfies the conditions in (35). Hence, the matrix is an M-matrix. erefore, the solution of (27) exists and the maximum principle easily follows. For more details the reader can refer [28,29]. □ Whenever, the conditions of the maximum principle are satisfied, we can take {D i } as a barrier function for {U i }.

Convergence Analysis of the Proposed Method
In this section, we establish the ε-uniform error estimate of the hybrid scheme (27). For this, we consider the two cases C ε > MN − 1 and C ε ≤ MN − 1 separately. For both the cases, analogous to the continuous solution u, we decompose the discrete solution U into a smooth component V and a singular component W, such that U≔V + W, where V is the solution of the nonhomogeneous problem given by and W the solution of the homogeneous problem en the error at each mesh point is which implies and so the error in the smooth and singular components of the solution can be estimated separately. First, to bound the errors we need to consider the truncation error of associated with the discrete operators in (27). For any smooth function y(x), the truncation errors L N M ± applied to y at x i±1/2 and L N C applied to y at y i , becomes T 1± : � L N M ± (y i ) − (Ly)(x i±1/2 ) and T 2 : � L N C (y i )− (Ly)(x i ) respectively, where y i : � y(x i ). us, the bounds are given in the following Lemma.

Lemma 4. Let y(x) be a smooth function defined on [− 1, 1]. en there exists a positive constant M such that
x i+1 x i y ‴ (t) dt, Proof. By repeated use of the fundamental theorem of calculus, one can obtain the proof as in Lemma 3.3 of [28].
To bound the truncation error of the scheme the comparison principle of Lemma 3 alone is not enough, so we will consider the following lemma which enables us to bound the error. □ Lemma 5. Let Z i � 2 + x i for 0 ≤ i ≤ N be the mesh function for (27). en there exists a positive constant M such that Proof. e proof is an easy computation. Sometimes, the truncation error contains a term of magnitude greater than the desired order of convergence, when this happens we shall combine Lemma 3 with the following results.
Whenever C ε ≤ MN − 1 , we define the auxiliary discrete function on the mesh Ω N by using the mesh elements given in (21) as 6 Abstract and Applied Analysis where c is a positive constant.

Lemma 6. For any c > 0 the discrete function {S i } from (45) is such that
Proof. e lower bound for S i follows from the inequality e − t ≤ (1 + t) − 1 which holds true for t ≥ 0. e upper bound for S i is obvious for i > 3N/4. For i ≤ 3N/4, it follows from the inequality (1 + t) − 1 ≤ e − t+t 2 , which holds true for t ≥ 0. Setting t ≔ ch/ε and using (23), we get Proof. e proof is similarly to Lemma 3.3 of [30].

Abstract and Applied Analysis
Proof. For the proof, one can follow similarly like Lemma 3.4 of [30].
□ Remark 1. Because of the symmetry of the mesh and the adaptive nature of the hybrid scheme, it is easy to derive a similar result like Lemmas 7 and 8 using the mesh function {S N− i } related to the layer function e − α(1+x)/C ε . Now we have assembled the tools for the proof of the ε-uniform convergence.

Theorem 4.
Assume that the conditions of (34) holds true. en the hybrid scheme (27) satisfies the following error estimates. (51) where U i the solution of the discrete problem (27) and u(x i ) is the solution of the continuous problem (9) at the mesh points in Ω N .
Proof. Here we estimate the error bounds by considering the two cases separately as follows: Case 1: when C ε > MN − 1 , we employed the central difference scheme on the entire domain. First let us compute the nodal error for the smooth part V i . Using Lemma 4 and the bound of v in eorem 3, the truncation error of the scheme becomes (53) Here, using h � 2N − 1 on the above inequality, we obtain the following estimate: Now, let us take D i ≔ MN − 2 (2 + x i ) as a barrier function for |V i − v(x i )|, then from (39), it is easy to see that (55) Next, we analyze the error bounds for the singular component W i . e local truncation error is bounded in standard way as done above. More precisely, Application of eorem 3 and using (22) on the above inequality and simplifying gives

(57)
Now, arguing similarly like the smooth part, we obtain (58) Using (55) and (58) in to (42) gives the required result of the first case (51). Case 2: for C ε ≤ MN − 1 , the mesh becomes piecewiseuniform and we employed combinations of midpoint upwind and central difference schemes. First, let us compute the nodal error for the smooth part V i .
Similarly like the smooth part of Case 1, the truncation error becomes then using the bound for the truncation error of Lemma 4 and the bound for v(x) from eorem 3, we get Since, C ε ≤ MN − 1 and h i ≤ 4N − 1 , then using these in the above inequality gives us Now, arguing similarly like the smooth part of the previous case, we obtain Next, we analyze the error bounds for the singular component W i . A different argument is used to bound |W − w| in the outer region and layer regions.
In the subinterval with no boundary layer Ω C , both W and w are small, and by the triangle inequality, we have so it suffices to bound W(x i ) and w(x i ) separately. eorem 3 for i � N/4, . . ., 3N/4 gives en using (23) in the above inequality, we get To obtain a similar bound for W(x i ), we set where S i and � S i are from Lemmas 7 and 8, respectively. Now for sufficiently large M 1 , using eorem 3 in (40) and Lemma 8, we get Further, for i � 1, . . ., N − 1, the property of the discrete operator from Lemmas 7 and 8 implies Abstract and Applied Analysis using this in the above inequality, we get From (66)-(69), we can easily observe that D i can be a barrier function for W i for M 1 sufficiently large. erefore, by the discrete maximum principle of Lemma 3 we get In particular, in the coarser region, (70) and Lemmas 7 and 8 together imply that erefore, combining (63), (65), and (71) we get It remains to prove the bound for the singular component in the layer regions Ω L and Ω R . First we estimate the bounds for the singular component in Ω R . Since we employ central difference scheme in Ω R , so as we did for the smooth component, we use the truncation error to bound the error.
us, for i � 3N/4 + 1, . . ., N − 1 the above inequality is reduced to Further, taking i � 3N/4 in (72), we get |W 3N/4 − w (3N/4) | ≤ MN − 2 , and for i � N, the boundary condition in (40) gives |W N − w(x N )| � 0. Now let be the mesh function, where S i is from Lemma 7. If M 2 is chosen large enough, our estimates shows that D i is a barrier function for |W i − w(x i )|. So by using the discrete maximum principle of Lemma 3 and Lemma 7, together with (23), we get Similarly the proof follows for the left boundary layer region, Ω L , i.e., Finally, properly using (62), (72), (76), and (77) in to (42) gives the required bound of the second case (52), which completes the proof.

Test Problems and Numerical Results
To demonstrate the applicability of the proposed method, we have implemented it on two boundary value problems of the form (1) and (2). Since the exact solution for the problems are not available, the pointwise error (e N i ) and maximum absolute error (Ê N ) are calculated by using the double mesh principle given by where U N and U 2N denotes the numerical solutions obtained using N and 2N meshes points, respectively. Further, we determine the corresponding rate of convergence bŷ Example 1. Consider the following homogeneous SPDDE with a turning point at x � 0: In Tables 1 and 2

Abstract and Applied Analysis
Example 2. Consider the following nonhomogeneous SPDDE with a turning point at x � 1/2: In Tables 3 and 4

Discussion
Singularly perturbed differential-difference turning point problems exhibiting twin boundary layers which contain small mixed shifts on the reaction coefficients are considered. For the numerical treatment of these problems, first we employ a second-order Taylor's series approximation on the terms containing shift parameters and obtain a modified singularly perturbed problem which approximates the original problem. And then an efficient hybrid difference scheme on an appropriate piecewise-uniform Shishkin mesh is developed for the modified problem. e proposed method is analyzed for stability and convergence, and it has been shown that the method is ε-uniformly convergent with   an almost second-order rate of convergence. Further, two numerical experiments are examined to support the theoretical results and to illustrate the effect of the small shifts on the layer behavior of the solutions. Tables 1-4 present the computed maximum pointwise error and the rate of convergence for the considered examples. e results demonstrate that the method is robust, i.e., converges for all ε. We also observed that the maximum pointwise error and the rate of convergence stabilizes as ε ⟶ 0 for each appropriate N. Further, the numerical results clearly support the theoretical error bounds and order of convergence.
In addition, to demonstrate the effect of the small shifts on the behavior of the solution graphs of the considered problems are plotted in Figures 1-4 for different values of δ and η. It is observed that the boundary layers are maintained but layer gets shifted as delay/advance argument changes. Shifts in the layer depend upon the type of shift as well as on the value of the coefficients of the term containing delay/ advance.

Data Availability
No data were used to support this study.

Conflicts of Interest
e authors declare that there are no conflicts of interest regarding to the publication of this paper. Abstract and Applied Analysis 13