A NON-ASYMPTOTIC METHOD FOR SINGULARLY PERTURBED DELAY DIFFERENTIAL EQUATIONS
A NON-ASYMPTOTIC METHOD FOR SINGULARLY PERTURBED DELAY DIFFERENTIAL EQUATIONS
Journal of Applied Mathematics & Informatics. 2014. Jan, 32(1_2): 39-53
• Received : April 30, 2012
• Accepted : August 13, 2013
• Published : January 28, 2014
PDF
e-PUB
PPT
Export by style
Share
Article
Author
Metrics
Cited by
TagCloud
GEMECHIS FILE
Y.N. REDDY

Abstract
In this paper, a non-asymptotic method is presented for solving singularly perturbed delay differential equations whose solution exhibits a boundary layer behavior. The second order singularly perturbed delay differential equation is replaced by an asymptotically equivalent first order neutral type delay differential equation. Then, Simpson’s integration formula and linear interpolation are employed to get three term recurrence relation which is solved easily by Discrete Invariant Imbedding Algorithm. Some numerical examples are given to validate the computational efficiency of the proposed numerical scheme for various values of the delay and perturbation parameters. AMS Mathematics Subject Classification : 65L10.
Keywords
1. Introduction
A singularly perturbed delay differential equation is an ordinary differential equation in which the highest derivative is multiplied by a small parameter and containing delay term. In recent years, Delay differential equations have attracted the attention of many researchers because of their applications in many scientific and technical fields such as in population ecology, control theory, viscoelasticity, and materials with thermal memory, El’sgolts [4] , and Bellman and Cooke [3] . Boundary value problems involving delay differential equations are also used to study signal transmission with time delays in control theory El’sgolts [4] , first-exit problems in neurobiology Stein [13] , and Tuckwell [14] . As a result, many researchers tried to develop and present numerical schemes for solving such problems. For example, Amiraliyev and Cimen [1] , presented numerical method comprising a fitted difference scheme on a uniform mesh to solve second order delay differential equations. Lange and Miura [10 , 11] gave an asymptotic approach for a class of boundary-value problems for linear second-order differential-difference equations. Kadalbajoo and Sharma [7 , 8 , 9] presented numerical approaches to solve singularly perturbed differential-difference equations, which contains negative shift in the convention term (i.e. in the derivative term). Furthermore, Phaneendra, et al. [12] , presented higher order numerical method to solve the problems of the same type.
The objective of this paper is, therefore, to present a non-asymptotic method for solving singularly perturbed delay differential equations. The method is named as non-asymptotic as it does not depend on asymptotic expansion as well as on the matching of the coefficients. In this method, we first replace the second order singularly perturbed delay differential equation by an asymptotically equivalent first order neutral type delay differential equation and employ Simpson’s rule for the numerical integration of the obtained first order delay differential equation. Then, linear interpolation is used to get the three term recurrence relation which is solved easily by discrete Invariant Imbedding algorithm. The method is demonstrated by implementing several model examples by taking various values for the delay and perturbation parameters.
2. Description of the Method
Consider singularly perturbed singular boundary value problems of the form
PPT Slide
Lager Image
with the interval and boundary conditions
PPT Slide
Lager Image
PPT Slide
Lager Image
where ε is small parameter, 0 < ε << 1, and the delay parameter, δ is of o ( ε ) satisfying ( ε δ a ( x )) > 0 for all [0, 1]; a ( x ), b ( x ), f ( x ) and ϕ ( x ) are smooth functions and β is a finite constant. For δ = 0, the problem (1)-(3) becomes a boundary value problem for singularly perturbed ordinary differential equations. The layer behavior of the problem under consideration is maintained only for δ ≠ 0 but sufficiently small (i.e δ is of o ( ε )). When the delay parameter, δ exceeds the perturbation parameter, ε (i.e δ is of O ( ε )), the layer behavior of the solution is no longer maintained, rather the solution exhibits an oscillatory behavior or diminished. In this paper, we consider the cases where the boundary layer is either on the left side or right side of the interval [0, 1].
2.1. Layer on the Left Side: In general, the solution of problem (1)-(3) exhibits boundary layer behavior at one end of the interval [0, 1] depending on the sign of a ( x ). Now, we assume that a ( x ) ≥ M > 0 throughout the interval [0, 1], where M is positive constant. This assumption merely implies that the boundary layer will be in the neighborhood of x = 0. Further, it is also assumed that b ( x ) satisfies the condition that b ( x ) ≤ − θ < 0 , ∀ [0, 1] where θ is a positive constant.
By using Taylor series expansion in the neighborhood of the point x , we have
PPT Slide
Lager Image
By using (4) in (1) we get an asymptotically equivalent first order delay neutral type differential equation:
PPT Slide
Lager Image
The transition from equation (1) to equation (5) is admitted, because of the condition that ε is small, 0 < ε << 1; and the truncated term ε 2 /2 is very small and approaches to 0 as ε is small. Further details on the validity of this transition can also be found in Elsgolt’s and Norkin [5] .
Now we divide the interval [0, 1] into N equal subintervals of mesh size h = 1/ N so that xi = ih , i = 0, 1, 2, · · · , N .
Integrating equation (5) with respect to x from xi to x i+1 for i = 1, 2, · · · , N −1, we get
PPT Slide
Lager Image
where yi = y ( xi ), ai = a ( xi ), bi = b ( xi ), fi = f ( xi ).
By using Simpson’s rule to evaluate the integral, we get:
PPT Slide
Lager Image
Since delay term is of the small order of ε , to tackle the term containing delay, we can use Taylor’s approximation (Kadalbajoo and Sharma [9] , Lange and Miura [11] ). Thus, by means of Taylor series expansion and then by approximating
PPT Slide
Lager Image
by linear interpolation, we get:
PPT Slide
Lager Image
PPT Slide
Lager Image
PPT Slide
Lager Image
PPT Slide
Lager Image
In similar way,
PPT Slide
Lager Image
Hence, by making use of equations (8)-(12) in (7) we obtain:
PPT Slide
Lager Image
To make equation (13) a three term recurrence relation, we can express yi +1/2 in terms of yi −1 , yi and yi +1 using Hermite’s interpolation as follows:
PPT Slide
Lager Image
In view of equations(5) and (14), we get:
PPT Slide
Lager Image
By making use equations (8)- (12) in (15) and finite difference approximations, we get
PPT Slide
Lager Image
Finally, making use of equations (16) in (13) and rearranging as three term recurrence relation, we get:
PPT Slide
Lager Image
for i = 1, 2, · · · , N − 1, where
PPT Slide
Lager Image
PPT Slide
Lager Image
PPT Slide
Lager Image
PPT Slide
Lager Image
We solve the above tridiagonal system by using method of invariant imbedding system described in section 3.
2.2. Right-End Boundary Layer Problem: We now assume that a ( x ) ≤ − M < 0 throughout the interval [0, 1], where M is positive constant. This assumption merely implies that the boundary layer will be in the neighborhood of x = 1. By using Taylor series expansion in the neighborhood of the point x , we have:
PPT Slide
Lager Image
and consequently, equation (1) is replaced by the following first order differential equation:
PPT Slide
Lager Image
Now we divide the interval [0, 1] into N equal subintervals of mesh size h = 1/ N so that xi = ih , i = 0, 1, 2, · · · , N . Integrating equation (23) with respect to x from xi −1 to xi , we get:
PPT Slide
Lager Image
By using Simpson’s rule to evaluate the integral, we get:
PPT Slide
Lager Image
Again, by means of Taylor series expansion and then by approximating
PPT Slide
Lager Image
by linear interpolation, we get:
PPT Slide
Lager Image
PPT Slide
Lager Image
PPT Slide
Lager Image
PPT Slide
Lager Image
In similar way,
PPT Slide
Lager Image
Hence, by making use of equations (26)-(30) in (25) we obtain:
PPT Slide
Lager Image
To make equation (31) a three term recurrence relation, we can express yi −1/2 in terms of yi −1 , yi and yi +1 using Hermite’s interpolation as follows:
PPT Slide
Lager Image
In view of equations (23) and (32), we get:
PPT Slide
Lager Image
By making use of equations (26)- (30) in (33) and finite difference approximations, we get:
PPT Slide
Lager Image
Finally, making use of equations (34) in (31) and rearranging as three term recurrence relation, we get:
PPT Slide
Lager Image
for i = 1, 2, · · · , N − 1, where
PPT Slide
Lager Image
PPT Slide
Lager Image
PPT Slide
Lager Image
PPT Slide
Lager Image
We solve the above tridiagonal system by using method of invariant imbedding described in the next section.
3. Discrete Invariant Imbedding Algorithm
We now describe the Thomas algorithm which is also called Discrete Invariant Imbedding (Angel and Bellman [2] ) to solve the three term recurrence relation:
PPT Slide
Lager Image
for i = 1, 2, · · · , N − 1. Let us set a difference relation of the form:
PPT Slide
Lager Image
where Wi = W ( xi ) and Ti = T ( xi ) are to be determined. From (41), we have:
PPT Slide
Lager Image
substituting (42) in (40), we get:
PPT Slide
Lager Image
By comparing (43) and (41), we get the recurrence relations:
PPT Slide
Lager Image
PPT Slide
Lager Image
To solve these recurrence relations for i = 1, 2, 3, · · · , N −1 , we need the initial conditions for W 0 and T 0 . If we choose W 0 = 0 , then we get T 0 = ϕ (0) . With these initial values, we compute Wi and Ti for i = 1, 2, 3, · · · , N − 1 from equations (44) and (45) in forward process, and then obtain yi in the backward process from (41).
The conditions for the stability of discrete invariant imbedding algorithm are (Angel and Bellman [2] , Kadalbajoo and Reddy [6] ):
PPT Slide
Lager Image
In our method, one can easily show that if the assumptions a ( x ) > 0, b ( x ) < 0 and ( ε δa ( x )) > 0 hold, then the above conditions (46) hold and thus the invariant imbedding algorithm is stable.
4. Numerical Examples
To demonstrate the applicability of the method, we have considered two numerical experiments with left-end boundary layer and two numerical experiments with right-end boundary layer. We compared the computed results with the exact solution of the problems where the exact solution of the problems is known; and we also have tested the effect of small delay parameter on computed solution of the problem for different values of δ of o ( ε ) which are presented by numerical experiments whose exact solutions are not known.
Example 4.1 . Consider an example of singularly perturbed delay differential equation with left layer:
PPT Slide
Lager Image
The exact solution is given by
PPT Slide
Lager Image
where
PPT Slide
Lager Image
The computational results are presented in the Tables 1 , 2 , 3 and 4 for ε = 0.01 and 0.001 for different values of δ .
Example 4.2. Now we consider an example of variable coefficient singularly perturbed delay differential equation with left layer:
PPT Slide
Lager Image
For which the exact solution is not known.This example is also considered to show the effect of the small shift on the boundary layer solution. The computational results are presented in the Tables 5 and 6 for ε = 0.01 and 0.001 for different values of δ .
Example 4.3. Consider a singularly perturbed delay differential equation with right layer:
PPT Slide
Lager Image
The exact solution is given by
PPT Slide
Lager Image
where
PPT Slide
Lager Image
The computational results are presented in the Tables 7 and 8 for ε = 0.01 and 0.001 for different values of δ .
Example 4.4. Now we consider an example of variable coefficient singularly perturbed delay differential equation with right layer:
PPT Slide
Lager Image
For which the exact solution is not known. The computational results are presented in the Tables 9 and 10 for ε = 0.01 and 0.001 for different values of δ .
5. Discussion and Conclusion
We have presented a numerical integration method to solve singularly perturbed delay differential equations whose solutions exhibits the layer behavior. The scheme is repeated for different choices of the delay parameter, δ and perturbation parameter ε . The choice of δ is not unique, but can assume any number of values satisfying the condition, δ ( ε ) = τε with τ = O (1) and τ is not too large Lange and Miura [10] .
To demonstrate the efficiency of the method, some numerical examples of boundary value problems with constant and variable coefficients are considered for different values of ε and δ . Besides solving for the numerical solutions, we have considered numerical results of boundary value problems, examples 4.2 & 4.4 , whose exact solutions are not known to show the effect of delay parameter on the boundary layer solutions.
From the tabulated solutions of the numerical experiments considered, it can be observed that the present method approximates the exact solution very well ( Tables 1 - 4 , 7 - 8 ) for h ε for which other classical finite difference methods fails to give good results; and the small shift, δ affects both the boundary layer solutions (left and right) in similar fashion but reversely. That is, as δ increases the thickness of the left boundary layer decreases ( Tables 5 - 6 ) while that of the right boundary layer increases ( Tables 9 - 10 ). This method does not depend on asymptotic expansion as well as on the matching of the coefficients. Thus, we have devised an alternative technique of solving boundary value problems for singularly perturbed delay differential equations, which is easily implemented on computer and is also practical.
Numerical results of Example 4.1 forε= 0.01,δ= 0.001, N = 100
PPT Slide
Lager Image
Numerical results of Example 4.1 for ε = 0.01, δ = 0.001, N = 100
Numerical results of Example 4.1 forε= 0.01,δ= 0.008, N = 100
PPT Slide
Lager Image
Numerical results of Example 4.1 for ε = 0.01, δ = 0.008, N = 100
Numerical results of Example 4.1 forε= 0.001,δ= 0.0001, N = 100
PPT Slide
Lager Image
Numerical results of Example 4.1 for ε = 0.001, δ = 0.0001, N = 100
Numerical results of Example 4.1 forε= 0.001,δ= 0.0008, N = 100
PPT Slide
Lager Image
Numerical results of Example 4.1 for ε = 0.001, δ = 0.0008, N = 100
Numerical results of Example 4.2 forε= 0.01,N= 100 and different values ofδ.
PPT Slide
Lager Image
Numerical results of Example 4.2 for ε = 0.01, N = 100 and different values of δ.
Numerical results of Example 4.2 forε= 0.001,N= 100 and different values ofδ.
PPT Slide
Lager Image
Numerical results of Example 4.2 for ε = 0.001, N = 100 and different values of δ.
Numerical results of Example 4.3 forε= 0.01,δ= 0.001, N = 100
PPT Slide
Lager Image
Numerical results of Example 4.3 for ε = 0.01, δ = 0.001, N = 100
Numerical results of Example 4.3 forε= 0.01,δ= 0.008, N = 100
PPT Slide
Lager Image
Numerical results of Example 4.3 for ε = 0.01, δ = 0.008, N = 100
Numerical results of Example 4.4 forε= 0.01,N= 100 and different values ofδ.
PPT Slide
Lager Image
Numerical results of Example 4.4 for ε = 0.01, N = 100 and different values of δ.
Numerical results of Example 4.4 forε= 0.001,N= 100 and different values ofδ.
PPT Slide
Lager Image
Numerical results of Example 4.4 for ε = 0.001, N = 100 and different values of δ.
BIO
Gemechis File received his B.Sc from Bahir Dar University, M.Sc from Addis Ababa University, and currently persuing Ph.D. from National Institute of Technology, Warangal. His research interest focuses on the numerical treatments of singular perturbation problems.
Department of Mathematics, National Institute of Technology, Warangal, 506 004 A.P, India.
e-mail: gammeef@yahoo.com
Y. N. Reddy received his M.Sc from NIT Warangal, and Ph.D. from IIT-Kampur. His research interest focuses on the numerical treatments of singular perturbation problems.
Department of Mathematics, National Institute of Technology, Warangal, 506 004 A.P, India.
e-mail: ynreddy@nitw.ac.in
References
Amiraliyev G. M. , Cimen E. (2010) Numerical Method For A Singularly Perturbed Convection Diffusion Problem with Delay Applied Mathematics and Computation 216 2351 - 2359    DOI : 10.1016/j.amc.2010.03.080
Angel E. , Bellman R. 1972 Dynamic Programming and Partial differential equations Academic Press New York
Bellman R. , Cooke K. L. 1963 Differential Difference Equations Academic Press New York
El’sgolts L.E. 1964 Qualitative Methods in Mathematical Analysis Trans. Math. Monogr. American Mathematical Society, Providence, R. I. 12
Elsgolt’s L. E. , Norkin S. B. 1973 Introduction to the Theory and Applications of Differential Equations with Deviating Arguments Academic Press New York
Kadalbajoo M. K. , Reddy Y. N. (1986) A non asymptotic method for general linear singular perturbation problems Journal of Optimization Theory and Applications 55 256 - 269
Kadalbajoo M. K. , Sharma K. K. (2005) Numerical treatment of boundary value problems for second-order singularly perturbed delay differential equations Comput. Appl. Math. 24 151 - 172    DOI : 10.1590/S0101-82052005000200001
Kadalbajoo M. K. , Sharma K. K. (2008) A numerical method on finite difference for boundary value problems for singularly perturbed delay differential equations Applied Mathematics and Computation 197 692 - 707    DOI : 10.1016/j.amc.2007.08.089
Kadalbajoo M. K. , Sharma K. K. (2004) Numerical analysis of singularly perturbed delay differential equations with layer behavior Applied Mathematics and Computation 157 11 - 28    DOI : 10.1016/j.amc.2003.06.012
Lange C. G. , Miura R. M. (1994) Singular perturbation analysis of boundary-value problems for differential-difference equations. v. small shifts with layer behavior SIAM J. Appl. Math. 54 249 - 272    DOI : 10.1137/S0036139992228120
Lange C. G. , Miura R. M. (1994) Singular perturbation analysis of boundary-value problems for differential-difference equations. vi. Small shifts with rapid oscillations SIAM J. Appl. Math. 54 273 - 283    DOI : 10.1137/S0036139992228119
Phaneendra K. , Reddy Y. N. , Soujanya GBSL. (2011) A seventh order numerical method for singularly perturbed differential difference equations with negative shift Non linear Analysis: Modeling and Control 16 206 - 219
Stein R. B. (1965) A theoretical analysis of neuronal variability Biophys. J. 5 173 - 194    DOI : 10.1016/S0006-3495(65)86709-1
Tuckwell H. C. (1976) On the first exit time problem for temporally homogeneous Markov processes J. Appl. Probab. 13 39 - 46    DOI : 10.2307/3212663