A finite difference method on layer-adapted mesh for singularly perturbed delay differential equations

Consider an initial value problem for the semilinear second order singularly perturbed delay differential equation [12] in the interval Ī = [0,T]: (1) $ε u^{″} (t) + a (t) u^{'} (t) + f (t, u (t), u (t - r)) = 0, t \in I,$ \varepsilon {u^{\prime\prime}}(t) + a(t){u^\prime}(t) + f(t,u(t),u(t - r)) = 0,\quad t \in I,(2) $u (t) = φ (t), t \in I_{0},$ u(t) = \varphi (t),\quad t \in {I_0},(3) $u^{'} (0) = A / ε,$ {u^\prime}(0) = A/\varepsilon , where $I = (0, T] = \cup_{p = 1}^{m} I_{p}$ I = (0,T] = \mathop \cup \limits_{p = 1}^m {I_p} , I_p = {t : r_p−1 <t ≤ r_p}, 1 ≤ p ≤ m and r_s = sr, for 0 ≤ s ≤ m and I₀ = (−r,0]. 0 < ɛ ≤ 1 is the perturbation parameter, a(t) ≥ α > 0, b(t), c(t), f (t) and ϕ(t) are given sufficiently smooth functions satisfying certain regularity conditions to be specified and r is a constant delay, which is independent of ɛ. Moreover (4) $| \frac{\partial f}{\partial u} | \leq b^{*} and | \frac{\partial f}{\partial v} | \leq c^{*} .$ \left| {{{\partial f} \over {\partial u\;}}} \right| \le {b^ *}\,{\rm{and}}\;\left| {{{\partial f} \over {\partial v\;}}} \right| \le {c^ * }. Delay differential equations play an important role in the mathematical modelling of various practical phenomena in the biosciences and control theory. Any system involving a feedback control will almost always involve time delays. These arise because a finite time is required to sense information and then react to it. A singularly perturbed delay differential equation is an ordinary differential equation in which the highest derivative is multiplied by a small parameter and involving at least one delay term [1,2,3,4]. Such problems arise frequently in the mathematical modelling of various practical phenomena, for example, in the modelling of several physical and biological phenomena like the optically bistable devices [5], description of the human pupil-light reflex [6], a variety of models for physiological processes or diseases and variational problems in control theory where they provide the best and in many cases the only realistic simulation of the observed phenomena [7]. An overview of numerical treatment for first and second order singularly perturbed delay differential equations, may be obtained in [8,9,10,11,12,13,14,15,16] (see, also references therein).

The numerical analysis of singular perturbation cases has always been far from trivial because of the boundary layer behavior of the solution. Such problems undergo rapid changes within very thin layers near the boundary or inside the problem domain. It is well known that standard numerical methods for solving singular perturbation problems do not give satisfactory result when the perturbation parameter is sufficiently small. Therefore, it is important to construct suitable numerical methods for these problems, whose accuracy does not depend on the perturbation parameter, i.e. methods that are uniformly convergent with respect to the perturbation parameter [17,18,19,20,21,22].

In a singularly perturbed delay differential equation, one encounters with two difficulties, one is because of occurrence of the delay term and another one is due to presence of perturbation parameter. To overcome the first difficulty, we employed the numerical method of steps [2] for the delay argument which converted the problem to a initial value problem for a singularly perturbed differential equation. Then we constructed a numerical scheme based on finite difference method on non uniform Shishkin mesh for the numerical solution.

In the present paper we discretize (1)–(2) using a numerical method, which is composed of an exponentially fitted difference scheme on piecewise uniform Shishkin mesh on each time-subinterval. In section 2, we state some important properties of the exact solution. In section 3, we describe the finite difference discretization and introduce the piecewise uniform mesh. In section 4, we present convergence analysis for approximate solution. Uniform convergence is proved in the discrete maximum norm. Some numerical results are being presented in section 5. The technique to construct discrete problem and error analysis for approximate solution is similar to those in [8,9,23,24]

Throughout the paper, C will denote a generic positive constant independent of ɛ and of the mesh parameter.

The Continuous Problem

Here we show some properties of the solution of (1)–(3), which are needed in later sections for the analysis of the appropriate numerical solution. For any continuous function g(t), we use ‖g‖_∞ for the continuous maximum norm on the corresponding interval.

Lemma 2.1

Let δ (t) be nonnegative and continuous function such that (5) $δ (t) \leq δ_{0} + \int_{0}^{t} {c_{0} δ (τ) + c_{1} δ (τ - r)} d τ, t > 0,$ \delta (t) \le {\delta _0} + \int\limits_0^t \left\{ {{c_0}\delta (\tau ) + {c_1}\delta (\tau - r)} \right\}d\tau ,{\kern 1pt} {\kern 1pt} t > 0,(6) $δ (t) = δ_{0}, - r \leq t \leq 0,$ \delta (t) = {\delta _0},\;{\kern 1pt} - r \le t \le 0, where δ₀, c₀ and c₁ are nonnegative constants. Then it holds that $δ (t) \leq γ_{0} \sum_{s = 0}^{p} {\frac{γ_{1}^{s} {(t - r_{p - 1})}^{s}}{s!}}, t \in I_{p}, p = 1, 2, ..., m,$ \delta (t) \le {\gamma _0}\sum\limits_{s = 0}^p \left\{ {{{\gamma _1^s{{(t - {r_{p - 1}})}^s}} \over {s!}}} \right\},{\kern 1pt} \;{\kern 1pt} t \in {I_p}{\kern 1pt} \;{\kern 1pt} ,\,p = 1,2,...,m, where γ₀ = δ₀exp(c₀T) and γ₁ = c₁exp(c₀T).

Proof

See [11].

Lemma 2.2

Let a,b,c, f ∈ C¹ (Ī),ϕ ∈ C¹ (Ī₀). Then for the solution u(t) of problem (1)–(3) the following estimates hold (7) ${‖ u (t) ‖}_{\infty, \bar{I}} \leq C,$ {\left\| {u(t)} \right\|_{\infty ,\bar I}} \le C,(8) $| u^{'} (t) | \leq C {1 + \frac{1}{ε} \exp (- \frac{α t}{ε})}, t \in I,$ \left| {{u^\prime}(t)} \right| \le C\left\{ {1 + {1 \over \varepsilon }exp\,\left( { - {{\alpha t} \over \varepsilon }} \right)} \right\},{\kern 1pt} \;{\kern 1pt} t \in I,(9) $| u^{″} (t) | \leq C {1 + \frac{{(t - r_{p - 1})}^{p - 1}}{ε^{2}} exp (- \frac{α (t - r_{p - 1})}{ε})}, t \in I_{p}, p = 1, 2,$ \left| {{u^{\prime\prime}}(t)} \right| \le C\left\{ {1 + {{{{(t - {r_{p - 1}})}^{p - 1}}} \over {{\varepsilon ^2}}}\exp \left( {\; - {{\alpha (t - {r_{p - 1}})} \over \varepsilon }} \right)} \right\},t \in {I_p},p = 1,2,(10) $| u^{″} (t) | \leq C, t \in I_{p}, 3 \leq p \leq m,$ \left| {{u^{\prime\prime}}(t)} \right| \le C{\kern 1pt} ,{\kern 1pt} \;{\kern 1pt} t \in {I_p},{\kern 1pt} \;{\kern 1pt} 3 \le p \le m,

Proof

The semilinear equation (1) can be written in the form (11) $ε u^{″} (t) + a (t) u^{'} + b (t) u (t) + c (t) u (t - r) = - f (t, 0, 0), t \in I,$ \varepsilon {u^{\prime\prime}}(t) + a(t){u^\prime} + b(t)u(t) + c(t)u(t - r) = - f(t,0,0),\quad t \in I, where $\begin{matrix} b (t) = \frac{\partial f}{\partial u} (t, \tilde{u}, \tilde{v}), c (t) = \frac{\partial f}{\partial v} (t, \tilde{u}, \tilde{v}), \\ \tilde{u} = γ u, \tilde{v} = γ u (t - r), (0 < γ < 1) - intermediate value . \end{matrix}$ \matrix{ {b(t) = {{\partial f} \over {\partial u}}(t,\tilde u,\tilde v),{\kern 1pt} \;{\kern 1pt} c(t) = {{\partial f} \over {\partial v}}(t,\tilde u,\tilde v),} \cr {\tilde u = \gamma u,{\kern 1pt} \;{\kern 1pt} \tilde v = \gamma u(t - r),{\kern 1pt} \;{\kern 1pt} (0 < \gamma < 1) - {\kern 1pt} {\rm intermediate}\;{\rm value}\ .}}

We rewrite (11) in the form (12) $ε u^{″} (t) + a (t) u^{'} (t) = F (t),$ \varepsilon {u^{\prime\prime}}(t) + a(t){u^\prime}(t) = F(t), where $F (t) = - f (t, 0, 0) - b (t) u (t) - c (t) u (t - r) .$ F(t) = - f(t,0,0) - b(t)u(t) - c(t)u(t - r). From (12), we have the following relation for u^′(t) (13) $u^{'} (t) = u^{'} (0) exp (- \frac{1}{ε} \int_{0}^{t} a (s) ds) + \frac{1}{ε} \int_{0}^{t} F (τ) exp (- \frac{1}{ε} \int_{τ}^{t} a (s) ds) d τ .$ {u^\prime}(t) = {u^\prime}(0)\exp \left( {\; - {1 \over \varepsilon }\int\limits_0^t a(s)ds} \right) + {1 \over \varepsilon }\int\limits_0^t F(\tau )\exp \left( {\; - {1 \over \varepsilon }\int\limits_\tau ^t a(s)ds} \right)d\tau . Integrating (13) from 0 to t, we have $\begin{array}{r} u (t) = φ (0) + A ε^{- 1} \int_{0}^{t} exp (- \frac{1}{ε} \int_{0}^{s} a (τ) d τ) ds + \\ + \frac{1}{ε} \int_{0}^{t} ds \int_{0}^{s} F (τ) exp (- \frac{1}{ε} \int_{τ}^{t} a (λ) d λ) d τ . \end{array}$ \matrix{ \hfill {u(t) = \varphi (0) + A{\varepsilon ^{ - 1}}\int\limits_0^t \exp \left( {\; - {1 \over \varepsilon }\int\limits_0^s a(\tau )d\tau } \right)ds + } \cr \hfill {{\kern 1pt} + {1 \over \varepsilon }\int\limits_0^t ds\int\limits_0^s F(\tau )\exp \left( {\; - {1 \over \varepsilon }\int\limits_\tau ^t a(\lambda )d\lambda } \right)d\tau .}} By the change of integral bounds, we get $\begin{array}{r} | u (t) | \leq | φ (0) | + | A ε^{- 1} \int_{0}^{t} exp (- \frac{1}{ε} \int_{0}^{s} a (τ) d τ) ds | + \\ + \frac{1}{ε} | \int_{0}^{t} d τ {F (τ)} \int_{τ}^{t} exp (- \frac{1}{ε} \int_{τ}^{s} a (λ) d λ) ds | \end{array}$ \matrix{ \hfill {|u(t)| \le |\varphi (0)| + \left| {A{\varepsilon ^{ - 1}}\int\limits_0^t \exp \left( {\; - {1 \over \varepsilon }\int\limits_0^s a(\tau )d\tau } \right)ds} \right| + } \cr \hfill {{\kern 1pt} + {1 \over \varepsilon }\left| {\int\limits_0^t d\tau \{ F(\tau )\} \int\limits_\tau ^t \exp \left( {\; - {1 \over \varepsilon }\int\limits_\tau ^s a(\lambda )d\lambda } \right)ds} \right|}} Then $| u (t) | \leq | φ (0) | + α^{- 1} (| A | + {‖ f ‖}_{1}) + α^{- 1} \int_{0}^{t} ({‖ b ‖}_{\infty} | u (τ) | + {‖ c ‖}_{\infty} | u (τ - r) |) d τ$ |u(t)| \le |\varphi (0)| + {\alpha ^{ - 1}}(|A| + {\left\| f \right\|_1}) + {\alpha ^{ - 1}}\int\limits_0^t \left( {{{\left\| b \right\|}_\infty }|u(\tau )| + {{\left\| c \right\|}_\infty }|u(\tau - r)|} \right)d\tau

Applying Lemma 2.1, we obtain (7).

Due to ‖F‖_∞ ≤ C, it now follows from (13) that $\begin{matrix} | u^{'} (t) | \leq | u^{'} (0) | exp (- \frac{1}{ε} \int_{0}^{t} a (s) ds) + \frac{1}{ε} \int_{0}^{t} | F (τ) | exp (- \frac{1}{ε} \int_{τ}^{t} a (s) ds) d τ . \\ \leq C {| A | ε^{- 1} \exp (- \frac{α t}{ε}) + α^{- 1} (1 - \exp (- \frac{α t}{ε})} \end{matrix}$ \matrix{ {|{u^\prime}(t)| \le |{u^\prime}(0)|\exp \left( {- {1 \over \varepsilon }\int\limits_0^t a(s)ds} \right) + {1 \over \varepsilon }\int\limits_0^t \left| {F(\tau )} \right|\exp \left( {- {1 \over \varepsilon }\int\limits_\tau ^t a(s)ds} \right)d\tau .} \cr { \le C\left\{ {|A|{\varepsilon ^{ - 1}}exp \left( { - {{\alpha t} \over \varepsilon }} \right) + {\alpha ^{ - 1}}(1 - exp \left( { - {{\alpha t} \over \varepsilon }} \right)} \right\}}} which proves (8)

Now, differentiating the Eq.(1) we have (14) $ε u^{'″} (t) + a (t) u^{″} (t) = Φ (t), t \in I_{k + 1},$ \varepsilon {u^{\prime\prime\prime}}(t) + a(t){u^{\prime\prime}}(t) = \Phi (t),\quad t \in {I_{k + 1}}, where $Φ (t) = f^{'} (t, 0, 0) - (a^{'} (t) + b (t)) (t) u^{'} (t) - b^{'} (t) u (t) - c^{'} (t) u (t - r) - c (t) u^{'} (t - r) .$ \Phi (t) = {f^\prime}(t,0,0) - ({a^\prime}(t) + b(t))(t){u^\prime}(t) - {b^\prime}(t)u(t) - {c^\prime}(t)u(t - r) - c(t){u^\prime}(t - r). From (14) we have the following relation for u^″(t) (15) $u^{″} (t) = u^{″} (r_{k}) exp (- \frac{1}{ε} \int_{r_{k}}^{t} a (s) ds) + \frac{1}{ε} \int_{r_{k}}^{t} Φ (τ) exp (- \frac{1}{ε} \int_{τ}^{t} a (s) ds) d τ .$ {u^{\prime\prime}}(t) = {u^{\prime\prime}}({r_k})\exp \left( {\; - {1 \over \varepsilon }\int\limits_{{r_k}}^t a(s)ds} \right) + {1 \over \varepsilon }\int\limits_{{r_k}}^t \Phi (\tau )\exp \left( {\; - {1 \over \varepsilon }\int\limits_\tau ^t a(s)ds} \right)d\tau . and it is easy to see that (16) $\begin{matrix} | Φ (t) | \leq | f^{'} (t) | + | (a^{'} (t) + b (t)) u^{'} (t) | + | b^{'} (t) u (t) | + | c^{'} (t) u (t - r) | + | c (t) u^{'} (t - r) | \\ \leq C (1 + | (a^{'} (t) + b (t)) u^{'} (t) | + | c (t) u^{'} (t - r) |) . \end{matrix}$ \matrix{ {\left| {\Phi (t)} \right| \le \left| {{f^\prime}(t)} \right| + \left| {({a^\prime}(t) + b(t)){u^\prime}(t)} \right| + \left| {{b^\prime}(t)u(t)} \right| + \left| {{c^\prime}(t)u(t - r)} \right| + \left| {c(t){u^\prime}(t - r)} \right|} \cr { \le C\left( {1 + \left| {({a^\prime}(t) + b(t)){u^\prime}(t)} \right| + \left| {c(t){u^\prime}(t - r)} \right|} \right).}}

It follows from (1.1) that (17) $| u^{″} (0) | \leq \frac{C}{ε^{2}} .$ \left| {{u^{\prime\prime}}(0)} \right| \le {C \over {{\varepsilon ^2}}}.

Using the estimates (8), (16) and (17) in (21) for k = 0, we have $\begin{matrix} | u^{″} (t) | \leq \frac{C}{ε^{2}} exp (\frac{- α t}{ε}) + \frac{1}{ε} C \int_{0}^{t} (1 + \frac{1}{ε} \exp (\frac{- α τ}{ε})) exp (\frac{- α (t - τ)}{ε}) d τ \\ \leq C (1 + \frac{1}{ε^{2}} \exp (\frac{- α t}{ε})), t \in I_{1}, \end{matrix}$ \matrix{ {\left| {{u^{\prime\prime}}(t)} \right| \le {C \over {{\varepsilon ^2}}}\exp \left( {{{ - \alpha t} \over \varepsilon }} \right) + {1 \over \varepsilon }C\int\limits_0^t \left( {1 + {1 \over \varepsilon }exp \left( {{{ - \alpha \tau } \over \varepsilon }} \right)} \right)\exp \left( {{{ - \alpha (t - \tau )} \over \varepsilon }} \right)d\tau } \cr { \le C\left( {1 + {1 \over {{\varepsilon ^2}}}exp \left( {{{ - \alpha t} \over \varepsilon }} \right)} \right),{\kern 1pt} \;{\kern 1pt} t \in {I_1},}} i.e., the inequality (9) is valid for p = 1.

Next, from the last inequality for t = r, we have $| u^{″} (r) | \leq C$ \left| {{u^{\prime\prime}}(r)} \right| \le C and thereby from (15) it follows that $\begin{matrix} | u^{″} (t) | \leq | u^{″} (r) | exp (\frac{- α t}{ε}) \\ + \frac{1}{ε} C \int_{0}^{t} (1 + \frac{1}{ε} \exp (\frac{- α (τ - r)}{ε})) exp (\frac{- α (t - τ)}{ε}) d τ \\ \leq C (1 + \frac{1}{ε^{2}} (t - r) \exp (\frac{- α t}{ε})), t \in I_{2} . \end{matrix}$ \matrix{ {\left| {{u^{\prime\prime}}(t)} \right| \le |{u^{\prime\prime}}(r)|\exp \left( {{{ - \alpha t} \over \varepsilon }} \right)} \cr { + {1 \over \varepsilon }C\int\limits_0^t \left( {1 + {1 \over \varepsilon }exp \left( {\;{{ - \alpha (\tau - r)} \over \varepsilon }} \right)} \right)\exp \left( {{{ - \alpha (t - \tau )} \over \varepsilon }} \right)d\tau } \cr { \le C\left( {1 + {1 \over {{\varepsilon ^2}}}(t - r)exp \left( {{{ - \alpha t} \over \varepsilon }} \right)} \right),{\kern 1pt} \;{\kern 1pt} t \in {I_2}.}} i.e. the inequality (9) is also valid for p = 2.

Since $| Φ^{'} (t) | \leq C, t \in I_{p}, p = 3, ..., m,$ \left| {{\Phi ^\prime}(t)} \right| \le C,{\kern 1pt} \;{\kern 1pt} t \in {I_p}{\kern 1pt} \;{\kern 1pt} ,{\kern 1pt} \;{\kern 1pt} p = 3,...,m, the estimate (10) follows immediately from (15).

Discretization and Mesh

In this section, we construct a numerical scheme for solving the initial value problem (1)–(2).

Let ${\bar{ω}}_{N_{0}}$ {\bar \omega _{{N_0}}} be any non-uniform mesh on Ī: ${\bar{ω}}_{N_{0}} = {0 = t_{0} < t_{1} < ... < t_{N_{0}} = T, h_{i} = t_{i} - t_{i - 1}}$ {\bar \omega _{{N_0}}} = \left\{ {0 = {t_0} < {t_1} < ... < {t_{{N_0}}} = T,{\kern 1pt} \;{\kern 1pt} {h_i} = {t_i} - {t_{i - 1}}} \right\} which contains by N mesh point at each subinterval I_p(1 ≤ p ≤ m) : $ω_{N, p} = {t_{i} : (p - 1) N + 1 \leq i \leq pN}, 1 \leq p \leq m,$ {\omega _{N,p}} = \left\{ {{t_i}:(p - 1)N + 1 \le i \le pN} \right\},{\kern 1pt} \;{\kern 1pt} 1 \le p \le m, and consequently $ω_{N_{0}} = \cup_{p = 1}^{m} ω_{N, p} .$ {\omega _{{N_0}}} = \bigcup\limits_{p = 1}^m {\omega _{N,p}}.

To simplify the notation we set w_i = w(t_i) for any function w(t), and moreover y_i denotes an approximation of u(t) at t_i. For any mesh function {w_i} defined on ${\bar{ω}}_{N_{0}}$ {\bar \omega _{{N_0}}}, we use $\begin{matrix} w_{\bar{t}, i} = (w_{i} - w_{i - 1}) / h_{i}, w_{t, i} = (w_{i + 1} - w_{i}) / h_{i + 1}, \\ w_{\overset{0}{t, i}} = (w_{\bar{t}, i} + w_{t, i}) / 2, w_{\bar{t} t, i} = (w_{t, i} - w_{\bar{t}, i}) / h_{i}, \\ {‖ w ‖}_{\infty, N, p} = {‖ w ‖}_{\infty, ω_{N, p}} : = max_{1 \leq i \leq N} | w_{i} | . \end{matrix}$ \matrix{ {{w_{\bar t,i}} = ({w_i} - {w_{i - 1}})/{h_i},{\kern 1pt} \;{\kern 1pt} {w_{t,i}} = ({w_{i + 1}} - {w_i})/{h_{i + 1}},} \cr {{w_{\mathop t\limits^0 ,i}} = ({w_{\bar t,i}} + {w_{t,i}})/2,{\kern 1pt} \;{\kern 1pt} {w_{\bar tt,i}} = ({w_{t,i}} - {w_{\bar t,i}})/{h_i},} \cr {{{\left\| w \right\|}_{\infty ,N,p}} = {{\left\| w \right\|}_{\infty ,{\omega _{N,p}}}}: = \mathop {\max }\limits_{1 \le i \le N} \left| {{w_i}} \right|.}}

For the difference approximation of (1), we are using the following identity (18) $h_{i}^{- 1} \int_{t_{i - 1}}^{t_{i + 1}} Lu (t) ψ_{i} (t) dt = 0, 1 \leq i \leq N_{0} - 1,$ h_i^{ - 1}\int\limits_{{t_{i - 1}}}^{{t_{i + 1}}} Lu(t){\psi _i}(t)dt = 0,1 \le i \le {N_0} - 1, where exponential basis functions $ψ_{i} (t) = {\begin{array}{l} ψ_{1 i} (t) \equiv \frac{e^{a_{i} (t - t_{i - 1}) / ε} - 1}{e^{a_{i} h / ε} - 1}, t_{i - 1} < t < t_{i} \\ ψ_{2 i} (t) \equiv \frac{1 - e^{- a_{i} (t_{i + 1} - t) / ε}}{1 - e^{- a_{i} h / ε}}, t_{i} < t < t_{i + 1}, \\ 0, t \notin (t_{i - 1}, t_{i + 1}), \end{array}$ {\psi _i}(t) = \left\{ {\matrix{ {{\psi _{1i}}(t) \equiv {{{e^{{a_i}(t - {t_{i - 1}})/\varepsilon }} - 1} \over {{e^{{a_i}h/\varepsilon }} - 1}},{t_{i - 1}} < t < {t_i}} \hfill \cr {{\psi _{2i}}(t) \equiv {{1 - {e^{ - {a_i}({t_{i + 1}} - t)/\varepsilon }}} \over {1 - {e^{ - {a_i}h/\varepsilon }}}},{t_i} < {\kern 1pt} t < {t_{i + 1}},} \hfill \cr {0\quad \quad \quad \quad \quad \quad \quad \quad \quad ,t \notin ({t_{i - 1}},{t_{i + 1}}),} \hfill \cr } } \right. and $h_{i}^{- 1} \int_{t_{i - 1}}^{t_{i + 1}} ψ_{i} (t) dt = 1 .$ h_i^{ - 1}\int_{{t_{i - 1}}}^{{t_{i + 1}}} {\psi _i}\left( t \right)dt = 1.

We note that functions ψ_1i(t) and ψ_2i(t) are the solutions of the following problems, respectively: $\begin{array}{l} \begin{array}{l} \begin{array}{l} ε ψ^{″} (t) - a_{i} ψ^{'} (t) = 0, t_{i - 1} < t < t_{i} \\ ψ (t_{i - 1}) = 0, ψ (t_{i}) = 1, \end{array} \end{array} \\ \begin{array}{l} ε ψ^{″} (t) - a_{i} ψ^{'} (t) = 0, t_{i} < t < t_{i + 1} \\ ψ (t_{i}) = 1, ψ (t_{i + 1}) = 0 . \end{array} \end{array}$ \matrix{ \matrix{ \matrix{ {\varepsilon {\psi ^{\prime\prime}}(t) - {a_i}{\psi ^\prime}(t) = 0{\kern 1pt} {\kern 1pt} ,{\kern 1pt} {\kern 1pt} {t_{i - 1}} < t < {t_i}} \hfill \cr {\psi ({t_{i - 1}}) = 0{\kern 1pt} {\kern 1pt} ,{\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} \psi ({t_i}) = 1,} \hfill \cr } \hfill \cr \hfill \cr} \hfill \cr {\matrix{ {\varepsilon {\psi ^{\prime\prime}}(t) - {a_i}{\psi ^\prime}(t) = 0{\kern 1pt} {\kern 1pt} ,{\kern 1pt} {\kern 1pt} {t_i} < t < {t_{i + 1}}} \hfill \cr {\psi ({t_i}) = 1{\kern 1pt} {\kern 1pt} ,{\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} \psi ({t_{i + 1}}) = 0.} \hfill} } \hfill}

Using interpolating quadrature rules with the weight and remainder terms in integral form on subintervals [t_i−1, t_i+1], consistent with [23,24], after a simple calculation, we have the following relation: (19) $ε θ_{i} u_{\bar{t} t, i} + a_{i} u_{\overset{0}{t}, i} + f (t_{i}, u_{i}, u_{i - N}) + R_{i} = 0, i = 1, 2, ..., N_{0} - 1,$ \varepsilon {\theta _i}{u_{\bar tt,i}} + {a_i}{u_{\mathop t\limits^0 ,i}} + f({t_i},{u_i},{u_{i - N}}) + {R_i} = 0,\quad i = 1,2,...,{N_0} - 1, where the factor (20) $θ_{i} = \frac{ρ a_{i}}{2} coth (a_{i} ρ / 2), ρ = h / ε$ {\theta _i} = {{\rho {a_i}} \over 2}\coth ({a_i}\rho /2){\kern 1pt} {\kern 1pt} ,{\kern 1pt} {\kern 1pt} \rho = h/\varepsilon and remainder term (21) $\begin{matrix} R_{i} = h_{i}^{- 1} \int_{t_{i - 1}}^{t_{i + 1}} [a (t) - a (t_{i})] u^{'} (t) ψ_{i} (t) dt \\ + h_{i}^{- 1} \int_{t_{i - 1}}^{t_{i + 1}} \frac{d}{dt} f (τ, u (τ), u (τ - r)) K_{0, i}^{*} (t, τ) d τ \\ K_{0, i}^{*} (t, τ) = T_{0} (t_{i} - τ) - T_{0} (t_{i} - τ), T_{0} (t) = 1, t > 0; T_{0} (t) = 0, t < 0 . \end{matrix}$ \matrix{ {{R_i} = h_i^{ - 1}\int\limits_{{t_{i - 1}}}^{{t_{i + 1}}} [a(t) - a({t_i})]{u^\prime}(t){\psi _i}(t)dt} \cr { + h_i^{ - 1}\int\limits_{{t_{i - 1}}}^{{t_{i + 1}}} {d \over {dt}}f(\tau ,u\left( \tau \right),u\left( {\tau - r} \right))K_{0,i}^ * \left( {t,\tau } \right)d\tau } \cr {K_{0,i}^ * \left( {t,\tau } \right) = {T_0}\left( {{t_i} - \tau } \right) - {T_0}\left( {{t_i} - \tau } \right),{T_0}(t) = 1,\;t > 0;{T_0}(t) = 0,\;t < 0.}}

To define an approximation for the boundary condition (3), we proceed our discretization process by $\int_{t_{0}}^{t_{1}} Lu (t) φ_{0} (t) dt = 0$ \int\limits_{{t_0}}^{{t_1}} Lu(t){\varphi _0}(t)dt = 0 where the exponential basis function $ψ_{0} (t) = {\begin{array}{l} \frac{1 - e^{- a_{0} (t_{1} - t) / ε}}{1 - e^{- a_{0} h / ε}}, & t_{0} < t < t_{1} \\ 0 & , t \notin (t_{0}, t_{1}) \end{array}$ {\psi _0}(t) = \left\{ {\matrix{ {{{1 - {e^{ - {a_0}({t_1} - t)/\varepsilon }}} \over {1 - {e^{ - {a_0}h/\varepsilon }}}},} \hfill & {{t_0} < {\kern 1pt} t < {t_1}} \hfill \cr {0\quad \quad \quad \quad\,\,\,\,\,,} \hfill & {t \notin ({t_0},{t_1})} \hfill \cr } } \right.

We note that functions ϕ₀(t) is the solution of the following problem: $\begin{array}{l} ε ψ^{″} (t) - a_{0} ψ^{'} (t) = 0, t_{0} < t < t_{1}, \\ ψ (t_{0}) = 0, ψ (t_{1}) = 1 . \end{array}$ \eqalign{& {\varepsilon {\psi ^{\prime\prime}}(t) - {a_0}{\psi ^\prime}(t) = 0{\kern 1pt} {\kern 1pt} ,{\kern 1pt} {\kern 1pt} {t_0} < t < {t_1},} \cr & {\psi ({t_0}) = 0{\kern 1pt} {\kern 1pt} ,{\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} \psi ({t_1}) = 1.}}

Whence, as similar above we can write the following difference relation: (22) $ε σ u_{t, 0} - A + r^{(0)} = 0$ \varepsilon \sigma {u_{t,0}} - A + {r^{(0)}} = 0 where the coefficient and the remainder term (23) $σ = 1 + a_{0} ε^{- 1} \int_{t_{0}}^{t_{1}} ψ_{0} (t) dt = \frac{a_{0} ρ}{1 - e^{- a_{0} ρ}}$ \sigma = 1 + {a_0}{\varepsilon ^{ - 1}}\int\limits_{{t_0}}^{{t_1}} {\psi _0}(t)dt = {{{a_0}\rho } \over {1 - {e^{ - {a_0}\rho }}}}(24) $r^{(0)} = \int_{t_{0}}^{t_{1}} [a (t) - a_{0}] u^{'} (t) ψ_{0} (t) dt + \int_{t_{0}}^{t_{1}} f (t, u (t), u (t - r)) ψ_{0} (t) dt .$ {r^{(0)}} = \int\limits_{{t_0}}^{{t_1}} \left[ {a(t) - {a_0}} \right]{u^\prime}(t){\psi _0}(t)dt + \int\limits_{{t_0}}^{{t_1}} f(t,u(t),u(t - r)){\psi _0}(t)dt.

Neglecting R_i and r⁽⁰⁾ in (19) and (22), we propose the following difference scheme for approximation (1)–(2): (25) $ε θ_{i} y_{\bar{t} t, i} + a_{i} y_{\overset{0}{t}, i} + f (t_{i}, y_{i}, y_{i - N}) = 0, i = 1, 2, ..., N_{0} - 1,$ \varepsilon {\theta _i}{y_{\bar tt,i}} + {a_i}{y_{\mathop t\limits^0 ,i}} + f({t_i},{y_i},{y_{i - N}}) = 0,\quad i = 1,2,...,{N_0} - 1,(26) $y_{i} = φ_{i}, - N \leq i \leq 0, ε σ y_{t, 0} - A = 0,$ {y_i} = {\varphi _i},\ - N \le i \le 0,\varepsilon \sigma {y_{t,0}} - A = 0, where θ_i and σ are defined by (20),(23), respectively.

The difference scheme (25)–(26), in order to be ɛ- uniform convergent, we will use the non uniform Shishkin mesh. For the even number N, the piecewise uniform mesh ω_{N, p} divides each of the interval [r_p−1, σ_p] and [σ_p, r_p] into N/2 equidistant subintervals, where the transition point σ_p, which separates the fine and coarse portions of the mesh is obtained by $σ_{p} = r_{p - 1} + min {r / 2, α^{- 1} β_{p} ε ln N},$ {\sigma _p} = {r_{p - 1}} + \min \left\{ {r/2,{\kern 1pt} \;{\kern 1pt} {\alpha ^{ - 1}}{\beta _p}\varepsilon \ln N} \right\}, where β₁ ≥ 1 and β_p > 1 (2 ≤ p ≤ m) are some constants. Hence, if $h_{p}^{(1)}$ h_p^{(1)} and $h_{p}^{(2)}$ h_p^{(2)} denote the stepsizes in [r_p−1, σ_p] and [σ_p, r_p] respectively, we have $h_{p}^{(1)} = 2 (σ_{p} - r_{p - 1}) N^{- 1}, h_{p}^{(2)} = 2 (r_{p} - σ_{p}) N^{- 1}, 1 \leq p \leq m,$ h_p^{(1)} = 2({\sigma _p} - {r_{p - 1}}){N^{ - 1}},{\kern 1pt} \;{\kern 1pt} h_p^{(2)} = 2({r_p} - {\sigma _p}){N^{ - 1}},{\kern 1pt} \;{\kern 1pt} 1 \le p \le m, so $\begin{matrix} {\bar{ω}}_{N, p} = {\begin{array}{l} t_{i} = r_{p - 1} + (i - (p - 1) N) h_{p}^{(1)}, i = (p - 1) N, ..., (p - 1 / 2) N \\ t_{i} = σ_{p} + (i - (p - 1 / 2) N) h_{p}^{(2)}, i = (p - 1 / 2) N + 1, ..., pN, \end{array} \\ 1 \leq p \leq m . \end{matrix}$ \matrix{ {{{\bar \omega }_{N,p}} = \left\{ {\matrix{ {{t_i} = {r_{p - 1}} + (i - (p - 1)N)h_p^{(1)},{\kern 1pt} \;{\kern 1pt} i = (p - 1)N,...,(p - 1/2)N} \hfill \cr {{t_i} = {\sigma _p} + (i - (p - 1/2)N)h_p^{(2)},{\kern 1pt} \;{\kern 1pt} i = (p - 1/2)N + 1,...,pN,} \hfill \cr } } \right.} \cr {1 \le p \le m.}}

In the rest of the paper we only consider this type mesh.

Convergence Analysis

We now estimate the approximate error z_i = y_i − u_i, 0 ≤ i ≤ N₀, which satisfies the discrete problem (27) $ε θ_{i} z_{\bar{t} t, i} + a_{i} z_{\overset{0}{t}, i} + f (t_{i}, y_{i}, y_{i - N}) - f (t_{i}, u_{i}, u_{i - N}) = R_{i}, i = 1, ..., N_{0} - 1,$ \varepsilon {\theta _i}{z_{\bar tt,i}} + {a_i}{z_{\mathop t\limits^0 ,i}} + f({t_i},{y_i},{y_{i - N}}) - f({t_i},{u_i},{u_{i - N}}) = {R_i},i = 1,...,{N_0} - 1,(28) $z_{i} = 0, - N \leq i \leq 0, ε σ z_{t, 0} - A + r^{(0)} = 0,$ {z_i} = 0,\ - N \le i \le 0,\varepsilon \sigma {z_{t,0}} - A + {r^{(0)}} = 0, where the truncation error R_i and r⁽⁰⁾ are given by (21) and (24).

Lemma 4.1

Let z_i be the solution of the problem (27)–(28). Then the following estimate holds (29) ${‖ z ‖}_{\infty, N, p} \leq γ ((θ_{*} / σ) | r^{(0)} | + C \sum_{k = 1}^{p} {‖ R ‖}_{\infty, ω_{N, k}}, 1 \leq p \leq m .$ {\left\| z \right\|_{\infty ,N,p}} \le \gamma (({\theta _ * }/\sigma )|{r^{(0)}}| + C\sum\limits_{k = 1}^p {\left\| R \right\|_{\infty ,{\omega _{N,k}}}},{\kern 1pt} \;{\kern 1pt} 1 \le p \le m. where $θ_{*} = \frac{ρ α}{2} coth (α ρ / 2)$ {\theta _ * } = {{\rho \alpha } \over 2}\coth (\alpha \rho /2) , γ = 4α⁻¹ exp(4α⁻¹ (b^* + c^*))

Proof

Let z_{t, i} = v_i. Then the relation (27) can be rewritten as $ε θ_{i} v_{\bar{t}, i} + \frac{a_{i}}{2} (v_{i} + v_{i - 1}) = F_{i},$ \varepsilon {\kern 1pt} {\theta _i}{v_{\bar t,i}} + {{{a_i}} \over 2}({v_i} + {v_{i - 1}}) = {F_i}, where $F_{i} = R_{i} - f (t_{i}, y_{i}, y_{i - N}) + f (t_{i}, u_{i}, u_{i - N}) .$ {F_i} = {R_i} - f({t_i},{y_i},{y_{i - N}}) + f({t_i},{u_i},{u_{i - N}}).

Solving the first order difference equation with respect to v_i, we get (30) $v_{i} = v_{0} Q_{i} + h \sum_{k = 1}^{i} \frac{F_{𝓁}}{ε θ_{k} + {ha}_{k} / 2} Q_{i - k},$ {v_i} = {v_0}{Q_i} + h\sum\limits_{k = 1}^i {{{F_\ell }} \over {\varepsilon {\theta _k} + h{a_k}/2}}{Q_{i - k}}, where $Q_{i - k} = {\begin{array}{l} 1, k = i, \\ \prod_{j = k + 1}^{i} (\frac{1 - a_{j} ρ / (2 θ_{j})}{1 + a_{j} ρ / (2 θ_{j})}), 0 \leq k \leq i - 1 . \end{array}$ {Q_{i - k}} = \left\{ {\matrix{ {1{\kern 40pt},{\kern 1pt} {\kern 1pt} k = i,} \hfill \cr {\prod\limits_{j = k + 1}^i \left( {{{1 - {a_j}\rho /(2{\theta _j})} \over {1 + {a_j}\rho /(2{\theta _j})}}} \right),{\kern 1pt} {\kern 1pt} 0 \le k \le i - 1.} \hfill \cr } } \right.

Then, since $z_{p} = h \sum_{i = 0}^{p - 1} v_{i} = h \sum_{i = 1}^{p} v_{i - 1}$ {z_p} = h\sum\limits_{i = 0}^{p - 1} {v_i} = h\sum\limits_{i = 1}^p {v_{i - 1}} and taking into consideration (5), from (30) after some manipulations we have the following inequality $| z_{p} | \leq 4 α^{- 1} (ε θ_{*} | z_{t, 0} | + h \sum_{i = 1}^{p - 1} (| R_{i} | + b^{*} | z_{i} | + c^{*} | z_{i - N} |)), 1 \leq p \leq N_{0} - 1 .$ \left| {{z_p}} \right|{\kern 1pt} \le 4{\alpha ^{ - 1}}(\varepsilon {\theta _ * }\left| {{z_{t,0}}} \right| + h\sum\limits_{i = 1}^{p - 1} (|{R_i}| + {b^ * }|{z_i}| + {c^ * }|{z_{i - N}}|)),{\kern 1pt} 1 \le p \le {N_0} - 1.

Replacing i − N = j, we get $\begin{matrix} h \sum_{i = 1}^{p} | z_{i - N} | = h \sum_{j = 1 - N}^{p - N} | z_{j} | \leq h \sum_{j = 1 - N}^{0} | z_{j} | + h \sum_{j = 1}^{p - 1} | z_{j} | \\ = {‖ ψ ‖}_{L_{1} (w_{N, 0})} + h \sum_{j = 1}^{p - 1} | z_{j} | . \end{matrix}$ \matrix{ {h\sum\limits_{i = 1}^p |{z_{i - N}}| = h\sum\limits_{j = 1 - N}^{p - N} |{z_j}| \le h\sum\limits_{j = 1 - N}^0 |{z_j}| + h\sum\limits_{j = 1}^{p - 1} |{z_j}|} \cr {\kern 40pt} { = {{\left\| \psi \right\|}_{{L_1}({w_{N,0}})}} + h\sum\limits_{j = 1}^{p - 1} |{z_j}|.}}

Therefore $| z_{p} | \leq C (ε θ_{*} | z_{t, 0} | + {‖ ψ ‖}_{1} + \sum_{i = 1}^{p} (| R_{i} | + b^{*} | z_{i} | + c^{*} | z_{i - 1} |)$ \left| {{z_p}} \right| \le C(\varepsilon {\theta _ * }\left| {{z_{t,0}}} \right| + {\left\| \psi \right\|_1} + \sum\limits_{i = 1}^p (|{R_i}| + {b^*}|{z_i}| + {c^*}|{z_{i - 1}}|)

Taking also into account (28), and using difference analogue of Gronwall’s inequality this leads to (29).

Lemma 4.2

Let a,b,c, f ∈ C¹ (Ī), ψ ∈ C¹ (Ī₀). Then for the truncation errors R_i and r⁽⁰⁾, the following estimates hold (31) $\begin{matrix} {‖ R ‖}_{\infty, ω_{N}, p} \leq {CN}^{- 1} ln N, 1 \leq p \leq m, \\ | r^{(0)} | \leq {CN}^{- 1} ln N . \end{matrix}$ \matrix{ {{{\left\| R \right\|}_{\infty ,{\omega _N},p}} \le C{N^{ - 1}}\ln N,{\kern 1pt} \;{\kern 1pt} 1 \le p \le m,} \cr {|{r^{(0)}}| \le C{N^{ - 1}}\ln N.}}

Proof

From explicit expression (21) for R_i on an arbitrary mesh, we have $\begin{array}{l} | R_{i} | = h_{i}^{- 1} \int_{t_{i - 1}}^{t_{i + 1}} | (a (t) - a (t_{i})) | | u^{'} (t) | ψ_{i} (t) dt \\ + h_{i}^{- 1} \int_{t_{i - 1}}^{t_{i + 1}} | \frac{d}{dt} f (t, u (t), u (t - r)) | dt, 1 \leq i \leq N_{0} . \end{array}$ \matrix{ {|{R_i}| = h_i^{ - 1}\int\limits_{{t_{i - 1}}}^{{t_{i + 1}}} |(a(t) - a({t_i}))||{u^\prime}(t)|{\psi _i}(t)dt} \hfill \cr { + h_i^{ - 1}\int\limits_{{t_{i - 1}}}^{{t_{i + 1}}} \left| {{d \over {dt}}f(t,u\left( t \right),u\left( {t - r} \right))} \right|dt,{\kern 1pt} \;{\kern 1pt} 1 \le i \le {N_0}.} \hfill}

This inequality together with (7), enable us to write $| R_{i} | \leq C {h_{i} + \int_{t_{i - 1}}^{t_{i + 1}} (| u^{'} (t) | + | u^{'} (t - r) |) dt}, 1 \leq i \leq N_{0} .$ \left| {{R_i}} \right| \le C\left\{ {{h_i} + \int\limits_{{t_{i - 1}}}^{{t_{i + 1}}} \left( {\left| {{u^\prime}(t)} \right| + \left| {{u^\prime}(t - r)} \right|} \right)dt} \right\},{\kern 1pt} \;{\kern 1pt} 1 \le i \le {N_0}.

From here, in view of (8), it follows that (32) $| R_{i} | \leq C {h_{i} + \frac{1}{ε} \int_{t_{i - 1}}^{t_{i + 1}} e^{- \frac{α t}{ε}} dt}, 1 \leq i \leq N_{0}$ \left| {{R_i}} \right| \le C\left\{ {{h_i} + {1 \over \varepsilon }\int\limits_{{t_{i - 1}}}^{{t_{i + 1}}} {e^{ - {{\alpha t} \over \varepsilon }}}dt} \right\},{\kern 1pt} \;{\kern 1pt} 1 \le i \le {N_0} in which $h_{i} = {\begin{matrix} h_{p}^{(1)}, (p - 1) N \leq i \leq (p - 1 / 2) N \\ h_{p}^{(2)}, (p - 1 / 2) N + 1 \leq i \leq pN . \end{matrix}$ {h_i} = \left\{ {\matrix{ {h_p^{(1)},{\kern 1pt} \;{\kern 1pt} (p - 1)N \le i \le (p - 1/2)N} \cr {h_p^{(2)},{\kern 1pt} \;{\kern 1pt} (p - 1/2)N + 1 \le i \le pN.} \cr } } \right.

At the each submesh ω_{N, p}, we estimate the truncation error R as follows. We consider first the case σ_p = r_p−1 + r/2, and so r/2 ≤ α⁻¹θ_pɛ lnN and $h_{p}^{(1)} = h_{p}^{(2)} = τ_{p} = r / N$ h_p^{(1)} = h_p^{(2)} = {\tau _p} = r/N .

Hereby, since $\begin{matrix} h_{i}^{- 1} ε^{- 1} \int_{t_{i - 1}}^{t_{i + 1}} e^{- \frac{α t}{ε}} dt \leq ε^{- 1} h_{p} \leq \frac{2 β_{p} ln N}{α r} \frac{r}{N} = 2 β^{- 1} β_{p} N^{- 1} ln N, \\ (p - 1) N \leq i \leq pN, 1 \leq p \leq m - 1 \end{matrix}$ \matrix{ {h_i^{ - 1}{\varepsilon ^{ - 1}}\int\limits_{{t_{i - 1}}}^{{t_{i + 1}}} {e^{ - {{\alpha t} \over \varepsilon }}}dt \le {\varepsilon ^{ - 1}}{h_p} \le {{2{\beta _p}\ln N} \over {\alpha r}}{r \over N} = 2{\beta ^{ - 1}}{\beta _p}{N^{ - 1}}\ln N,} \cr {(p - 1)N \le i \le pN,{\kern 1pt} \;{\kern 1pt} 1 \le p \le m - 1}} it follows from (32) that $| R_{i} | \leq {CN}^{- 1} ln N, (p - 1) N \leq i \leq pN, 1 \leq p \leq m - 1 .$ \left| {{R_i}} \right| \le C{N^{ - 1}}\ln N,{\kern 1pt} \;{\kern 1pt} (p - 1)N \le i \le pN,{\kern 1pt} \;{\kern 1pt} 1 \le p \le m - 1.

We now consider the case σ_p = r_p−1 + α⁻¹β_pɛ lnN and estimate R_i on [r_p−1, σ_p] and [σ_p, r_p] separately. In the layer region [r_p−1, σ_p], the inequality (32) reduces to $\begin{matrix} | R_{i} | \leq C (1 + ε^{- 1}) {hp}^{(1)} = C (1 + ε^{- 1}) \frac{α^{- 1} β_{p} ε ln N}{N / 2}, (p - 1) N \leq i \leq (p - 1 / 2) N, \\ 1 \leq p \leq m - 1 . \end{matrix}$ \matrix{ {\left| {{R_i}} \right| \le C(1 + {\varepsilon ^{ - 1}})h{p^{(1)}} = C(1 + {\varepsilon ^{ - 1}}){{{\alpha ^{ - 1}}{\beta _p}\varepsilon \ln N} \over {N/2}},(p - 1)N \le i \le (p - 1/2)N,} \cr {1 \le p \le m - 1.}}

Hence $| R_{i} | \leq {CN}^{- 1} ln N, (p - 1) N \leq i \leq (p - 1 / 2) N, 1 \leq p \leq m - 1 .$ \left| {{R_i}} \right| \le C{N^{ - 1}}\ln N,{\kern 1pt} \;{\kern 1pt} (p - 1)N \le i \le (p - 1/2)N,{\kern 1pt} \;{\kern 1pt} 1 \le p \le m - 1.

Next we estimate R_i for (p−1/2)N +1 ≤ i ≤ pN. In this case, recalling that $t_{i} = σ_{p} + (i - (p - 1 / 2) N) h_{p}^{(2)} = r_{p - 1} + α^{- 1} β_{p} ε ln N + (i - (p - 1 / 2) N) h_{p}^{(2)}$ {t_i} = {\sigma _p} + (i - (p - 1/2)N)h_p^{(2)} = {r_{p - 1}} + {\alpha ^{ - 1}}{\beta _p}\varepsilon \ln N + (i - (p - 1/2)N)h_p^{(2)} , we obtain from (32) $| R_{i} | \leq C {h_{p}^{(2)} + α^{- 1} β_{p} (e^{- \frac{α (t_{i - 1})}{β_{p} ε}} - e^{- \frac{α (t_{i})}{β_{p} ε}})}$ \left| {{R_i}} \right| \le C\left\{ {h_p^{(2)} + {\alpha ^{ - 1}}{\beta _p}({e^{ - {{\alpha ({t_{i - 1}})} \over {{\beta _p}\varepsilon }}}} - {e^{ - {{\alpha ({t_i})} \over {{\beta _p}\varepsilon }}}})} \right\} and this implies that $| R_{i} | \leq {CN}^{- 1} .$ \left| {{R_i}} \right| \le C{N^{ - 1}}.

Then, we estimate the remainder term r⁽⁰⁾. From the explicit expression (32) and |ϕ₀ (t)| ≤ 1, after similar calculation as above, it follows that $| r^{(0)} | \leq {CN}^{- 1} ln N .$ \left| {{r^{\left( 0 \right)}}} \right| \le C{N^{ - 1}}\ln N.

Combining the two previous lemmas gives us the following main convergence result.

Theorem 4.1

Let a ∈ C¹ (Ī), ϕ ∈ C¹ (I₀). The continuously differentiable function f (t,u,v) satisfies the conditions (1.4) and the derivative $\frac{\partial}{\partial t} f (t, u, v)$ {\partial \over {\partial t}}f(t,u,v) is bounded for 0 ≤ t ≤ T and |u|, |v| ≤ C. Then the following estimate holds $| y_{i} - u_{i} | \leq {CN}^{- 1} ln N, 0 \leq i \leq N_{0} .$ \left| {{y_i} - {u_i}} \right| \le C{N^{ - 1}}\ln N,\quad 0 \le i \le {N_0}.

Numerical Results

In this section, a simple numerical example is devised to verify the validity for the proposed method. Consider the test problem with $a = 1, f = - u (t - 1) + t^{2}, T = 2, φ = 1 + t, A = - 1 .$ a = 1,f = - u(t - 1) + {t^2},T = 2,\varphi = 1 + t,A = - 1.

We use the double mesh principle to estimate the errors and compute the experimental rates of convergence in our computed solution, i.e. We compare the computed solution with the solution on a mesh that is twice as fine [17]. The error estimate $e_{ε}^{N, p}$ e_\varepsilon ^{N,p} and the computed convergence rate $r_{ε}^{N, p}$ r_\varepsilon ^{N,p} obtained in this way are denoted by $e_{ε}^{N, p} = max_{ω_{N, p}} | y^{ε, N} - y^{ε, 2 N} |, p = 1, 2; r_{ε}^{N, p} = ln (e_{ε}^{N, p} / e_{ε}^{2 N, p}) / ln 2 .$ e_\varepsilon ^{N,p} = \mathop {\max }\limits_{{\omega _{N,p}}} \left| {{y^{\varepsilon ,N}} - {y^{\varepsilon ,2N}}} \right|,p = 1,2;r_\varepsilon ^{N,p} = \ln \left( {e_\varepsilon ^{N,p}/e_\varepsilon ^{2N,p}} \right)/\ln 2. The resulting errors $e_{ε}^{N, p}$ e_\varepsilon ^{N,p} and the corresponding numbers $r_{ε}^{N, p}$ r_\varepsilon ^{N,p} for particular values of ɛ, N, for first and second subinterval are listed in the Tables 1 and 2.

Table 1

Maximum errors and rates of convergence on ω_N,1

ɛ	N = 64	N = 128	N = 256	N = 512	N = 1024
2⁻¹	0.0085116	0.0042832	0.0021485	0.0010760	0.0005384
2⁻¹	0.99	1.00	1.00	1.00
2⁻²	0.0140007	0.0070930	0.0035693	0.0017090	0.0008906
2⁻²	0.98	0.98	0.99	1.00
2⁻⁴	0.0241181	0.0143603	0.0083166	0.0047146	0.0026331
2⁻⁴	0.75	0.78	0.81	0.84
2⁻⁶	0.0230541	0.0137267	0.0079497	0.0045066	0.0025149
2⁻⁶	0.75	0.78	0.81	0.84
2⁻⁸	0.0227881	0.0135684	0.0078580	0.0044546	0.0024859
2⁻⁸	0.75	0.78	0.81	0.84
2⁻¹⁰	0.0227216	0.0135288	0.0078350	0.0044416	0.0024786
2⁻¹⁰	0.75	0.78	0.81	0.84
2⁻¹²	0.0227050	0.0135189	0.0078293	0.0044382	0.0024768
2⁻¹²	0.75	0.78	0.81	0.84
2⁻¹⁴	0.0227008	0.0135164	0.0078279	0.0044361	0.0024763
2⁻¹⁴	0.75	0.78	0.81	0.84
2⁻¹⁶	0.0226998	0.0135158	0.00782756	0.0044374	0.0024762
2⁻¹⁶	0.75	0.78	0.81	0.84

Table 2

Maximum errors and rates of convergence on ω_N,2

ɛ	N = 64	N = 128	N = 256	N = 512	N = 1024
2⁻¹	0.0079814	0.0040103	0.0020101	0.0010063	0.0005033
2⁻¹	0.98	0.99	1.00	1.00
2⁻²	0.0012138	0.0061392	0.0030875	0.0015482	0.0007775
2⁻²	0.98	0.99	1.00	1.00
2⁻⁴	0.0202600	0.0120754	0.0069885	0.0039609	0.0022110
2⁻⁴	0.75	0.79	0.81	0.84
2⁻⁶	0.0194243	0.0115426	0.0066802	0.0037862	0.0021126
2⁻⁶	0.75	0.78	0.82	0.84
2⁻⁸	0.0191427	0.0114094	0.0066031	0.0037425	0.0020882
2⁻⁸	0.75	0.78	0.82	0.84
2⁻¹⁰	0.0190868	0.0113761	0.0065838	0.0037315	0.0020821
2⁻¹⁰	0.75	0.78	0.82	0.84
2⁻¹²	0.0190729	0.0113678	0.0065790	0.0037288	0.0020806
2⁻¹²	0.75	0.78	0.82	0.84
2⁻¹⁴	0.0190696	0.0113657	0.0065778	0.0037281	0.0020668
2⁻¹⁴	0.75	0.78	0.82	0.84
2⁻¹⁶	0.0190685	0.0113652	0.0065754	0.0037280	0.0020801
2⁻¹⁶	0.75	0.78	0.82	0.84

The values of ɛ for which we solve the test problem are ɛ = 2⁻ⁱ,i = 1,2,4,...,16. These convergence rates are increasing as N increases for any fixed ɛ. Tables 1 and 2 verify the ɛ-uniform convergence of the numerical solution on both subintervals and computed rates are essentially in agreement with our theoretical analysis.

Conclusions

A numerical method is developed to solve singularly perturbed initial value problem for semilinear second-order delay differential equation. This method is based on exponentially fitted difference scheme on an equidistant mesh, which gives first order uniform convergence in the discrete maximum norm. The method is shown to be uniformly convergent to the continuous solution i.e., independent of mesh parameter and perturbation parameter. Efficiency of the present method is demonsrated by a numerical example and also by comparing the results with exact solution of the problem. Presented numerical results are essentially in agreement with our theoretical analysis.

eISSN:: 2444-8656
Language:: English

Publication timeframe:: Volume Open
Journal Subjects:: Life Sciences, other, Mathematics, Applied Mathematics, General Mathematics, Physics

Journal RSS Feed

A finite difference method on layer-adapted mesh for singularly perturbed delay differential equations

Published Online: Mar 31, 2020

Page range: 425 - 436

Received: Jun 11, 2019

Accepted: Aug 13, 2019

DOI: https://doi.org/10.2478/amns.2020.1.00040

KeywordsDelay differential equation, Singular perturbation, Finite difference scheme, Piecewise-uniform mesh, Error estimates

© 2019 Fevzi Erdogan et al., published by Sciendo

This work is licensed under the Creative Commons Attribution 4.0 International License.

Keywords
Delay differential equation, Singular perturbation, Finite difference scheme, Piecewise-uniform mesh, Error estimates