A UNIFORMLY CONVERGENT NUMERICAL METHOD FOR A WEAKLY COUPLED SYSTEM OF SINGULARLY PERTURBED CONVECTION-DIFFUSION PROBLEMS WITH BOUNDARY AND WEAK INTERIOR LAYERS†
A UNIFORMLY CONVERGENT NUMERICAL METHOD FOR A WEAKLY COUPLED SYSTEM OF SINGULARLY PERTURBED CONVECTION-DIFFUSION PROBLEMS WITH BOUNDARY AND WEAK INTERIOR LAYERS†
Journal of Applied Mathematics & Informatics. 2015. Sep, 33(5_6): 635-648
• Received : February 05, 2015
• Accepted : June 01, 2015
• Published : September 30, 2015
PDF
e-PUB
PPT
Export by style
Share
Article
Author
Metrics
Cited by
TagCloud
SHEETAL CHAWLA
S. CHANDRA SEKHARA RAO

Abstract
We consider a weakly coupled system of singularly perturbed convection-diffusion equations with discontinuous source term. The diffusion term of each equation is associated with a small positive parameter of different magnitude. Presence of discontinuity and different parameters creates boundary and weak interior layers that overlap and interact. A numerical method is constructed for this problem which involves an appropriate piecewise uniform Shishkin mesh. The numerical approximations are proved to converge to the continuous solutions uniformly with respect to the singular perturbation parameters. Numerical results are presented which illustrates the theoretical results. AMS Mathematics Subject Classification : 65L10, 65L11, 65L12.
Keywords
1. Introduction
An extensive research had been done on numerical methods for a single singularly perturbed convection-diffusion differential equation [1] - [4] , but for system of equations very few works had been done. The classical numerical methods fail to produce good approximations for singularly perturbed problems. Various non-classical approaches produce better approximations and converge uniformly with respect to the small perturbation parameter. In the literature [7] - [14] methods were available to obtain numerical approximation for system of singularly perturbed convection- diffusion differential equations the source term are smooth on the whole domain. Farrell et.al [5] - [6] considered scalar singularly perturbed convection-diffusion equation with discontinuous source term. The interior layers in [6] were strong, in the sense that the solution was bounded but the magnitude of the first derivative grew unboundedly as ε → 0, but in [5] they were weak, in the sense that the solution and the first derivative were bounded but the magnitude of the second derivative grows unboundedly as ε → 0. In this work, we present a uniformly convergent numerical method for a weakly coupled system of singularly perturbed convection-diffusion equations having discontinuous source term with different diffusion parameters. The solution to such equations has overlapping and interacting boundary and interior layers which makes the construction of numerical methods and analysis quite difficult. Tamilselvan and Ramanujam [15] considered the same problem but with equal diffusion parameters.
Consider a weakly coupled system of singularly perturbed convection-diffusion equations with discontinuous source term on the unit interval Ω = (0, 1), having a single discontinuity in the source term at a point d ∈ Ω. Let Ω 1 = (0, d ) and Ω 2 = ( d , 1). Let the jump in a function ω at a point d ∈ Ω given as [ ω ]( d ) = ω ( d +) − ω ( d −). The corresponding boundary value problem is to find u 1 , u 2 C 0 (
PPT Slide
Lager Image
) ∩ C 1 (Ω) ∩ C 2 1 ∪ Ω 2 ), such that
PPT Slide
Lager Image
PPT Slide
Lager Image
where E = diag( ε 1 , ε 2 ), the coupling matrix A = diag( a 1 , a 2 ) and B = ( bij ) 2×2 with 0 < ε 1 ε 2 ≤ 1, f = ( f 1 , f 2 ) T , and u = ( u 1 , u 2 ) T . Assume for each i = 1, 2 and x
PPT Slide
Lager Image
, the matrices A and B satisfy
PPT Slide
Lager Image
PPT Slide
Lager Image
Let α = min{ α 1 , α 2 }. Further assume that the source terms f 1 , f 2 are sufficiently smooth on
PPT Slide
Lager Image
and their derivatives have jump discontinuity at the same point.
Notations. Throughout the paper, C denotes a generic positive constant and C = ( C , C ) T denotes a generic positive constant vector, both are independent of perturbation parameters ε 1 , ε 2 and the discretization parameter N, but may not be same at each occurrence. Define v w if vi wi , i = 1, 2, and | v | = (| v 1 |, | v 2 |) T . We consider the maximum norm and denote it by║.║ S , where S is a closed and bounded subset in
PPT Slide
Lager Image
. For a real valued function v C ( S ) and for a vector valued function v = ( v 1 , v 2 ) T C ( S ) 2 , we define
PPT Slide
Lager Image
and ║ v S = max{║ v 1 S ,║ v 2 S }. Now let a mesh
PPT Slide
Lager Image
be a set of points satisfying x 0 < x 1 <⋯< xN = 1. A mesh function
PPT Slide
Lager Image
is a real-valued function defined on Ω N . Define the discrete maximum norm for such functions by
PPT Slide
Lager Image
and for vector mesh functions V = ( V 1 , V 2 ) T =
PPT Slide
Lager Image
are used and define║ V ΩN = max{║ V 1 ΩN ,║ V 2 ΩN }.
- 2. The continuous problem
Theorem 2.1 (Continuous maximum principle). Suppose u 1 , u 2 C 0 (
PPT Slide
Lager Image
) ∩ C 1 (Ω) ∩ C 2 1 ∪ Ω 2 ). Further suppose that u = ( u 1 , u 2 ) T satisfies u (0) ≥ 0, u (1) ≥ 0, Lu ( x ) ≥ 0 in Ω 1 ∪ Ω 2 and [ u′ ]( d ) ≤ 0. Then u ( x ) ≥ 0, for all x
PPT Slide
Lager Image
.
Proof . Let u :=
PPT Slide
Lager Image
Now let p and q be the points at which θ 1 ( p ) :=
PPT Slide
Lager Image
and θ 2 ( q ) :=
PPT Slide
Lager Image
. Assume without loss of generality θ 1 ( p ) ≤ θ 2 ( q ). If θ 1 ( p ) ≥ 0, then there is nothing to prove. Suppose that θ 1 ( p ) < 0, then proof is completed by showing that this leads to a contradiction. Note that p ≠ {0, 1}. So either p ∈ Ω 1 ∪ Ω 2 or p = d .
In the first case for x ∈ Ω 2 ,
PPT Slide
Lager Image
In the second case, that is, p = d , we have [ u′ ]( d ) =
PPT Slide
Lager Image
, and at a negative minimum [Θ′]( d ) ≥ 0, which gives a contradiction. □
Lemma 2.2 (Stability Result). Let u = ( u 1 , u 2 ) T be the solution of (1) − (2), then,
PPT Slide
Lager Image
Proof . Define the function Ψ ± ( x ) := max {║ u (0)║,║ u (1)║,║ Lu ║Ω 1 ∪Ω 2 } (2 − x , 2 − x ) T ± u ( x ). Then Ψ ± (0) ≥ 0 , Ψ ± (1) ≥ 0 , L Ψ ± ( x ) ≥ 0 for each x ∈ Ω 1 ∪ Ω 2 , and [Ψ ± ′ ]( d ) = ±[ u′ ]( d ) = 0 , since u C 1 (Ω) 2 . It follows from the maximum principle that Ψ ± ( x ) ≥ 0 for all x ∈ Ω, which leads to the required bound on u . Consequently, the problem (1) − (2) has a unique and stable solution. □
To derive sharper bounds on the derivatives of solution, the solution is decomposed into a sum, composed of a regular component v and a singular component w . That is, u = v + w . The regular component v , can be written in the form v =
PPT Slide
Lager Image
, where v 0 = ( v 01 , v 02 ) T , v 1 = ( v 11 , v 12 ) T and v 2 = ( v 21 , v 22 ) T are defined respectively to be the solutions of the problems
PPT Slide
Lager Image
and
PPT Slide
Lager Image
Thus the regular component v is the solution of
PPT Slide
Lager Image
Further we decompose w as w = w 1 + w 2 where w 1 = ( w 11 , w 12 ) T , w 2 = ( w 21 , w 22 ) T . Thus w 1 = w 11 + w 21 and w 2 = w 12 + w 22 , where w 1 is the solution of
PPT Slide
Lager Image
and w 2 is the solution of
PPT Slide
Lager Image
Lemma 2.3. For each integer k, satisfying 0 ≤ k ≤ 3 , the regular component v and its derivatives satisfy the bounds given by
PPT Slide
Lager Image
Proof . The proof follows from [7] and [2] . □
Lemma 2.4. For each integer k, satisfying 0 ≤ k ≤ 3 , the singular component w 1 and its derivatives satisfy the bounds given by
PPT Slide
Lager Image
Proof . The proof follows from [7] and [2] . □
Lemma 2.5. For each integer k, satisfying 0 ≤ k ≤ 3 , the singular component w2 and its derivatives satisfy the bounds given by
PPT Slide
Lager Image
Proof . Consider the barrier function ϕ (1, 1) T ± w 2 , where
PPT Slide
Lager Image
to bound w 2 . To bound derivatives of w 2 , use the technique used in [7] and bound on w 2 on the domain Ω 1 and Ω 2 . □
3. Discretization of the Problem
We use piecewise uniform Shishkin mesh which uses these transition parameters:
PPT Slide
Lager Image
The interior points of the mesh are denoted by
PPT Slide
Lager Image
Let hi = xi x i −1 be the ith mesh step and
PPT Slide
Lager Image
, clearly
PPT Slide
Lager Image
= { xi : i = 0, 1,..., N }. Let N = 2 l , l ≥ 5 be any positive integer.
We divide
PPT Slide
Lager Image
into three sub-intervals [0, σєl1 ], [ σєl1 , σєl2 ] and [ σєl2 , d ] for some 0 < σєl1 σєl2
PPT Slide
Lager Image
. The sub-intervals [0, σєl1 ] and [ σєl1 , σєl2 ] are divided into N /8 equidistant elements and the sub-interval [ σєl2 , d ] is divided into N /4 equidistant elements. Similarly, in
PPT Slide
Lager Image
the sub-intervals [ d , d + σєr1 ] and [ d + σєr1 , d + σєr2 ] are divided into N /8 equidistant elements and the sub-interval [ d + σєr2 , 1] is divided into N /4 equidistant elements, for some 0 < σєr1 σєr2
PPT Slide
Lager Image
.
Define the discrete finite difference operator LN as follows
PPT Slide
Lager Image
with boundary conditions
PPT Slide
Lager Image
where
PPT Slide
Lager Image
and at x N /2 = d the scheme is given by
PPT Slide
Lager Image
where
PPT Slide
Lager Image
Lemma 3.1. Suppose that a mesh function Z ( xi ) satisfies Z ( x 0 ) ≥ 0,Z ( xN ) ≥ 0, LNZ ( xi ) ≥ 0 for all xi ∈ Ω N and
PPT Slide
Lager Image
0, then Z ( xi ) ≥ 0 for all xi
PPT Slide
Lager Image
.
Lemma 3.2. If Z ( xi ) is any mesh function, then,
PPT Slide
Lager Image
The discrete solution U can be decomposed into the sum U = V + W . The function V , is defined as the solution of the following problem:
PPT Slide
Lager Image
PPT Slide
Lager Image
The function W , is defined as the solution of the following problem:
PPT Slide
Lager Image
PPT Slide
Lager Image
where the jump in the discrete derivative of a mesh function Z at the point xi = d is given by:
PPT Slide
Lager Image
Further decompose W as W = W 1 + W 2 , where the function W 1 is defined as the solution of the following problem:
PPT Slide
Lager Image
PPT Slide
Lager Image
and the function W 2 is defined as the solution of the following problem:
PPT Slide
Lager Image
PPT Slide
Lager Image
4. Convergence analysis
By Taylor’s expansion and bounds on regular components defined in lemma 2.3 gives
PPT Slide
Lager Image
Define the mesh function Ψ ± ( xi ) as
PPT Slide
Lager Image
Using discrete maximum principle, the error of the regular component satisfies the estimate
PPT Slide
Lager Image
As in [7] , the error of the singular component satisfies the estimate
PPT Slide
Lager Image
Lemma 4.1. The following ε 1 , ε 2 uniform bound
PPT Slide
Lager Image
where W 2 is the solution of (14) − (15).
Proof . At the point x = d we know that
PPT Slide
Lager Image
First consider
PPT Slide
Lager Image
From lemma 2.3 we have
PPT Slide
Lager Image
Therefore, | D V ( d )| ≤
PPT Slide
Lager Image
Similarly, consider
PPT Slide
Lager Image
Again from lemma 2.3 we have
PPT Slide
Lager Image
Therefore, | D + V ( d )| ≤
PPT Slide
Lager Image
On Ω 1 , | W 1 ( xi )| ≤
PPT Slide
Lager Image
implies that | D W 1 ( d )| ≤
PPT Slide
Lager Image
. On Ω 2 , D + W 1 ( d ) = D + ( W 1 w 1 )( d ) + D + w 1 ( d ). From lemma 2.4
PPT Slide
Lager Image
and | D + ( W 1 w 1 )( xi )| ≤
PPT Slide
Lager Image
.
Therefore,
PPT Slide
Lager Image
Lemma 4.2. The following ε 1 , ε 2 uniform bound
PPT Slide
Lager Image
is valid, where W 2 is the solution of (14) − (15).
Proof . Consider the following function
PPT Slide
Lager Image
= 1, 2 where
PPT Slide
Lager Image
where Ψ = ( ψ 1 , ψ 2 ) T is the solution of
PPT Slide
Lager Image
Using the discrete maximum principle we get the required result.
Lemma 4.3. The error of the singular component satisfies the estimate
PPT Slide
Lager Image
Proof .
PPT Slide
Lager Image
Now
PPT Slide
Lager Image
From lemma 4.1 we have
PPT Slide
Lager Image
and
PPT Slide
Lager Image
Hence
PPT Slide
Lager Image
Likewise,
PPT Slide
Lager Image
Using the bounds on derivatives of w 2 given in lemma 2.5, we have
PPT Slide
Lager Image
Using the bounds on derivatives of w 1 given in lemma 2.4, we have
PPT Slide
Lager Image
Also,
PPT Slide
Lager Image
Collecting all the previous inequalities we get that
PPT Slide
Lager Image
By the Taylor’s expansion and bounds on the derivatives of w 2 , given in lemma 2.5 we have
Case (i) For xi ∈ [ d + σєr2 , 1].
PPT Slide
Lager Image
Similar arguments prove a similar result for the sub-interval [ σєl1 , d ).
Case (ii) For xi ∈ (0, σєl1 ).
PPT Slide
Lager Image
Similar arguments prove a similar result for the sub-interval ( d , d + σєr1 ).
Case (iii) [ σєl1 , σєl2 ).
PPT Slide
Lager Image
Using the inequality
PPT Slide
Lager Image
Likewise, |( LN ( W 2 w 2 )) 2 ( xi )| ≤ CN −1 ln N .
Similar arguments prove a similar result for the sub-interval [ d + σєr1 , d+ σєr2 ).
Combining all these gives,
PPT Slide
Lager Image
Consider the following function
PPT Slide
Lager Image
, j = 1, 2 where
PPT Slide
Lager Image
where Ψ = ( ψ 1 , ψ 2 ) T is the solution of
PPT Slide
Lager Image
Using the discrete maximum principle we get the required result. □
Theorem 4.4. Let u be the solution of given problem (1) − (2) and U is the solution of discrete problem on the Shishkin mesh defined in section 3, then
PPT Slide
Lager Image
Proof . Using the equation (16), (17) and lemma 4.3 we get the required result. □
5. Numerical Results
To illustrate the theoretical results the scheme in Section 3 is implemented on these test examples.
Example 5.1 Consider the following singularly perturbed convection-diffusion problem with discontinuous source term:
PPT Slide
Lager Image
where
PPT Slide
Lager Image
For the construction of piecewise-uniform Shishkin mesh
PPT Slide
Lager Image
, we take α = 0.8. The Exact solution of the examples are not known. Therefore we estimate the error for U by comparing it to the numerical solution
PPT Slide
Lager Image
obtained on the mesh
PPT Slide
Lager Image
that contains the mesh points of the original mesh and their midpoints, that is,
PPT Slide
Lager Image
= xj , j=0,...,N,
PPT Slide
Lager Image
= ( xj + x j +1 )/2, j=0,...,N-1.
For different values of N and ε 1 , ε 2 , we compute
PPT Slide
Lager Image
Maximum point-wise errorsandε1,ε2−uniform rate of convergencepNfor Example 5.1.
PPT Slide
Lager Image
Maximum point-wise errors and ε1, ε2−uniform rate of convergence pN for Example 5.1.
If ε 1 = 10 j for some non-negative integer j , set
PPT Slide
Lager Image
Then the parameter-uniform error is computed as DN :=
PPT Slide
Lager Image
and the order of convergence is calculated using the formula pN :=
PPT Slide
Lager Image
Finally, we want to show that similar results can be obtained for coupled system of M (> 2) singularly perturbed convection diffusion problem with discontinuous source term. Letting N = 2 M × τ , where τ is some positive power of 2, the mesh is defined using the following transition points
PPT Slide
Lager Image
Then we divide the interval [0, d ] into M+1 subintervals [0, σєl1 ], [ σєl1 , σєl2 ],..., [ σєlM−1 , σєlM ], [ σєlM , d ]. On the subinterval [0, σєl1 ] a uniform mesh of N /2 M+1 mesh intervals, on [ σєlk , σєlk+1 ], 1 ≤ k M − 1, a uniform mesh of N /2 M−k+2 mesh intervals, and on [ σєlM , d ] a uniform mesh of N/4 mesh intervals are placed. Similarly, we divide the interval [ d , 1] into subintervals [ d , d + σєr1 ], [ d + σєr1 , d + σєr2 ],..., [ d + σєrM−1 , d + σєrM ], [ d + σєrM , 1]. On the subinterval [ d , d + σєr1 ] a uniform mesh of N /2 M+1 mesh intervals, on [ d + σєrk , d + σєrk +1], 1 ≤ k ≤ M − 1, a uniform mesh of N /2 M−k+2 mesh intervals, and on [ d + σєrM , 1] a uniform mesh of N /4 mesh intervals are placed. Let h єl1 and h єr1 be the mesh lengths on [0, σєl1 ] and on [ d , d + σєr1 ] respectively. Let H 1 and H 2 be the mesh lengths on [ σєlM , d ] and on [ d + σєrM , 1] respectively; hєlk and hєrk be the mesh lengths on [ σєlk , σєlk+1 ] and on [ d + σєrk , d + σєrk+1 ], k = 2,..., M respectively.
In this case also, we obtain the scheme similar to (3.1), with u = ( u 1 ; u 2 ,..., uM ) T
PPT Slide
Lager Image
C 1 (Ω) M C 2 1 ∪ Ω 2 ) M and also expect the error bound
PPT Slide
Lager Image
C ( N −1 ln N ) to hold, although attempts of the authors have failed so far to provide a proof. To illustrate the order of uniform convergence of this method we consider the following test example.
Example 5.2 Consider the following singularly perturbed convection-diffusion problem with discontinuous source term:
PPT Slide
Lager Image
where f 1 ( x ) =
PPT Slide
Lager Image
and
PPT Slide
Lager Image
Maximum point-wise errors, withε2= 10−4,ε3= 10−1andε1,ε2,ε3−uniform rate of convergencepNfor Example 5.2.
PPT Slide
Lager Image
Maximum point-wise errors , with ε2 = 10−4, ε3 = 10−1 and ε1, ε2, ε3−uniform rate of convergence pNfor Example 5.2.
BIO
Sheetal Chawla is doing her Ph.D. at Indian Institute of Technology, Delhi under the direction of Dr. S. Chandra Sekhara Rao. Her area of interests are Numerical Analysis and Differential Equations.
Department of Mathematics, Indian Institute of Technology Delhi, New Delhi-110016.
e-mail: chawlaasheetal@gmail.com
S. Chandra Sekhara Rao is working as a Professor, Department of Mathematics, Indian Institute of Technology, Delhi. His area of interests are Numerical Analysis and Parallel Computing.
Department of Mathematics, Indian Institute of Technology Delhi, New Delhi-110016.
e-mail: scsr@maths.iitd.ac.in
References
Miller J.J.H. , Riordan E.O. , Shishkin G.I. 1996 Fitted Numerical Methods for Singular Perturbation Problems: Error Estimates in the Maximum Norm for Linear Problems in One and Two Dimensions World Scientific Singapore
Farrell P.A. , Hegarty A.F. , Miller J.J.H. , Riordan E.O. , Shishkin G.I. 2000 Robust Computational Techniques for Boundary Layers Chapman & Hall, CRC Press Boca Raton, Florida
Roos H.G. , Stynes M. , Tobiska L. 1996 Numerical Methods for Singularly Perturbed Differential Equations Springer Berlin
Roos H.G. , Stynes M. , Tobiska L. 2008 Robust Numerical Methods for Singularly Perturbed Differential Equations Springer Berlin
Farrell P.A. , Hegarty A.F. , Miller J.J.H. , Riordan E.O , Shishkin G.I. (2004) Singularly perturbed convection-diffusion problems with boundary and weak interior layers J. Comput. Appl. Math. 166 133 - 151    DOI : 10.1016/j.cam.2003.09.033
Farrell P.A. , Hegarty A.F. , Miller J.J.H. , Riordan E.O. , Shishkin G.I. (2004) Global Maximum Norm Parameter-Uniform Numerical Method for a Singularly Perturbed Convection-Diffusion Problem with Discontinuous Convection Coefficient Math. Comput. Modelling 40 1375 - 1392    DOI : 10.1016/j.mcm.2005.01.025
Cen Z. (2005) Parameter-uniform finite difference scheme for a system of coupled singularly perturbed convection-diffusion equations Intern. J. Comput. Math. 82 177 - 192    DOI : 10.1080/0020716042000301798
Linss T. (2009) Analysis of a system of singularly perturbed convection-diffusion equations with strong coupling SIAM J. Numer. Anal. 47 1847 - 1862    DOI : 10.1137/070683970
Clavero C. , Gracia J.L. , Lisbona F. (1999) High order methods on Shishkin meshes for singular perturbation problems of convection-diffusion type Numer. Algorithms 22 73 - 97    DOI : 10.1023/A:1019150606200
Riordan E.O. , Stynes M. (2009) Numerical analysis of a strongly coupled system of two singularly preturbed convection-diffusion problems Adv. Comput. Math. 30 101 - 121    DOI : 10.1007/s10444-007-9058-z
Bellew S. , Riordan E.O. (2004) A parameter robust numerical method for a system of two singularly perturbed convection-diffusion equations Appl. Numer. Math. 51 171 - 186    DOI : 10.1016/j.apnum.2004.05.006
Munyakazi J.B. (2014) A uniformly convergent nonstandard finite difference scheme for a system of convection-diffusion equations Comp. Appl. Math.    DOI : 10.1007/s40314-014-0171-6
Priyadharshini R.M. , Ramanaujam N. (2013) Uniformly-convergent numerical methods for a system of coupled singularly perturbed convection-diffusion equations with mixed type boundary conditions Math. Model. Anal. 18 577 - 598    DOI : 10.3846/13926292.2013.851629
Tamilselva A. , Ramanujam N. , Priyadharshini R.M. , Valanarasu T. (2010) Parameter - uniform numerical method for a system of coupled singularly perturbed convection-diffusion equations with mixed type boundary conditions J. Appl. Math. Informatics 28 109 - 130
Tamilselvan A. , Ramanujam N. (2009) A numerical method for singularly perturbed system of second order ordinary differential equations of convection diffusion type with a discontinuous source term J. Appl. Math. Informatics 27 1279 - 1292