Fitted mesh B-spline collocation method for solving self-adjoint singularly perturbed boundary value problems

Fitted mesh B-spline collocation method for solving self-adjoint singularly perturbed boundary value problems

Applied Mathematics and Computation 161 (2005) 973–987 www.elsevier.com/locate/amc Fitted mesh B-spline collocation method for solving self-adjoint s...

243KB Sizes 0 Downloads 74 Views

Applied Mathematics and Computation 161 (2005) 973–987 www.elsevier.com/locate/amc

Fitted mesh B-spline collocation method for solving self-adjoint singularly perturbed boundary value problems Mohan K. Kadalbajoo, Vivek K. Aggarwal

*

Department of Mathematics, IIT-Kanpur, 208016 Kanpur, India

Abstract In this paper we develop B-spline collocation method for solving a class of selfadjoint singularly perturbed equations given by ðaðxÞy 0 Þ0 þ bðxÞyðxÞ ¼ gðxÞ; aðxÞ P a > 0; yð0Þ ¼ a; yð1Þ ¼ b:

a0 ðxÞ P 0;

bðxÞ P b > 0; ð0:1Þ

We use the fitted mesh technique to generate piecewise uniform mesh, and use B-spline method which leads to a tridiagonal linear system. The convergence analysis is given and the method is shown to have uniform convergence of second order. Numerical illustrations are given in the end to demonstrate the efficiency of our method. Ó 2004 Elsevier Inc. All rights reserved. Keywords: Self-adjoint singularly perturbed boundary value problems; B-spline collocation method; Fitted mesh method; Boundary layers; Uniform convergence

1. Introduction We consider the following class of self-adjoint singularly perturbed twopoint boundary value problems 0

Ly  ðaðxÞy 0 Þ ðxÞ þ bðxÞyðxÞ ¼ gðxÞ where 0 6 x 6 1;

*

ð1:1Þ

Corresponding author. E-mail addresses: [email protected] (M.K. Kadalbajoo), [email protected] (V.K. Aggarwal).

0096-3003/$ - see front matter Ó 2004 Elsevier Inc. All rights reserved. doi:10.1016/j.amc.2003.12.078

974

M.K. Kadalbajoo, V.K. Aggarwal / Appl. Math. Comput. 161 (2005) 973–987

subject to yð0Þ ¼ a

yð1Þ ¼ b;

a; b 2 R;

ð1:2Þ

where  is a small positive parameter and aðxÞ, bðxÞ and gðxÞ are smooth functions and satisfy aðxÞ P a > 0;

a0 ðxÞ P 0;

bðxÞ P b > 0:

Under these conditions the operator L admits a maximum principle [1]. Earlier this type of problems have been solved by numerous researchers. Boglave [2] and Schatz and Wahlbin [3] used finite element techniques to solve such problems. Niijima [4,5] gave uniformly second-order accurate difference schemes whereas Miller [6] gave sufficient conditions for the uniform first-order convergence of a general three-point difference scheme. These class of problems arises in various fields of science and engineering, for instance, fluid mechanics, quantum mechanics, optimal control, chemical-reactor theory, aerodynamics, reaction-diffusion process, geophysics etc. Surla and Jerkovic [7] considered the singularly perturbed boundary value problem using spline collocation method. Sakai and Usmani [8] gave a new concept of B-spline in terms of hyperbolic and trigonometric splines which are different from earlier ones. It is proved that the hyperbolic and trigonometric B-spline are characterized by a convolution of some special exponential functions and a characteristic function on the interval ½0; 1 . In this paper we use B-Spline collocation method [9–11], to solve the problems of above type (1.1). We reduce the original problem (1.1) to the normal form. Then we solve the normalized equation by B-spline collocation method using fitted mesh technique. In Section 2 we have given the derivation for the B-spline method and mesh strategy. Uniform convergence for the method has been discussed in Section 3. In Sections 4 and 5 we have solved some problems using the method and the results and graphs have also been shown.

2. B-spline collocation method Rewrite (1.1) as y 00 þ P ðxÞy 0 þ QðxÞy ¼ RðxÞ;

ð2:1Þ

where P ðxÞ ¼ a0 ðxÞ=aðxÞ;

QðxÞ ¼ bðxÞ=aðxÞ

and

RðxÞ ¼ gðxÞ=aðxÞ:

Let yðxÞ ¼ U ðxÞV ðxÞ;

ð2:2Þ

M.K. Kadalbajoo, V.K. Aggarwal / Appl. Math. Comput. 161 (2005) 973–987

975

and transform Eq. (2.1) into the normal form, as V 00 ðxÞ þ AðxÞV ðxÞ ¼ BðxÞ;

ð2:3Þ

where 2

AðxÞ ¼ QðxÞ  1=2P 0 ðxÞ  1=4ðP ðxÞÞ ; Z x BðxÞ ¼ RðxÞ expð1=2 P ðtÞ dtÞ; 0 Z x U ðxÞ ¼ expð1=2 P ðtÞ dtÞ;

ð2:4Þ

0

with V ð0Þ ¼ Uyð0Þ ¼ a0 , V ð1Þ ¼ Uyð1Þ ¼ a1 , x 2 ð0; 1Þ and a0 ; a1 2 R. Multiplying ð0Þ ð1Þ (2.3) throughout by , we get L1 V ðxÞ ¼ V 00 ðxÞ þ W ðxÞV ðxÞ ¼ ZðxÞ

ð2:5Þ

with boundary conditions V ð0Þ ¼ a0 ;

V ð1Þ ¼ a1 ;

ð2:6Þ

where W ðxÞ ¼ AðxÞ;

ZðxÞ ¼ BðxÞ

and

W ðxÞ P W  > 0:

The approximate solution of the problem (2.5) and (2.6) is obtained by using B-spline collocation method as described below. We subdivide the interval ½0; 1 , and we choose piecewise uniform mesh points represented by p ¼ fx0 ; x1 ; x2 ; . . . ; xN g, such that x0 ¼ 0 and xN ¼ 1 and ~h (discussed in the next section) is the piecewise uniform spacing. We define L2 ½0; 1 as a vector space of all the square integrable functions on ½0; 1 , and X be the linear subspace of L2 ½0; 1 . Now we define for i ¼ 0; 1; 2; . . . ; N . 8 > ðx  xi2 Þ3 ; > > > > if x 2 ½xi2 ; xi1 ; > > > 2 3 3 ~ > þ 3~ h2 ðx  xi1 Þ þ 3~ hðx  xi1 Þ  3ðx  xi1 Þ ; h > > > > 1 < 3 if x 22 ½xi1 ; xi ; 2 3 ~ ð2:7Þ Bi ðxÞ ¼ h ðxiþ1  xÞ þ 3~ hðxiþ1  xÞ  3ðxiþ1  xÞ ; h þ 3~ ~ h3 > > > if x 2 ½xi ; xiþ1 ; > > > 3 > ðxiþ2  xÞ ; > > > > > if x 2 ½xiþ1 ; xiþ2 ; > : 0; otherwise: We introduce four additional knots as x2 < x1 < x0 and xN þ2 > xN þ1 > xN . From the above Eq. (2.7) we can simply check that each of the functions Bi ðxÞ is twice continuously differentiable on the entire real line, Also

976

M.K. Kadalbajoo, V.K. Aggarwal / Appl. Math. Comput. 161 (2005) 973–987

8 < 4; Bi ðxj Þ ¼ 1; : 0;

if i ¼ j; if i  j ¼ 1; if i  j ¼ 2:

ð2:8Þ

and that Bi ðxÞ ¼ 0 for x P xiþ2 and x 6 xi2 . Similarly we can show that 8 if i ¼ j; < 0; B0i ðxj Þ ¼  3~h ; if i  j ¼ 1; : 0; if i  j ¼ 2:

ð2:9Þ

and 8 > ; < 12 ~ h2 00 Bi ðxj Þ ¼ ~62 ; h > : 0;

if i ¼ j; if i  j ¼ 1; if i  j ¼ 2:

ð2:10Þ

Each Bi ðxÞ is also a piece-wise cubic with knots at p, and Bi ðxÞ 2 X . Let X ¼ fB1 ; B0 ; B1 ; . . . ; BN þ1 g and let U3 ðpÞ ¼ span X. The functions X are linearly independent on ½0; 1 , thus U3 ðpÞ is ðN þ 3Þ-dimensional. Even one can show that U3 ðpÞ sp X . Let L be a linear operator whose domain is X and whose range is also in X . Now we define SðxÞ ¼ c1 B1 ðxÞ þ c0 B0 ðxÞ þ c1 B1 ðxÞ þ    þ cN BN ðxÞ þ cN þ1 BN þ1 ðxÞ: ð2:11Þ Then force SðxÞ to satisfy the collocation equations plus the boundary conditions. We have L1 Sðxi Þ ¼ Zðxi Þ

ð2:12Þ

0 6 xi 6 N ;

and Sð0Þ ¼ a0 ;

Sð1Þ ¼ a1 :

ð2:13Þ

On solving Eq. (2.12) we get ci1 ðB00i1 ðxi Þ þ Wi Bi1 ðxi ÞÞ þ ci ðB00i ðxi Þ þ Wi Bi ðxi ÞÞ þ ciþ1 ðB00iþ1 ðxi Þ þ Wi Biþ1 ðxi ÞÞ ¼ Zi

8i ¼ 0; 1; 2; . . . ; N ;

ð2:14Þ

where W ðxi Þ ¼ Wi and Zðxi Þ ¼ Zi . Simplifying Eq. (2.14) we get ð6 þ Wi ~ h2 Þci1 þ ð12 þ 4Wi ~ h2 Þci þ ð6 þ Wi ~h2 Þciþ1 ¼~ h2 Zi ;

8i ¼ 0; 1; . . . ; N :

ð2:15Þ

The given boundary conditions (2.13) becomes c1 þ 4c0 þ c1 ¼ a0 ;

ð2:16Þ

M.K. Kadalbajoo, V.K. Aggarwal / Appl. Math. Comput. 161 (2005) 973–987

977

and cN 1 þ 4cN þ cN þ1 ¼ a1 :

ð2:17Þ

Eqs. (2.15), (2.16) and (2.17) lead to a ðN þ 3Þ  ðN þ 3Þ tridiagonal system t with ðN þ 3Þ unknowns CN ¼ ðc1 ; c0 ; . . . ; cN þ1 Þ (where t stands for transpose). Now eliminating c1 from first equation of (2.15) and (2.16): we find 36 c0 ¼ Z0 ~ h2  a0 ð6 þ W0 ~ h2 Þ:

ð2:18Þ

Similarly, eliminating cN þ1 from the last equation of (2.15) and from (2.17), we find 36 cN ¼ ZN ~ h2  a1 ð6 þ WN ~ h2 Þ:

ð2:19Þ

Coupling Eqs. (2.18) and (2.19) with the second through ðN  1Þst equations of (2.15), we are lead to the system of ðN þ 1Þ linear equations TxN ¼ dN in the t ðN þ 1Þ unknowns xN ¼ ðc0 ; . . . ; cN Þ with right-hand side  2  h  a0 ð6 þ W0 ~h2 Þ; ~h2 Z1 ; ~h2 Z2 ; . . . ; ~ h2  a1 ð6 þ WN ~ h2 Þ ; dN ¼ Z0 ~ h2 ZN 1 ; ZN ~ and the coefficient matrix T given by 2

36 6 6 þ W1 ~h2 6 6 6 0 6 6 0 6 6 .. 6 . 6 4 0 0

0 12 þ 4W1 ~ h2 .. . 0 .. . ... ...

0 6 þ W1 ~ h2 .. . 6 þ Wi ~ h2 .. . 0 0

0 0 .. . 12 þ 4Wi ~ h2 .. . 6 þ WN1 ~ h2 0

... ... 0 6 þ Wi ~ h2 .. . 12 þ 4WN1 ~ h2 0

0 0

3

7 7 7 7 ... 7 7: ... 7 7 .. 7 . 7 25 ~ 6 þ WN1 h 36

ð2:20Þ Since W ðxÞ > 0, it is easily seen that the matrix T is strictly diagonally dominant and hence nonsingular. Since T is nonsingular, we can solve the system TxN ¼ dN for c0 ; c1 ; . . . ; cN and substitute into the boundary equations (2.16) and (2.17) to obtain c1 and cN þ1 . Hence this method of collocation applied to (2.5) and (2.6) using a basis of cubic B-spline has a unique solution SðxÞ. 2.1. Mesh selection strategy We form the piecewise-uniform grid in such a way that more points are generated in the boundary layers region than outside of it. We divide the interval ½0; 1 into three sub-intervals ð0; jÞ, ðj; 1  jÞ and ð1  j; 1Þ, where   pffiffi  j ¼ min 1=4;  ln N : ð2:21Þ W

978

M.K. Kadalbajoo, V.K. Aggarwal / Appl. Math. Comput. 161 (2005) 973–987

Assuming N ¼ 2r with r P 3 be the total no. of subintervals in the partitions of ½0; 1 . We divide the intervals ð0; jÞ and ð1  j; 1Þ each into N =4 equal mesh elements, while the interval ðj; 1  jÞ is divided into N =2 equal mesh elements. The resulting piecewise mesh depends upon just one parameter j. Now we consider ~ h ¼ h1 in the interval ½0; j ; 4j ; N xi ¼ xi1 þ h1

ð2:22Þ

h1 ¼

for i ¼ 1; 2; 3; . . . ; N =4;

where x0 ¼ 0. Also ~ ¼ h2 in the interval ½j; 1  j ; h 2ð1  2jÞ ; h2 ¼ N xi ¼ xi1 þ h2 for i ¼ N =4 þ 1; . . . ; 3N =4: Similarly ~ h ¼ h3 in the interval ½1  j; 1 ; 4j ; h3 ¼ N xi ¼ xi1 þ h3 for i ¼ 3N =4 þ 1; . . . ; N :

ð2:23Þ

ð2:24Þ

3. Derivation for uniform convergence First we prove the following lemma Lemma 1. The B-splines fB1 ; B0 ; . . . ; BN þ1 g defined in equation (2.7), satisfy the inequality N þ1 X

jBi ðxÞj 6 10;

0 6 x 6 1:

i¼1

Proof. We know that    X X N þ1 N þ1   Bi ðxÞ 6 jB ðxÞj:   i¼1 i  i¼1 At any ith nodal point xi , we have N þ1 X jBi j ¼ jBi1 j þ jBi j þ jBiþ1 j ¼ 6 < 10: i¼1

M.K. Kadalbajoo, V.K. Aggarwal / Appl. Math. Comput. 161 (2005) 973–987

979

Also we have jBi ðxÞj 6 4

and

jBi1 ðxÞj 6 4

for x 2 ½xi1 ; xi

similarly jBi2 ðxÞj 6 1

and

jBiþ1 ðxÞj 6 1

for x 2 ½xi1 ; xi :

Now for any point x 2 ½xi1 ; xi we have N þ1 X

jBi ðxÞj ¼ jBi2 j þ jBi1 j þ jBi j þ jBiþ1 j 6 10:

i¼1

Hence this proves the lemma.

h

Now to estimate the error kV ðxÞ  SðxÞk1 , let Yn be the unique spline interpolate from U3 ðpÞ to the solution V ðxÞ of our boundary value problem (2.5) and (2.6). If ZðxÞ 2 C 2 ½0; 1 , and V ðxÞ 2 C 4 ½0; 1 , and it follows from the De Boor-Hall error estimates that kDj ðV ðxÞ  Yn Þk1 6 cj h4j c ;

j ¼ 0; 1; 2;

ð3:1Þ

where hc ¼ maxfh1 ; h2 ; h3 g and cj Õs are independent of hc and N . Let Yn ðxÞ ¼

N þ1 X

bi Bi ðxÞ;

ð3:2Þ

i¼1

while SðxÞ ¼

N þ1 X

ci Bi ðxÞ;

i¼1

is our collocation solution. It is clear that from Eq. (3.1) that jL1 Sðxi Þ  L1 Yn ðxi Þj ¼ jZðxi Þ  L1 Yn ðxi Þ þ L1 V ðxi Þ  L1 V ðxi Þj 6 kh2c ;

ð3:3Þ

where k ¼ ½c2 þ c0 kW ðxÞk1 h2c . Also L1 Sðxi Þ ¼ L1 V ðxi Þ ¼ Zðxi Þ. Let L1 Yn ðxi Þ ¼ t Z^n ðxi Þ 8i, and Z^ n ¼ ðZ^n ðx0 Þ; Z^n ðx1 Þ; . . . ; Z^n ðxN ÞÞ . Now from the system TxN ¼ dN and equation (3.3) that the ith coordinate ½T ðxN  y n Þ i of T ðxN  y n Þ, y n ¼ ðb0 ; b1 ; . . . ; bN Þt , satisfy the inequality j½T ðxN  y n Þ i j ¼ h2c jZi  Z^i j 6 kh4c ;

ð3:4Þ

since ðTxN Þi ¼ h2c Zi and ðTy n Þi ¼ h2c Z^n ðxi Þ for i ¼ 1; 2; 3; . . . ; N  1. Also h2c Þ; and ðTxN Þ0 ¼ h2c Z0  a0 ð6 þ W0 ~ ðTy n Þ ¼ h2 Z^n ðx0 Þ  a0 ð6 þ W0 ~ h2 Þ 0

c

c

ð3:5Þ

980

M.K. Kadalbajoo, V.K. Aggarwal / Appl. Math. Comput. 161 (2005) 973–987

similarly ðTxN ÞN ¼ h2c ZN  a1 ð6 þ WN ~ h2c Þ; and ðTy n Þ ¼ h2 Z^n ðxN Þ  a1 ð6 þ WN ~ h2 Þ: N

c

ð3:6Þ

c

But the ith coordinate of ½T ðxN  y n Þ is the ith equation ð6 þ Wi h2c Þdi1 þ ð12 þ 4Wi h2c Þdi þ ð6 þ Wi h2c Þdiþ1 ¼ ni

8i ¼ 1; 2; . . . ; N  1;

ð3:7Þ

where di ¼ bi  ci , 1 6 i 6 N þ 1, and ni ¼ h2c ½Zðxi Þ  Z^n ðxi Þ

8i ¼ 1; 2; 3; . . . ; N  1:

ð3:8Þ

Obviously from Eq. (3.3) jni j 6 kh4c :

ð3:9Þ t

Let n ¼ max1 6 x 6 N 1 jni j. Also consider d ¼ ðd1 ; d0 ; . . . ; dN þ1 Þ , then we define ei ¼ jdi j, and ~e ¼ max1 6 i 6 N 1 jei j. Now Eq. (3.7) becomes ð12 þ 4Wi h2c Þdi ¼ ni þ ð6  Wi h2c Þðdi1 þ diþ1 Þ 8i ¼ 1; 2; . . . ; N  1: ð3:10Þ Taking absolute values with sufficiently small hc . We have ð12 þ 4Wi h2c Þei 6 n þ 2~eð6  Wi h2c Þ:

ð3:11Þ



Also 0 < W 6 W ðxÞ for all x. We get ð12 þ 4W  h2c Þei 6 n þ 2~eð6  W  h2c Þ 6 n þ 2~eð6 þ W  h2c Þ:

ð3:12Þ

In particular ð12 þ 4W  h2c Þ~e 6 n þ 2~eð6 þ W  h2c Þ;

ð3:13Þ

which gives 2W  h2c ~e 6 n 6 kh4c

implies ~e 6

kh2c : 2W 

ð3:14Þ

To estimate e1 , e0 , eN and eN þ1 , first we observe that the first equation of the system T ðxN  y n Þ ¼ h2c ½Z N  Z^ n (where Z N ¼ ðZ0 ; Z1 ; . . . ; ZN ÞÞ gives 36 d0 ¼ h2c ½Z0  Z^0 ;

ð3:15Þ

which gives e0 6

kh4c 36

ð3:16Þ

M.K. Kadalbajoo, V.K. Aggarwal / Appl. Math. Comput. 161 (2005) 973–987

981

similarly we can calculate eN 6

kh4c : 36

ð3:17Þ

Now e1 and eN þ1 can be calculated using Eqs. (2.16) and (2.17) as d1 ¼ ð4d0  d1 Þ and dN 1 ¼ ð4dN  dN 1 Þ e1 6

kh4c kh2 þ c ; 9 2W

ð3:18Þ

also eN þ1 6

kh4c kh2 þ c ; 9 2W

ð3:19Þ

using value k ¼ ½c2 þ c0 kW ðxÞk1 h2c , we get e¼

max fei g ) e 6 xh2c ;

ð3:20Þ

1 6 i 6 N þ1 2

6

where x ¼ c29hc þ 2Wk  where ðhc  0Þ. The above inequality enables us to estimate kSðxÞ  Yn ðxÞk1 , and hence kV ðxÞ  SðxÞk1 . In particular SðxÞ  Yn ðxÞ ¼

N þ1 X

ðci  bi ÞBi ðxÞ;

ð3:21Þ

i¼1

thus jSðxÞ  Yn ðxÞj 6 max jci  bi j

N þ1 X

jBi ðxÞj:

ð3:22Þ

i¼1

But N þ1 X

jBi ðxÞj 6 10;

06x61

ðusing Lemma 1Þ:

ð3:23Þ

i¼1

Combining Eqs. (3.20), (3.22) and (3.23), we see that kS  Yn k1 6 10xh2c :

ð3:24Þ

But kV  Yn k1 6 c0 h4c . Since kV  Sk1 6 kV  Yn k1 þ kYn  Sk1 , we see that kV  Sk1 6 Mh2c ; where M ¼ 10x þ c0 h2c . Combining the results we have proved.

ð3:25Þ

982

M.K. Kadalbajoo, V.K. Aggarwal / Appl. Math. Comput. 161 (2005) 973–987

Theorem 1. Let W ðxÞ, ZðxÞ 2 C 2 ½0; 1 , and W ðxÞ P W  > 0. Then the collocation approximation SðxÞ from the space U3 ðpÞ to the solution V ðxÞ of the boundary value problem (2.5) with (2.6) exists, and uniform error is given by kV ðxÞ  SðxÞk1 6 Mh2c ; for h2c sufficiently small and M is a positive constant (independent of hc and ). After knowing the value of V ðxÞ in the given domain, we can calculate the value of yðxÞ using Eq. (2.2).

4. Numerical results In this section the numerical results of some model problems are presented (Figs. 1 and 2). Problem 1. We solve Eq. (1.1) where aðxÞ ¼ 1;

bðxÞ ¼ 1;

and

gðxÞ ¼  cos2 ðpxÞ  2p2 cosð2pxÞ;

ð4:1Þ

having boundary conditions yð0Þ ¼ 0

and

yð1Þ ¼ 0;

ð4:2Þ

which has the exact solution [12]  pffiffi pffiffi   pffiffi  yðxÞ ¼ expð1ð1  xÞ= Þ þ expðx= Þ = 1 þ expð1= Þ  cos2 ðpxÞ:

N=64

0.4

S(x)=Appr. sol. u(x)=Exact sol.

0.2

0

Y– axis

Y–axis

S u

S(x)=Appr. . sol u(x)=Exact sol.

0.2

0 –0.2

–0.2

–0.4

–0.4

–0.6

–0.6

–0.8

–0.8

–1

N=128

0.4

S u

0

0.2

0.4

0.6

X–axis

0.8

1

1.2

–1

0

0.2

0.4

0.6

X– axis

Fig. 1. Problem 1 with uniform mesh for  ¼ 104 .

0.8

1

1.2

M.K. Kadalbajoo, V.K. Aggarwal / Appl. Math. Comput. 161 (2005) 973–987 N=64

N=128

1

0.8

0.8

Y–axis

Y–axis

1

0.6

0.6

0.4

0.4

0.2

0.2

0.2

0.4

0.6

0.8

1

S u

. S(x)=Appr. sol u(x)=Exact sol.

S u

. S(x)=Appr. sol u(x)=Exact sol.

0 0

983

1.2

0 0

0.2

0.4

X–axis

0.6

0.8

1

1.2

X–axis

Fig. 2. Problem 3 with uniform mesh for  ¼ 103 .

Problem 2. We solve Eq. (1.1) where aðxÞ ¼ 1;

bðxÞ ¼ 1 þ x;

  gðxÞ ¼ 40 xðx2  1Þ  2 ;

ð4:3Þ

having boundary conditions yð0Þ ¼ 0

and

yð1Þ ¼ 0;

ð4:4Þ

which has the exact solution [4] yðxÞ ¼ 40xð1  xÞ: Problem 3. We solve Eq. (1.1) where aðxÞ ¼ 1;

bðxÞ ¼ 1 þ xð1  xÞ; pffiffi pffiffi gðxÞ ¼ 1 þ xð1  xÞ þ ½2   x2 ð1  xÞ exp½ð1  xÞ=  pffiffi pffiffi 2 þ ½2   xð1  xÞ exp½x=  ;

having boundary conditions yð0Þ ¼ 0

and

yð1Þ ¼ 0;

ð4:5Þ

which has the exact solution [13] pffiffi pffiffi yðxÞ ¼ 1 þ ðx  1Þ exp½x=   x exp½ð1  xÞ=  : The maximum error for Problem 1 for different values of  and grid points N is presented in Table 1. The comparison in maximum error for Problem 2 with the existing methods for different values of  and grid points N is presented in Tables 2 and 3. The maximum error for Problem 3 for different values of  and grid points N is presented in Table 4.

984

 ¼ 2k

N ¼ 16

N ¼ 32

N ¼ 64

N ¼ 128

N ¼ 256

N ¼ 512

N ¼ 1024

k¼4 6 8 10 12 14 16 18 20 25 30

7.098E)03 4.073E)03 6.305E)02 7.752E)02 6.341E)02 6.237E)02 6.251E)02 6.250E)02 6.250E)02 6.250E)02 6.250E)02

1.779E)03 1.016E)03 2.441E)03 5.022E)02 3.066E)02 3.174E)02 3.119E)02 3.124E)02 3.125E)02 3.125E)02 3.125E)02

4.450E)04 2.541E)04 9.336E)04 3.928E)03 2.021E)02 1.576E)02 1.581E)02 1.560E)02 1.562E)02 1.562E)02 1.562E)02

1.112E)04 6.353E)05 2.323E)04 9.612E)04 3.803E)03 6.303E)03 7.909E)03 7.871E)03 7.804E)03 7.812E)03 7.812E)03

2.782E)05 1.588E)05 5.803E)05 2.392E)04 9.633E)04 5.367E)03 3.468E)03 3.940E)03 3.921E)03 3.906E)03 3.906E)03

6.955E)06 3.970E)06 1.450E)05 5.028E)05 2.768E)04 9.917E)04 9.731E)04 1.826E)03 1.963E)03 1.952E)03 1.953E)03

1.738E)06 9.926E)07 3.626E)06 1.276E)05 5.988E)05 2.398E)04 9.635E)04 6.840E)04 9.404E)04 9.759E)04 9.765E)04

M.K. Kadalbajoo, V.K. Aggarwal / Appl. Math. Comput. 161 (2005) 973–987

Table 1 Maximum error for Problem 1 using fitted mesh

M.K. Kadalbajoo, V.K. Aggarwal / Appl. Math. Comput. 161 (2005) 973–987

985

Table 2 Comparison of maximum errors for Problem 2 fitted mesh: N ¼ 16 

Miller [6]

Niijima [5]

Niijima [4]

Our method (uniform mesh)

0.1E)03 0.1E)04 0.1E)05 0.1E)06 0.1E)07 0.1E)08 0.1E)09

0.25E)01 0.21E)01 0.70E)02 0.75E)03 0.74E)04 0.67E)05 0.00E+00

0.26E)01 0.24E)01 0.17E)01 0.69E)02 0.23E)02 0.76E)03 0.24E)03

0.65E)04 0.36E)04 0.33E)04 0.26E)04 0.20E)04 0.20E)04 0.11E)04

1.776E)15 1.776E)15 1.776E)15 1.776E)15 1.783E)15 1.783E)15 1.784E)15

Table 3 Comparison of maximum errors for Problem 2 Fitted Mesh: N ¼ 32 

Miller [6]

Niijima [5]

Niijima [4]

Our method (uniform mesh)

0.1E)03 0.1E)04 0.1E)05 0.1E)06 0.1E)07 0.1E)08 0.1E)09

0.64E)02 0.61E)02 0.41E)02 0.77E)03 0.76E)04 0.67E)05 0.00E+00

0.65E)02 0.64E)02 0.56E)02 0.31E)02 0.12E)02 0.38E)03 0.13E)03

0.59E)04 0.21E)04 0.35E)04 0.39E)04 0.21E)04 0.21E)04 0.14E)04

1.776E)15 1.776E)15 1.776E)15 1.776E)15 1.776E)15 1.788E)15 1.789E)15

Table 4 Maximum error for Problem 3 using fitted mesh  ¼ 2k

N ¼ 16

N ¼ 32

N ¼ 64

N ¼ 128

N ¼ 256

N ¼ 512

N ¼ 1024

k¼4 6 8 10 12 14 16 18 20 25 30

2.000E)03 5.101E)03 1.950E)02 4.891E)02 5.731E)02 6.040E)02 6.153E)02 6.202E)02 6.230E)02 6.251E)02 6.251E)02

4.908E)04 1.302E)03 5.210E)03 2.012E)02 2.691E)02 2.950E)02 3.050E)02 3.091E)02 3.112E)02 3.123E)02 3.124E)02

1.225E)04 3.137E)04 1.120E)03 4.219E)03 1.234E)02 1.432E)02 1.510E)02 1.542E)02 1.551E)02 1.560E)02 1.562E)02

3.063E)05 7.842E)05 2.750E)04 1.102E)03 4.120E)03 6.890E)03 7.401E)03 7.720E)03 7.723E)03 7.810E)03 7.810E)03

7.657E)06 1.960E)05 6.877E)05 2.579E)04 9.988E)04 3.001E)03 2.901E)03 2.912E)03 2.912E)03 2.912E)03 2.192E)03

1.884E)06 4.824E)06 1.692E)05 6.329E)05 2.447E)04 9.196E)04 9.112E)04 9.070E)04 9.049E)04 9.032E)04 9.029E)04

4.786E)07 1.225E)06 4.297E)06 1.607E)05 6.208E)05 2.442E)04 2.836E)04 2.823E)04 2.816E)04 2.811E)04 2.810E)04

5. Discussion We have described a B-spline collocation method for solving self-adjoint singularly perturbed problems using fitted mesh technique. We have first transformed the original problem into normal form and then solved using the method. If we do not transform problem (1.1) to the normal form and if we discretize (1.1) directly then after getting the tridiagonal system, one can ob-

986

M.K. Kadalbajoo, V.K. Aggarwal / Appl. Math. Comput. 161 (2005) 973–987

serve that out of the two diagonal elements, viz., super and sub-diagonal elements one is exponentially large whereas the other is exponentially small as the parameter  ! 0 which cause the illconditioning in the system. This happens when aðxÞ is not a constant. These type of problems are very important physically and difficult to solve because of two boundary layers at both of the end points. In Tables 1 and 4, we have shown the maximum error for the problems 1 and 3 using B-spline collocation method with fitted mesh respectively. As we decrease the value of , to study the behaviour of the solution at the boundary layer, we need large number of mesh point in that region. So we conclude that uniform mesh is not a good technique for such kind of problems, because number of mesh points in the boundary layer region should be much more than that from the outer region. In Tables 2 and 3 we have compared our results with the results of [4,5] and it is clear that our results are much better. It can be seen from the respective Tables that if we take more points in the boundary layer region then we obtain better results which implies that the use of piecewise-uniform mesh is quite advantageous. To further corroborate the applicability of the proposed method, graphs have been plotted in Figs. 1–2 for Problems 1 and 3 for x 2 ½0; 1 versus the computed solution obtained at different values of x for a fixed . 6. Conclusion As is evident from the numerical results, this method gives Oðh2 Þ accuracy. The results obtained using this method are better than using the stated existing methods with the same no. of knots. Also this method produce a spline function which may be used to obtain the solution at any point in the range, whereas the finite-difference methods [4,5], gives the solution only at the chosen knots. References [1] M.H. Potter, H.F. Weinberger, Maximum principles in differential equations, Prentice-Hall, Inc., Englewood Cliffs, New Jersey, 1967. [2] I.P. Boglaev, A variational difference scheme for a boundary value problems with a small parameter in the highest derivative, USSR, Comput. Math. Math. Phys. 21 (4) (1981) 71–81. [3] A.H. Schatz, L.B. Wahlbin, On the finite element method for singularly perturbed reaction diffusion problems in two and one dimension, Math. Comput. 40 (1983) 47–89. [4] K. Niijima, On a three-point difference scheme for a singular perturbation problem without a first derivative term II, Mem. Numer. Math. 7 (1980) 11–27. [5] K. Niijima, On a three-point difference scheme for a singular perturbation problem without a first derivative term I, Mem. Numer. Math. 7 (1980) 1–10. [6] J.J.H. Miller, On the convergence, uniformly in , of difference schemes for a two-point boundary singular perturbation problem, in: Numerical Analysis of Singular Perturbation Problems, Academic Press, New York, 1979, pp. 467–474.

M.K. Kadalbajoo, V.K. Aggarwal / Appl. Math. Comput. 161 (2005) 973–987

987

[7] K. Surla, V. Jerkovic, Some possibilities of applying spline collocations to singular perturbation problems, Numerical methods and Approximation Theory, vol. II, Novisad, 1985, pp 19–25. [8] M. Sakai, R.A. Usmani, On exponential splines, J. Approx. Theory 47 (1986) 122–131. [9] L.L. Schumaker, Spline Functions: Basic Theory, Krieger Publishing Company, Florida, 1981. [10] J.M. Ahlberg, E.N. Nilson, J.L. Walsh, The Theory of Splines and Their Applications, Academic Press, New York, 1967. [11] G. Micula, Handbook of Splines, Kluwer Academic Publishers, Dordrecht, The Netherlands, 1999. [12] E.P. Doolan, J.J.H. Miller, W.H.A. Schilders, Uniform Numerical Methods for Problems with Initial and Boundary Layers, Boole Press, Dublin, 1980. [13] V. Vukoslavcevic, K. Surla, Finite element method for solving self-adjoint singularly perturbed boundary value problems, Math. Montisnigri. VII (1996) 69–86.