CC BY-NC-ND
0Int. J. Appl. Math. Comput. Sci, 2014, Vol. 24, No. 2, 387-395
DOI: 10.2478/amcs-2014-0029 cITHCS
A ROBUST COMPUTATIONAL TECHNIQUE FOR A SYSTEM OF SINGULARLY PERTURBED REACTION-DIFFUSION EQUATIONS
Vinod KUMAR *, RajeshK. BAWA **, Arvind K. LAL *
* School of Mathematics and Computer Applications Thapar University, Patiala, 147004, India e-mail: vinod.patiala@gmail.com
"Department of Computer Science Punjabi University, Patiala 147002, India
In this paper, a singularly perturbed system of reaction-diffusion Boundary Value Problems (BVPs) is examined. To solve such a type of problems, a Modified Initial Value Technique (MIVT) is proposed on an appropriate piecewise uniform Shishkin mesh. The MIVT is shown to be of second order convergent (up to a logarithmic factor). Numerical results are presented which are in agreement with the theoretical results.
Keywords: asymptotic expansion approximation, backward difference operator, trapezoidal method, piecewise uniform Shishkin mesh.
1. Introduction
In many fields of applied mathematics we often come across initial/boundary value problems with small positive parameters. If, in a problem arising in this manner, the role of the perturbation is played by the leading terms of the differential operator (or part of them), then the problem is called a Singularly Perturbed Problem (SPP). Applications of SPPs include boundary layer problems, WKB theory, the modeling of steady and unsteady viscous flow problems with a large Reynolds number and convective-heat transport problems with large Peclet numbers, etc.
The numerical analysis of singularly perturbed cases has always been far from trivial because of the boundary layer behavior of the solution. These problems depend on a perturbation parameter e in such a way that the solutions behave non-uniformly as e tends towards some limiting value of interest. Therefore, it is important to develop some suitable numerical methods whose accuracy does not depend on e, i.e., which are convergent e-uniformly. There are a wide variety of techniques to solve these types of problems (see the books of Doolan et al. (1980) and Roos et al. (1996) for further details). Parameter-uniform numerical methods for a scalar reaction-diffusion equation have been examined
extensively in the literature (see the works of Roos et al. (1996), Farrell et al. (2000), Miller et al. (1996) and the references therein), whereas for a system of singularly perturbed reaction-diffusion equations only few results (Madden and Stynes, 2003; Matthews et al., 2000; 2002; Natesan and Briti, 2007; Valanarasu and Ramanujam, 2004) have been reported.
In this paper, we treat the following system of two singularly perturbed reaction-diffusion equations:
L\u = -eu"(x) + an(x)ui (x) + a^(x)u2(x)
= fi(x), (1)
L2U = -¡u'^x) + a2i(x)ui (x) + a22 (x)u2(x)
= f2(x), (2)
where u = (u1,u2)T, x e Q = (0,1), with the boundary conditions
u(0) = Q , u(1) = Q . (3)
Without loss of generality, we shall assume that 0 <
e < l < 1. The functions a11(x), a12 (x), a21(x), a22(x), f1(x), f2(x) are sufficiently smooth and satisfy the following set of inequalities:
(i) on(x) > |oi2(x)|, a22(x) > |o2i(x)|,
x € Q = [0,1],
(ii) oi2(x) < 0, «21 (x) < 0, x G Q.
Shishkin (1995) classifies three separate cases for a system of two singularly perturbed reaction-diffusion problems with diffusion coefficients e, f: (i) 0 < e = f C 1, (ii) 0 < e C f = 1 and (iii) e, f arbitrary. Matthews et al. (2000) consider case (i), showing that a standard finite difference scheme is uniformly convergent on a fitted piecewise uniform mesh. They establish first-order convergence up to a logarithmic factor in the discrete maximum norm. The same authors have also obtained a similar result for case (ii), which they have strengthened to show almost second-order convergence (Matthews et al, 2002). Madden and Stynes (2003) obtained almost first-order convergence for the general case (iii). For case (ii), Natesan and Briti (2007) developed a numerical method which is a combination of a cubic spline and a finite difference scheme.
Das and Natesan (2013) obtained almost second-order convergence for the general case (iii) in which they used central difference approximation for an outer region with cubic spline approximation for an inner region of boundary layers. Melenk et al. (2013) have constructed full asymptotic expansions together with error bounds that cover the complete range of 0 < e C f C 1. Rao et al. (2011) proposed a hybrid difference scheme on a piecewise-uniform Shishkin mesh and showed that the scheme generates better approximations to the exact solution than the classical central difference one. Valanarasu and Ramanujam (2004) proposed an Asymptotic Initial Value Method (AIVM) to solve (1)—(3), whose theoretical order of convergence is 1. Bawa et al. (2011) used a hybrid scheme for a singularly perturbed delay differential equation, which is of second order convergent.
We construct a Modified Initial Value Technique (MIVT) for (1)—(3) which is based on the underlying idea of the AIVM (Valanarasu and Ramanujam, 2004). The aim of the present study is to improve the order of convergence to almost second order (up to a logarithmic factor) for case (i), i.e., for 0 < e = f C 1.
First, in this technique, an asymptotic expansion approximation for the solution of the Boundary Value Problem (BVP) (1)—(3) has been constructed. Then, Initial Value Problems (IVPs) and Terminal Value Problems (TVPs) are formulated whose solutions are the terms of this asymptotic expansion. The IVPs and TVPs are happened to be SPPs, and therefore they are solved by a hybrid scheme similar to that by Bawa et al. (2011). The scheme is a combination of the trapezoidal scheme and a backward difference operator. It not only retains the oscillation free behavior of the backward difference
operator but also retains the second order of convergence of the trapezoidal method.
The paper is organized as follows. Section 2 presents an asymptotic expansion approximation of (1)—(3). The initial value problem is discussed in Section 3. Section 4 deals with the error estimates of the proposed hybrid scheme. The Shishkin mesh and the MIVT are given in Section 5. Finally, numerical examples are presented in Section 6 to illustrate the applicability of the method. The paper ends with some conclusions.
Note. Throughout this paper, we let C denote a generic positive constant that may take different values in the different formulas, but is always independent of N and e. Here || • || denotes the maximum norm over Q.
2. Preliminaries
2.1. Maximum principle and the stability result.
Lemma 1. (Matthews et al., 2002) Consider the BVP system (l)-(3). IfLry > 0,_L2y > 0 in Q andy{0) > 0, y{ 1) > 0, then y{x) > 0 in Q.
Lemma 2. (Matthews et al., 2002) If y(x) is the solution ofBVP(1)-(3), then
II y(x) ll<- II /11 + 11 №) 11 + 11^1)11,
where 7 = min{an(x) + ai2(x), a2\(x) + a22(x)}.
2.2. Asymptotic expansion approximation. It is well known that, by using the fundamental idea of WKB (Valanarasu and Ramanujam, 2004; Nayfeh, 1981), an asymptotic expansion approximation for the solution of the BVP (1)-(3) can be constructed as
uas{x) = ur{x) + v(x) + 0{yfc), Url(x)
Ur(X) =
yUR2(x)
is the solution of the reduced problem of (1)—(3) and is given by
- ai2(x)uR2(x) = fi(x), (4)
- a22(x)uR2(x) = f2(x), (5)
an (x)uri(x) a2i (x)uri (x) x e [0,1), and
is given by
vi(x) = [p - uri(0)]
+ [q - uri(1)]
vi (x) V2 (x)
aii (0) + ai2 (0)
aii (x) + ai2 (x)_ aii(1) + ai2(1)
aii(x) + ai2(x)_
VLi(x)
Wri(x),
V2(x) =[r - UR2(0)]
[s - Ur2(l)]
o21(0) + o22(0) 0'21 (x) + a22 (x) a2i(1) + a22(1)
a2i (x) + a,22 (x)
Vli(x)
Vl2(x)
Wr2(x).
\VL2(x)
is a "left boundary layer correction" and
W (x) = fwR1(x)\ W R(x)= {wR2(x))
is a "right boundary layer correction" defined as
VLi{x) = exp j-^ \f-vL2(x)=exp j-^
ivR1(x) = expj-^ \j
wm{x) = expj-^ \j
[aii(s) + ai2(s)]
[a2i(s) + a22(s)]
[aii(s) + ai2(s)]
[a2i(s) + a22(s)]
d^, (6)
ds}, (7)
ds}, (8)
ds}. (9)
It is easy to verify that vL1(x), vL2(x), wR1(x) and wR2(x) satisfy the following IVPs and TVPs, respectively:
yfëv'L1{x) + ^[anix) + a,12{x)]vL1(x) = 0, (10) Vli(0) = 1, (11)
\fiv'L2(x) + \/[02i (x) + a,22{x)}vL2(x) = 0, (12) Vl2(0) = 1, (13)
^w'R1(x) - y/[au{x) + oi2(x)]u>_Ri(x)
WRi(1) = 1,
0, (14) (15)
3. Initial value problem
In this section, we describe a hybrid scheme for the following singularly perturbed initial value problem of the first order:
LeV(x) = ey'(x) + b(x)y(x) = g(x), (18)
y(0) = A, (19)
where A is a constant, x e Q = (0,1) and 0 < e C 1 is a small parameter, b and g are sufficiently smooth functions, such that b(x) > /3 > 0 on SI = [0,1]. Under these assumptions, (18)-(19) possesses a unique solution y(x) (Doolan etal., 1980).
On Q, a piecewise uniform mesh of N mesh intervals is constructed as follows. The domain Cl is sub-divided into two subintervals [0, a] U [a, 1] for some a that satisfy 0 < a < 1/2. On each sub-interval, a uniform mesh with N/ 2 mesh intervals is placed. The interior points of the mesh are denoted by
x0 = 0, xi = ^2 hk, hk = xk+1 - xk,
xN = 1, i = 1, 2,...,N - 1.
Clearly, xn = 0.5 and Q^ = . It is fitted to
(18)-(19) by choosing a to be the following functions of N and e:
a = min | -, <Jne In N
where a0 > 2/p. Note that this is a uniform mesh when a = 1/2. Further, we denote the mesh size in the regions [0, a] by h = 2a/N and in [a, 1] by H = 2(1 - a)/N.
We define the following hybrid scheme for the approximation of (18)-(19):
Ln Yi =
l> V: + ■
_ 9i-1 + 9i
bi-iYi_1 + biYi
0 < i < 2
r> Y. ■ h :Y: ,J:. % < i < N,
Wr2(1) = 1.
Theorem 1. (Valanarasu and Ramanujam, 2004) The ze-roth order asymptotic expansion approximation uas satisfies the inequality
II (u-uas)(x) II < C^fe, where u(x) is the solution of the BVP (1)-(3).
Yo = A,
yßw'm{x) - yj[a2i{x) +o22(x)]u>ñ2(x) = 0, (16) where
D-Yi =
xi - xi i
and bi = b(xi), gi = g(xH).
4. Error estimate
Theorem 2. Let y(x) and Yi be respectively the solutions of(18)-(19) and (20)-(21). Then the local truncation er-
\LN (Yi - v(xí))\<{
ror satisfies the following bounds:
'CN-2ol ln2 N for 0 <i < N/2, C (N-1e + N) for N/2 < i < Nand H < e, C (N-2 + N) for N/2 < i < Nand H > e.
Proof. We distinguish several cases depending on the location of the mesh points. Firstly, we state the bound for the derivatives of the continuous solution, i.e., the solution y(x) of the IVP (18)-(19) satisfies the following bound (Doolan et al., 1980):
We discuss the following two cases. First, if H < e, from (27), we obtain
\LN(Yi - y(xi))\<CeH\y»(Z)\,
< C[He + He~i exp(-fX/e)]
< C[N-ie + N].
Secondly, if H > e, then using the bounds of the derivatives of y(x) from (23), one can obtain the following:
\LN(Yi - y(xi))\ < C ( He + JX i(xi - £)e_2 exp(-3£/e)d£^ . (30)
|y(k)(x)| < C [l + e-k exp(-px/e)] . (23) Integrating by parts, we get
For xi G (0, a], by using the usual Taylor series expansion, we get
\LN (Yi - y(xi))\<Ceh2\y"'(0|
for 0 < i < N/2 and some point £, xi_1 < £ < xi.
First we consider the case when the mesh is uniform. Then, a =1/2 and e_1 < Ca0 ln N. Using the above bound, we have
\LN(Yi - y(xi))\ < Ceh2 [1+e-3 exp(-££/e)]
< CN-2a¡2 ln2 N
for 0 <i < N/2.
Secondly, we consider the case when the mesh is non-uniform. Using h = 2N_1a0e ln N on the above bound and bounding the exponential function by a constant, we have
\LN(Yi - y(xi))\<CN-2a2 ln2 N
for 0 <i < N/2.
For xi G (a, 1], by using the Taylor series expansion, we get
\LN (Yi - y(xi))\<CeH\y"(Z)\,
for N/2 < i < N. Note that the above expression for the truncation error in the interval [a, 1] can also be represented as
|Lf (1-, - y(x,))| = j^Riixi, Xi-uy), (28) hi_i
Rn(a,p,g) = - / (p-Ong{n+1HOd(i
n- J a
denotes the remainder obtained from Taylor expansion in an integral form.
(xi - £)e-2 exp(-/3Ç/e)dÇ
< C [He + e-1 exp(-
< C [He + N.
Assuming that H < 2N 1 and e < H, we get ILN(Yi - y(xi))l < C (N-2 + N) .
Combining all the previous results, we obtain the required truncation error. Hence, we arrive at the desired result. ■
Theorem 3. Let y(x) be the solution of the IVP (18)-(19) and Yi be the numerical solution obtained from the hybrid scheme (20)-(21). Then, for sufficiently large N, and N_iao ln Np* < 1, where
3* = max 6(xi),
we have
\Yi - y(xi)\ < C [N_2 ln2 N + N_ie + N_l3a°] ,
Mxi G (32)
Proof. Let B_ = (2 - Pibi), B+ = (2 + pA) and b+ = (1 + pibi), where pi = hi/e.
The solution of the scheme (20)—(21) can be expressed as follows: For 0 < i < N/2,
"5=2 Bt
-(gi +92) H-----h +9i)
and for N/2 <i < N, 1
j=N/2 + i j
Yn/2 +
j=N/2 + ibj
9N/2+i
j=N/2+2 j
9N/2+2 + ••• +
Clearly, B+ 's and b+ 's are non-negative. For B- > 0, 0 < i < N/2, we have
2 - pibi = 2--^.
Since hi = 2N-1a0e ln N and bi < (3*, we have Bi- > 0. Consequently, the solution satisfies the discrete maximum principle and hence there are no oscillations.
Let us define the discrete barrier function:
4>i = C [N-2 ln2 N + N-1e + N-l3a°] .
Now, choosing C sufficiently large and using the discrete maximum principle, it is easier to see that
LN (h ± (Yi - y(xi))) > 0
or, equivalently,
LN (&) > \Yi - y(xi)\. Therefore, it follows that
\Yi ~ y{xi)\ < \4>i\, Vxj G f2. Thus, we have the required e-uniform error bound. ■
Remark 1. In Theorem 2, one can notice that the truncation error is of order Nfor H > e. It is assumed that pao > 2 and we are interested in the case of e < N_i. Also, we obtain the error bound of order N_ie only in the interval [a, 1] for the case H < e, which is not the practical case. With these points, we conclude that the order of convergence is almost 2 (up to a logarithmic factor).
5. Mesh and the scheme
A fitted mesh method for the problem (l)-(3) is now introduced. On Q, a piecewise uniform mesh of N mesh intervals is constructed as follows. The domain Cl is subdivided into the three subintervals as
Q = [0, a] U (a, 1 - a] U (1 - a, 1]
for some a that satisfies 0 < a < 1/4. On [0,a] and [1 - a, 1], a uniform mesh with N/4 mesh-intervals is placed, while [a, 1 - a] has a uniform mesh with N/ 2 mesh intervals. It is obvious that mesh is uniform when
a = 1/4. It is fitted to the problem by choosing a to be the function of N and e and
a = min {1/4, erov^ln N} ,
where a0 > 2/v/T?- Then, the hybrid scheme (20)-(21) for (10)—(11)becomes
£D_VLhi +
+ V/°11 ,i + Ol2,iVLl,i) = 0 = < for 0 <i < N/4 and 3N/4 <i < N,
eD~VL^i + y/au,i + «12, iVLi,i = 0
for N/4 < i < 3N/4,
VLi,o = 1. (34)
Similarly, we can define the hybrid scheme for (12)—(13), (14)—(15) and (16)—(17).
5.1. Description of the method. In this subsection, we describe the MIVT to solve (1)—(3):
Step 1. Solve the IVP (10)—(11) by using the hybrid scheme described on the Shishkin mesh. Let VLiii be its solution.
Step 2. Solve the IVP (12)—(13) by using the hybrid scheme. Let VL2ii be its solution.
Step 3. Solve the TVP (14)—(15) by using the hybrid scheme. Let WRiii be its solution.
Step 4. Solve the TVP (16)—(17) by using the hybrid scheme. Let WR2}i be its solution.
Step 5. Define mesh function Ui as
U ■ = (Un U i = \U2,,
= uRi,i
\UR2,i (
[p - uRi (0)]
[r - u r2 (0)]
[q - Uri (1)] [s - Ur2 (1)]
fflll(0) + ffll2(0) aii(xi) + ai2 (xi) _ a2i(0) + a22(0)
a2i (xi) + a22(xi )_ aii(1) + ai2(1)
aii(xi ) + ai2 (xi) a2i(1) + a22(1)
a2i (xí ) + a22 (xí)_
Vli, Vl2,
Wr2, (35)
Theorem 4. Let u(x) be the solution of the BVP (1)-(3) and Ui be the numerical solution obtained by the MIVT. Then we have
Ui - u(xi)
N~2 In2 N + N^e + 0 + sft
Proof. Theorem 3, when applied to the IVPs (10)-(11), (12)-(13) and the TVPs (14)-(15), (16)-(17), yields
\Vl1,i - vL1(xi)\
< C [/V-2 In2 N + N~1e + 0
for 0 < xi < 1,
\VL2,i - VL2(xi)\
NIn2 N + N^e + 0
for 0 < xi < 1,
Wrií - Wrl(Xi)
N-2 In2 N + N^s + N-^170
for 0 < xi < 1,
\Wr2í - WR2(xi)
N-2 In2 N + N^s +
for 0 < xi < 1.
From the definitions of uas(x), Ui and the above inequalities, we have
uas(xi) Ui
N-2 In2 N + N^s + 0
for xi e . From Theorem 1, we have
II u{xi) - uas{xi) II< cy/ë, x £ q.
, (36)
The desired estimate follows from the inequalities (36) and (37). ■
The exact solution of this example is not available. Therefore, to obtain the maximum pointwise errors and rates of convergence, we use the double mesh principle. By following the idea of Sun and Stynes (1995), we modify the Shishkin mesh. We calculate the numerical solution UN on Qe and the numerical solution UNon the mesh QN, where the transition parameter a is altered slightly to a = m\n{1/4,a0 e ln(N/2)}. Note that this slightly altered value of a will ensure that the positions
of transition points remain the same in meshes Qe and Q2eN. Hence, the use of interpolation for the double mesh principle can be avoided. The double mesh difference is defined as
E? = mzxJP? -U?«\},
where UN and IJ2N respectively denote the numerical solutions obtained by using N and 2N mesh intervals. The rates of convergence are calculated as
ln£f -InK
Tables 1 and 2 display respectively the maximum pointwise errors for u1 and u2 for several values of e and N taking ao = 2. 4
Example 2. Consider the following problem:
-eu'{(x) + 2(x + 1)2u1(x) - (x3 + 1)u2(x) = 2ex, -eu'2(x) - 2cos(nx/4)u1(x) + 2.2e-x+1u2(x)
x e (0,1], u(0) =
10x + 1, 0
1(1) =
Maximum pointwise errors and rate of convergence for u1 and u2 are given in Tables 3 and 4, respectively. From the rates of convergence one can conclude that the present method has second-order convergence up to a logarithmic factor. ♦
6. Numerical experiments and discussions
To show the applicability and efficiency of the present technique, two examples are provided. The computational results are given in the form of tables. The results are presented with the maximum point-wise errors for various values of e and N. We have also computed the computational order of convergence, which is shown in the same table along with the maximum errors.
Example 1. Consider the following problem:
-eu"(x) + 3u1(x) - u2(x) = 2, -eu'^x) - u1(x) + 3u2(x) = 3, x e (0, 1],
7. Conclusions
In this article, a robust computational technique is proposed for solving the system of two singularly perturbed reaction-diffusion problems. It is observed that, although the backward difference operator satisfies the discrete maximum principle in the whole domain [0,1], the its order is 1 (up to a logarithmic factor). We can get the order 2 (up to a logarithmic factor) by applying the trapezoidal scheme in [0,1], but it results in small oscillations, hence the solution is not stable unless the mesh size is very small even in the outer region [a, 1], where a coarse mesh is enough to give satisfactory results.
In order to retain the second-order convergence of the implicit trapezoidal scheme together with the
Table 1. Maximum pointwise errors and rates of convergence of m for Example 1.
e Number of mesh points
16 32 64 128 256 512 1024
1CT4 2.94547E-01 1.358 1.14899E-01 1.689 3.56373E-02 1.593 1.18109E-02 1.625 3.82798E-03 1.669 1.20389E-03 1.698 3.71125E-04
icrb 2.94547E-01 1.358 1.14899E-01 1.689 3.56373E-02 1.593 1.18109E-02 1.625 3.82798E-03 1.669 1.20389E-03 1.698 3.71125E-04
1CT8 2.94547E-01 1.358 1.14899E-01 1.689 3.56373E-02 1.593 1.18109E-02 1.625 3.82798E-03 1.669 1.20389E-03 1.698 3.71125E-04
1CT1U 2.94547E-01 1.358 1.14899E-01 1.689 3.56373E-02 1.593 1.18109E-02 1.625 3.82798E-03 1.669 1.20389E-03 1.698 3.71125E-04
icráS 2.94547E-01 1.358 1.14899E-01 1.689 3.56373E-02 1.593 1.18109E-02 1.625 3.82798E-03 1.669 1.20389E-03 1.698 3.71125E-04
10-4U 2.94547E-01 1.358 1.14899E-01 1.689 3.56373E-02 1.593 1.18109E-02 1.625 3.82798E-03 1.669 1.20389E-03 1.698 3.71125E-04
Table 2. Maximum pointwise errors and rates of convergence of u2 for Example 1.
e Number of mesh points
16 32 64 128 256 512 1024
1er4 2.94547E-01 1.358 1.14899E-01 1.689 3.56373E-02 1.593 1.18109E-02 1.625 3.82798E-03 1.669 1.20389E-03 1.698 3.71125E-04
icrb 2.94547E-01 1.358 1.14899E-01 1.689 3.56373E-02 1.593 1.18109E-02 1.625 3.82798E-03 1.669 1.20389E-03 1.698 3.71125E-04
10-« 2.94547E-01 1.358 1.14899E-01 1.689 3.56373E-02 1.593 1.18109E-02 1.625 3.82798E-03 1.669 1.20389E-03 1.698 3.71125E-04
10-1U 2.94547E-01 1.358 1.14899E-01 1.689 3.56373E-02 1.593 1.18109E-02 1.625 3.82798E-03 1.669 1.20389E-03 1.698 3.71125E-04
10-3s 2.94547E-01 1.358 1.14899E-01 1.689 3.56373E-02 1.593 1.18109E-02 1.625 3.82798E-03 1.669 1.20389E-03 1.698 3.71125E-04
10-4u 2.94547E-01 1.358 1.14899E-01 1.689 3.56373E-02 1.593 1.18109E-02 1.625 3.82798E-03 1.669 1.20389E-03 1.698 3.71125E-04
non-oscillating behavior of the backward difference operator, we proposed the hybrid scheme. This paper demonstrates the effectiveness of the Shishkin mesh by modifying the initial value technique (Valanarasu and Ramanujam, 2004) in a very simple way so that a higher order (nearly the second order) of convergence can be achieved with no restrictions on the values of h and e. The nonlinear system of equations has been handled by the present technique after linearization.
References
Bawa, R.K., Lal, A.K. and Kumar, V. (2011). An e-uniform hybrid scheme for singularly perturbed delay differential equations, Applied Mathematics and Computation 217(21): 8216-8222.
Das, P. and Natesan, S. (2013). A uniformly convergent hybrid scheme for singularly perturbed system of reaction-diffusion Robin type boundary-value problems,
Journal of Applied Mathematics and Computing 41(1): 447-471.
Doolan, E.P., Miller, J.J.H. and Schilders, W.H.A. (1980). Uniform Numerical Methods for Problems with Initial and Boundary Layers, Boole Press, Dublin.
Farrell, P.E., Hegarty, A.F., Miller, J.J.H., O'Riordan, E. and Shishkin, G.I. (2000). Robust Computational Techniques for Boundary Layers, Chapman & Hall/CRC Press, New York, NY.
Madden, N. and Stynes, M. (2003). A uniformly convergent numerical method for a coupled system of two singularly perturbed linear reaction-diffusion problems, IMAJournal of Numerical Analysis 23(4): 627-644.
Matthews, S., Miller, J.J.H., O'Riordan, E. and Shishkin, G.I. (2000). Parameter-robust numerical methods for a system of reaction-diffusion problems with boundary layers, in G.I. Shishkin, J.J.H. Miller and L. Vulkov (Eds.), Analytical and Numerical Methods for Convection-Dominated and Singularly Perturbed Problems, Nova Science Publishers, New York, NY, pp. 219-224.
Table 3. Maximum pointwise errors and rates of convergence of m for Example 2.
e Number of mesh points
16 32 64 128 256 512 1024
io-4 6.57432E-01 1.128 3.00718E-01 1.471 1.08497E-01 1.799 3.11750E-02 1.674 9.77281E-03 1.693 3.02225E-03 1.749 8.99222E-04
10~b 6.55359E-01 1.125 3.00433E-01 1.465 1.08804E-01 1.784 3.15870E-02 1.658 1.00113E-02 1.668 3.15045E-03 1.706 9.65973E-04
6.55157E-01 1.125 3.00405E-01 1.465 1.08835E-01 1.783 3.16281E-02 1.656 1.00351E-02 1.666 3.16325E-03 1.701 9.72645E-04
10-1U 6.55137E-01 1.125 3.00402E-01 1.465 1.08838E-01 1.783 3.16322E-02 1.656 1.00375E-02 1.666 3.16453E-03 1.701 9.73312E-04
10^ 6.55135E-01 1.125 3.00402E-01 1.465 1.08838E-01 1.783 3.16326E-02 1.656 1.00377E-02 1.666 3.16466E-03 1.701 9.73379E-04
10-I4 6.55135E-01 1.125 3.00402E-01 1.465 1.08838E-01 1.783 3.16326E-02 1.656 1.00377E-02 1.666 3.16466E-03 1.701 9.73379E-04
6.55135E-01 1.125 3.00402E-01 1.465 1.08838E-01 1.783 3.16326E-02 1.656 1.00377E-02 1.666 3.16466E-03 1.701 9.73379E-04
10-4U 6.55135E-01 1.125 3.00402E-01 1.465 1.08838E-01 1.783 3.16326E-02 1.656 1.00377E-02 1.666 3.16466E-03 1.701 9.73379E-04
Table 4. Maximum pointwise errors and rates of convergence of u2 for Example 2.
e Number of mesh points
16 32 64 128 256 512 1024
io-4 4.60874E-01 1.368 1.78595E-01 1.631 5.76521E-02 1.547 1.97354E-02 1.530 6.83615E-03 1.469 2.46891E-03 1.370 9.55156E-04
10~b 4.61430E-01 1.393 1.75721E-01 1.692 5.43734E-02 1.602 1.79112E-02 1.613 5.85455E-03 1.647 1.86935E-03 1.659 5.91967E-04
10-« 4.61487E-01 1.395 1.75436E-01 1.699 5.40469E-02 1.608 1.77292E-02 1.623 5.75642E-03 1.667 1.81307E-03 1.693 5.60807E-04
10-1U 4.61493E-01 1.396 1.75408E-01 1.699 5.40143E-02 1.609 1.77110E-02 1.624 5.74661E-03 1.668 1.80791E-03 1.696 5.57848E-04
10^ 4.61494E-01 1.396 1.75405E-01 1.699 5.40110E-02 1.609 1.77092E-02 1.624 5.74563E-03 1.669 1.80739E-03 1.697 5.57552E-04
io-14 4.61494E-01 1.396 1.75405E-01 1.699 5.40110E-02 1.609 1.77092E-02 1.624 5.74563E-03 1.669 1.80739E-03 1.697 5.57552E-04
10-Ib 4.61494E-01 1.396 1.75405E-01 1.699 5.40110E-02 1.609 1.77092E-02 1.624 5.74563E-03 1.669 1.80739E-03 1.697 5.57552E-04
4.61494E-01 1.396 1.75405E-01 1.699 5.40110E-02 1.609 1.77092E-02 1.624 5.74563E-03 1.669 1.80739E-03 1.697 5.57552E-04
10-4U 4.61494E-01 1.396 1.75405E-01 1.699 5.40110E-02 1.609 1.77092E-02 1.624 5.74563E-03 1.669 1.80739E-03 1.697 5.57552E-04
World Scientific, Singapore.
Natesan, S. and Briti, S.D. (2007). A robust computational method for singularly perturbed coupled system of reaction-diffusion boundary value problems, Applied Mathematics and Computation 188(1): 353-364.
Nayfeh, A.H. (1981). Introduction to Perturbation Methods, Wiley, New York, NY.
Rao, S.C.S., Kumar, S. and Kumar, M. (2011). Uniform global convergence of a hybrid scheme for singularly perturbed
Matthews, S., O'Riordan, E. and Shishkin, G.I. (2002). A numerical method for a system of singularly perturbed reaction-diffusion equations, Journal of Computational and Applied Mathematics 145(1): 151-166.
Melenk, J.M., Xenophontos, C. and Oberbroeckling, L. (2013). Analytic regularity for a singularly perturbed system of reaction-diffusion equations with multiple scales, Advances in Computational Mathematics 39(2): 367-394.
Miller, J.J.H., O'Riordan, E. and Shishkin, G.I. (1996). Fitted Numerical Methods for Singular Perturbation Problems,
reaction-diffusion systems, Journal of Optimization Theory and Applications 151(2): 338-352.
Roos, H.-G., Stynes, M. and Tobiska, L. (1996). Numerical Methods for Singularly Perturbed Differential Equations, Springer, Berlin.
Shishkin, G.I. (1995). Mesh approximation of singularly perturbed boundary-value problems for systems of elliptic and parabolic equations, Computational Mathematics and Mathematical Physics 35(4): 429-446.
Sun, G. and Stynes, M. (1995). An almost fourth order uniformly convergent difference scheme for a semilinear singularly perturbed reaction-diffusion problem, Numerische Mathematik 70(4): 487-500.
Valanarasu, T. and Ramanujam, N. (2004). An asymptotic initial-value method for boundary value problems for a system of singularly perturbed second-order ordinary differential equations, Applied Mathematics and Computation 147(1): 227-240.
Vinod Kumar received his Ph.D. (mathematics) from Thapar University, Patiala, India, in 2013. He worked at SLIET, Longowal (Deemed University) and Chitkara University. His present area of interest includes parallel and scientific computations.
Rajesh K. Bawa obtained his Ph.D. in numerical computing from IIT Kanpur in 1994. Presently, he is a professor and the head of the Department of Computer Science, Punajbi University, Patiala. Before, he also served at SLIET, Longowal (Deemed University) and Thapar University at Patiala. His present area of interest is parallel and scientific computation. He has published numerous research papers in learned journals of reputed publishers such as Elsevier, Springer, Taylor and Francis, etc.
Arvind K. Lal received his M.Sc. (mathematics) in 1987 from the University of Bihar, Muzaffarpur, India, and a Ph.D. (mathematics) from the University of Roorkee, Roorkee, India (presently IIT, Roorkee) in
1995. He has been working in Thapar University, Patiala, India, since
1996. Currently, he is an associate professor at that university. His research interests include numerical analysis, theoretical astrophysics and reliability analysis
Received: 7 March 2013 Revised: 29 November 2013
