- Email: [email protected]
A computational procedure for solving a chemical flow reactor problem using shooting method
NORTH- HOLLAND
A Computational Procedure for Solving a Chemical Flow Reactor Problem Using Shooting Method K. Selvakumar
Department o] Mathematics Adhiparasakthi Engineering College Melmaruvathur-603 301 Chengai M.G.R. District Tamil Nadu, India Transmitted by John Cacti
ABSTRACT This paper presents a new computational procedure for solving a chemical flow reactor problem using shooting method. The convergence of the method is discussed and the numerical results are given to illustrate the computational procedure.
1.
INTRODUCTION Consider the following chemical flow reactor problem [1] L~(t) - ~ " + a(t)u' - b(t)u = / ( t ) ,
0 < t < 1, Bou(O) - -u'(0) = ¢I, Blu(1) - u(1) + eu'(1) = ¢~,
(la) (lb)
where e > 0 is a small parameter and a, b, and f are smooth functions satisfying a(t) >_ a > 0 and b(t) >_0 with respect to t. The solution of the problem (la, b) converges uniformly to the solution uo(t) of the reduced problem of (la, b) in the interval [0, 1] when ~ goes to zero: A P P L I E D M A T H E M A T I C S A N D C O M P U T A T I O N 68:27-40 (1995)
(~) Elsevier Science Inc., 1995 655 Avenue of the Americas, New York, NY 10010
0096-3003/95/$9.50 SSDI 0096-3003(94)00082-F
28
K. SELVAKUMAR
a(t)Uro(t)
-
b(t)uo(t) = f ( t ) ,
0
(2a) (2b)
Uo(1) = ¢2.
But their respective derivatives, in general, do not converge uniformly as goes to zero at t = 0. Because of this, classical schemes may not produce good approximations to the derivatives especially when e is small as in Table 1. T h e m a x i m u m principle and the stability result [2, 3] for the solution of the problem (la, b) are given in the following lemmas: LEMMA 1. Let v be a smooth function satisfying Lv(t) < 0, 0 <_ t < 1, Boy(O) >_ O, and BlV(1) > 0. Then v(t) >_ 0 for all 0 < t < 1. LEMMA 2. Let v be a smooth function. uniform stability estimate:
Then we have the following
{v(t)l <_ C(tBov(O)I + [Blv(1)] + max ILv(t)O t~[o,y] for all t in [0, 1]. The main aim of this paper is to solve the problem (la, b) using shooting method, t h a t is, the given problem (la, b) can be transformed into an initial value problem (IVP) then the resultant IVP is solved numerically using exponentially fitted finite difference schemes given in [3].
2.
SCHEMES F O R INITIAL VALUE P R O B L E M S
In this section, two finite difference schemes with variable and constant fitting factors for the IVP Ui - - U 2 = 0 ,
U 1(0) = ¢3
C ~ -~- a ( t ) u 2 -- b(~)Ul = f ( t ) ,
u2(0) = -(~1,
(3a) (3b)
where Ul = u and u2 ---- U I is given. The scheme with variable fitting factor is D+Ul,i - u2,~ = ((UlR(ti+l) - Uln(ti))/h) - u2R(ti) sa~(p)D+u% i + a(t~)u2, ~ - b(ti)ul, i = f h
(4a) (45)
Chemical Flow Reactor
29
TABLE
TEST
PROBLEM* --U'(0)
t
1
: Z U tr -F U r -- U ---- 0 , 0 < t < 1, ---- 0 , U ( 1 ) + ~ U ' ( 1 )
ui
---- 1
u(ti)
u ( t i ) - ui
0.00000E
+ 00
4.05826E
- 01
3.71469E
- 01
3.43570E
- 02
6.25000E
- 02
4.05826E
- 01
3.91354E
- 01
1.44721E
- 02
1.25000E
- 01
4.27691E
- 01
4.16332E
- 01
1.13598E
- 02
1.87500E
- 01
4.53751E
- 01
4.42911E
- 01
1.08399E
- 02
2.50000E
- 01
4.81793E
- 01
4.71188E
- 01
1.06058E
- 02
3.12500E
- 01
5.11620E
- 01
5.01269E
- 01
1.03506E
- 02
3.75000E
- 01
5.43300E
- 01
5.33271E
- 01
1.00282E
- 02
4.37500E
- 01
5.76942E
- 01
5.67317E
- 01
9.62520E
- 03
5.00000E
- 01
6.12668E
- 01
6.03536E
- 01
9,13203E
- 03
5.62500E
- 01
6.50605E
- 01
6.42067E
- 01
8.53872E
- 03
6.25000E
- 01
6.90892E
- 01
6.83058E
- 01
7.83479E
- 03
6.87500E
- 01
7.33674E
- 01
7.26666E
- 01
7.00849E
- 03
7.50000E
- 01
7.79105E-
01
7.73058E
- 01
6.04725E
- 03
8.12500E
- 01
8.27349E
- 01
8.22411E
- 01
4.93759E
- 03
8.75000E
- 01
8.78580E
- 01
8.74916E
- 01
3.66431E
- 03
9.37500E
- 01
9.32984E
- 01
9.30772E
- 01
2.21139E
- 03
1.00000E
+ 00
9.90756E
- 01
9.90195E
- 01
5.61357E
- 04
0.00000E
+ 00
0.00000E
+ 00
0,00000E
+ 00
0.00000E
+ 00
6.25000E
- 02
2.53641E
+ 00
3.86849E
- 01
2.14956E
+ 00
1.25000E
- 01
- 1.07797E
+ 01
4.12248E
- 01
1.11920E
+ 01
1.87500E
- 01
5.92668E
+ 01
4.38568E
- 01
5.88282E
+ 01
2.50000E
- 01
+ 02
4.66567E
- 01
3.08781E
+ 02
3.12500E
- 01
+ 03
4.96354E
- 01
1.62117E
+ 03
3.75000E
- 01
+ 03
5.28042E
- 01
8.51106E
+ 03
4.37500E
- 01
+ 04
5.61754E
- 01
4.46831E
+ 04
5.00000E
- 01
+ 05
5.97617E
- 01
2.34586E
+ 05
5.62500E
- 01
+ 06
6.35771E
- 01
1.23158E
+ 06
6.25000E
- 01
+ 06
6.76360E
- 01
6A6579E
+ 06
-3.08315E 1.62166E -8.51053E 4.46837E -2.34586E 1.23158E -6A6579E
6,87500E
- 01
3.39454E
+ 07
7.19540E
- 01
3.39454E
+ 07
7.50000E
- 01
- 1.78213E
+ 08
7.65477E
- 01
1.78213E
+ 08
8.12500E
- 01
9.35619E
+ 08
8.14347E
- 01
9.35619E
+ 08
8.75000E
-- 01
--4.91200E
+ 09
8.66337E
- 01
4.91200E
-4- 0 9
9.37500E
- 01
2.57880E
+ 10
9.21645E
- 01
2.57880E
+ 10
1.00000E
+ 00
9.24361E
- 01
9.80496E
- 01
5.613575E-02
*h =
1/16,
~ = 0.01
30
K. SELVAKUMAR Ul,0 : ¢3,
U2,0 ---- - - ¢ 1 ,
(4C)
where ci(P) -- pa(ti)/[1 - e x p ( - p a ( t i ) ) ] , f h = a( ti )u2R ( ti+ l ) -- b( ti)UlR ( ti )
p ---- h i e
(4d) (4e)
and Uln and u2a are the solution of the reduced problem of (3a, b) u~R - u2R = 0,
u l R ( 0 ) = ¢3
(ba)
a(t)u2R - b(t)ulR = f(t).
(55)
To reduce the computational time of the scheme (4a-e), a scheme with constant fitting factor is defined as follows: D+Ul# - u2# = ((UlR(~i+l) -- UlR(ti))/h) - u2R(ti) e a ( p ) D + u 2 # + a(~i)?L2,i+l - b(ti)l$1,i -- f h Ui,0 : ¢3,
U2,0 ~- - - ¢ 1 ,
(6a) (65) (6C)
where G(p) = pa(O)/[exp(pa(O))
- 1],
p = h/e
and UlR, u2R, and f h are defined as in (ba, b) and (4e). The above schemes are proposed in [3] for the IVP schemes are consistent, stable, and uniformly convergent reflect the asymptotic properties of the solution of (la, zero. The application of the schemes (4a-e) and (6a-d) next section.
3.
(6d)
(3a, b). These of order one and b) as e goes to are given in the
D E S C R I P T I O N OF T H E M E T H O D
In this section a new computational procedure to solve the problem (la, b) is presented using shooting method. For the problem (la, b), we make an initial guess ¢3 for u(0) from the zeroth order asymptotic expansion for the solution of the problem (la, b): u(t, c) -= Uo($) + c[(¢1 - Ulo(O))/a(O)]exp(-a(O)t/c),
(7)
where uo(t) is the solution of the problem (2a, b). Now, we denote u(t, ¢3) as the solution of the IVP ~ u " + a ( 0 u ' - b(t)u = f ( t ) ,
0 < t < 1
(8a)
31
Chemical Flow Reactor
~t(0) ~-~ (~3,
U'(0) ---~ --¢1
(8b)
such that [4] u(1, ¢3) + eu'(1, ¢3) - ¢2 = O(~).
(9)
As ~ goes to zero, u(1, ¢3) +eu'(1, ¢ 3 ) - ¢ 2 will tend to zero. The IVP (8a, b) can be reduced into a system of IVPs for first-order ordinary differential equations in which the first equation does not contain the parameter e and the second equation does contain a small parameter multiplied at the first derivative,
(10a)
u~ - u2 = 0,
eu~ + a(t)u2 - b(t)ul --- f(t), ~1(0) = ¢3,
where Ul : U and derivative but not given problem (la, solved numerically for (10a-c) is
0< t < 1
(10b)
(lOc)
us(0) = - ¢ 1 ,
u2 = u r. Here the nonuniformity occurs only at the in the solution. This is the similarity between the b) and the IVP (10a-c). Now the IVP (10a-c) can be using the scheme (4a~e). The finite difference scheme
D + u l # - u2,~ = ((UlR(t~+l)
-
-
ulR(t~))/h) - U2R(t,)
ea~(p)D+u2# + a(ti)u2,i - b(ti)ul,i = f h
ul,0 = ¢3,
i :> 0
u2,0 = --¢1,
(lla) (115) (Ilc)
where ul# = ui, cry(p), ulR, u2R, and fh are defined as in (4d), (5a-b), and (4e). Hence the given problem (la, b) can be solved numerically using shooting method. In the following the convergence of the solution u(t, ¢3) of the IVP (8a, b) to the solution u(t) of the problem (la, b) is given. THEOREM 1. Let u(t) and u(t,¢3) are the solutions of the problems (la, 5) and (8a, b), respectively. Then, for all t in [0, 1], we have
(12)
lu(t) - u(t, ¢3)t - c~.
PROOF. The proof of this theorem follows from the stability result (Lemma 2). From (la, b) and (8a, b) we have L[u(t) - u(t, ¢3)] : Lu(t) - Lu(t, 63) = f(t) - f(t) Bo[u(0)
= O,
- ~t(0, ¢ 3 ) ] ~-~ S o u ( O ) - S o u ( O , ¢ 3 ) m_ ¢1 - ¢ 1 = 0,
32
K. S E L V A K U M A R
and Bl[U(1) - u(1, ¢3)] - Blu(1) - Bltt(1, ¢3) = ¢2 - BlU(1, Ca) = ¢2 -- [U(1, ¢ 3 ) -t- eu'(1, ¢3)1.
From (9) we have IB1[u(1) - u(1, ~b3)]l _< 1¢2 - [u(1, ¢3) -~- eu'(1, ¢3)]I <_Ce. Now, using stability result, we have lu(t) - u(t, ¢3)1 --- Ce.
REMARK. W h e n e goes to zero, the solution u(t, ¢3) of (83, b) converges to the solution u(t) of (la, b). T h e following t h e o r e m gives the main result of this section, t h a t is, the convergence of the c o m p u t a t i o n a l m e t h o d . THEOREM 2. Let u(t) and Hi be the solutions of (la, b) and (lla~c), respectively. Then, we have, tu(td - u~l _< C(e + h).
PROOF.
Prom [3] we have the estimate for the solution of the scheme
(1 la-c) as
lu(t~, ¢3) - u~l <_ Ch. Prom T h e o r e m 1 we have, pu(t~) - u(t,, ~3)1 --- c e .
Therefore, I~(td - u~l < lu(td - u(t~, ¢3)1 + Ju(t~, ¢3) - ~ l <_ Ce + Ch = C(e + h).
33
Chemical Flow Reactor
REMARK. When the scheme (6a-d) is applied to the IVP (10a-c) the computational time will be reduced at each nodal points. From [3] we have lu(t~, ¢ 3 ) -
~1 < Ch,
where u(ti, ¢3) is the solution of the IVP (Sa, b) and us is the solution of the scheme (6a-d). Therefore, using the estimate (12), lu(t~) - u~l <_ l u ( t j
u(t~, ¢z)l
-
+
lu(t~, ¢3)
- u,I
< C(c + h).
4,
T E S T EXAMPLES AND NUMERICAL RESULTS EXAMPLE 1.
Consider the following homogeneous problem ~u" + u' = 0,
0
- u ' ( 0 ) = ¢1,
(13a)
u(1) + cu'(1) = ¢2,
(135)
where ¢1 = 0 and ¢2 = 0. The initial guess ¢3 for u(0) is given by u(0) = ¢2 +e¢1 from the zeroth order asymptotic expansion u(t,
~) = ¢2 + c ¢ l e x p ( - t / ~ ) .
The corresponding IVP is e u Ir + u ~ = 0,
(14a)
O< t < 1
u(0) = ¢3 = ¢2 + ~¢1,
~'(0) = - ¢ 1 ,
(14b)
whose solution is u(t,
¢3) = (¢3 - 6¢1) + c ¢ l e x p ( - t / ~ )
such that u(1, ¢3) + ~u'(1, ¢3) - ¢2 = 0.
The IVP (14a-b) can be reduced into a system of the form u~ - u2 = 0,
cu~ + u2 = 0,
ul(0)
= ¢3
u2(0) = -¢1,
(15~)
(15b)
34
K. SELVAKUMAR
where ul = u and u2 = u ~. The IVP (15a, b) can be solved numerically using the schemes given in Section 2. It is observed that the scheme (4a-d) reduces to the scheme (6a-d) when applied to the IVP (13a, b). Numerical results are given in Table 2. It is noted that this computational procedure solves exactly the sample problem (13a, b). EXAMPLE 2. able coefficient
Consider the following homogeneous problem with vari-
: u " + (1 + t)u' = 0,
0 < t < 1,
-ut(0) ~- (/)1,
(16a)
U(1) + eU'(1) = ¢2,
(16b)
where ¢1 = 0 and ¢2 = 1.5. The initial guess ¢3 for u(0) is given by u(0) = ¢2 + e e l , from the zeroth order asymptotic expansion
u(t, :) = ¢2 + : ¢ : e x p ( - t / e ) . The corresponding IVP is : u " + (1 + t)u' = 0,
0 < t < 1,
(17a)
u(0) = ¢3 = ¢2 + :¢1,
u'(0) = - ¢ 1 ,
(17b)
whose exact solution is
u(t, ¢3)
---- (¢3 -- : ¢ 1 ) q- : ¢ l e x p ( - - ( t
q-
t2/2)/C)
such that
u(1, Ca) + Eu'(1, Ca) - ¢2 = - : ¢ i e x p ( - 3 / 2 ~ ) = O(:). The IVP (17a, b) can be reduced into a system of the form u~ - u2 = 0, :U~ -~- (1 -[- t)~t2 ~-- 0,
u,(0) = ¢3
(18a)
~2(0) ---- --(~1,
(lSb)
where ul = u and u2 = u ~. Now the IVP (18a, b) can be solved numerically using the schemes (4a-d) and (6a-d). Table 3 gives the numerical result using the schemes (4a-d) and (6a-d), respectively, for the IVP (18a, b). EXAMPLE 3.
Consider the following nonhomogeneous problem
:u jt + u j = 1 + 2t,
0< t < 1
(19a)
Chemical Flow Reactor
35
TABLE 2* *
u~
u(t~)
u ( t i ) - u~
0 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + O0
0 . 0 0 0 0 0 E + 00
6 . 2 5 0 0 0 E - 02
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + O0
0 . 0 0 0 0 0 E + 00
1 . 2 5 0 0 0 E - 01
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00
0 . 0 0 0 0 0 E + 00
1 . 8 7 5 0 0 E - 01
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00
0 . 0 0 0 0 0 E + 00
2 . 5 0 0 0 0 E - O1
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00
0.0000OE + 0 0
3 . 1 2 5 0 0 E - 01
2 . 0 0 0 0 0 E + 00
2.O0000E + O0
0.0O000E + 00
3 . 7 5 0 0 0 E - 01
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00
0 . 0 0 0 0 0 E + 0O
4 . 3 7 5 0 0 E - 01
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00
0.0O000E + 0O
5.0000OE - 01
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00
0 . 0 0 0 0 0 E + 00
5 . 6 2 5 0 0 E - 01
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00
0 . 0 0 0 0 0 E + 00
6 . 2 5 0 0 0 E - 01 6 . 8 7 5 0 0 E - 01
2 . 0 0 0 0 0 E + 00 2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00 2 . 0 0 0 0 0 E + 00
0.0O000E + 00 0 . 0 0 0 0 0 E + 00
7 . 5 0 0 0 0 E - 01
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00
0 . 0 0 0 0 0 E + 00
8 . 1 2 5 0 0 E - Ol
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + O0
0 . 0 0 0 0 0 E + 00
8 . 7 5 0 0 0 E - 01
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00
0 . 0 0 0 0 0 E + 00
9 . 3 7 5 0 0 E - 01
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00
0 . 0 0 0 0 0 E + 00
1 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + O0
0 . 0 0 0 0 0 E + 00
*h = ¢ = 1/16, = O.O0000E + O0
¢3 = 2 . 0 0 0 0 0 E + 0 0 ,
u ( 1 , ¢ 3 ) + ¢ u ' ( 1 , ¢ 3 ) - ¢2
t
u~
u(t~)
u ( ~ ) - u~
0.0000OE + O0
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + O0
0 . 0 0 0 0 0 E + 00
6 . 2 5 0 0 0 E - 02
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00
0 . 0 0 0 0 0 E + 00
1 . 2 5 0 0 0 E - 01
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00
0 . 0 0 0 0 0 E + 00
1 . 8 7 5 0 0 E - 01
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00
0.00000E + 00
2 . 5 0 0 0 0 E - 01
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00
0 . 0 0 0 0 0 E + 00
3 . 1 2 5 0 0 E - 01
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00
0 . 0 0 0 0 0 E + 00
3 . 7 5 0 0 0 E - 01
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00
0 . 0 0 0 0 0 E + 00
4 . 3 7 5 0 0 E - 01
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00
0.00000E + 00
5 . 0 0 0 0 0 E - 01
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00
0 . 0 0 0 0 0 E + 00
5 . 6 2 5 0 0 E - 01
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00
0 . 0 0 0 0 0 E + 00
6 . 2 5 0 0 0 E - 01
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + O0
0 . 0 0 0 0 0 E + 00
6 . 8 7 5 0 0 E - 01
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00
0 . 0 0 0 0 0 E + 00
7 . 5 0 0 0 0 E - 01
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + O0
0.0000OE + 00
8 . 1 2 5 0 0 E - 01
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + O0
0 . 0 0 0 0 0 E + 00
8 . 7 5 0 0 0 E - 01
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00
0 . 0 0 0 0 0 E + 00
9 . 3 7 5 0 0 E - 01 1 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00 2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00 2 . 0 0 0 0 0 E + 00
0 . 0 0 0 0 0 E + 00 0 . 0 0 0 0 0 E + 00
*h = 1/16, e = 1.00000E-05, + ~ u ' ( 1 , ¢ 3 ) - ¢2 = O.00000E + 00
¢3 =
2.0000OE+00,
u(1,¢a)
36
K. S E L V A K U M A R
T A B L E 3*
0.00000E + 00
1.50000E + 00
1. 50000E + 00
0.00000E + 00
6. 25000E - 02
1.50000E + 00
1.50000E + 00
0.00000E + 00
1. 25000E -- 01
1.50000E q- 00
1.50000E -I- 00
0.00000E q- 00
1.87500E - 01
1.50000E + 00
1.50080E + 00
0.00000E + 00 0.00000E + 00
2.50000E - 01
1.50000E + 00
1.50000E + 00
3.12500E - 01
1.50000E + 00
1.50000E + 00
0.00000E + 00
3. 75000E - 01
1.50000E + 00
1,50000E + 00
0.00000E + 00
4 . 3 7 5 0 0 E - 01
1.50000E + 00
1.50000E + 00
0.00000E + 00
5 , 0 0 0 0 0 E - 01
1.50000E + 00
1,50000E + 00
0.00000E + 00
5.62500E - 01
1.50000E + 00
1.50000E + 00
0.00000E + 00
6.25000E - 01
1.50000E + 00
1.50000E ÷ 00
0.00000E + 00
6.87500E - 01
1.50000E + 00
1,50000E + 00
0.00000E + 00
7.50000E - 01
1.50000E + 00
1.50000E + 00
0.00000E + 00
8.12500E - 01
1.50000E + 00
1.50080E + 80
0.00000E + 00
8.75000E - 01
1.50000E + 00
1.50000E + 00
0.00000E + 00
9.37500E - 01
1.50000E + 00
1.50000E + 00
0.00000E + 00
1.00000E -b 00
1.50000E + 00
1.50000E + 00
0.00000E -t- 00
*h -- ~ --= 1/18, = 0.00000E + O0
¢3
=
1.50000E+00,
u(1, q~3) + Eu'(1,~b3) -- ¢2
t
ui
u(t~)
u ( t i ) - u~
0.00000E + 00
1.50000E + 00
1.50000E + 00
0.00000E + 00
6.25000E - 02
1.50000E + 00
1,50000E + 00
0,00000E + 00
1.25000E - 01
1.50000E + 00
1.50000E + 00
0.00000E + 00
1.87500E - 01
1.50000E + 00
1.50000E + 00
0.00000E -b 00
2.50000E - 01
1.50000E + 00
1.50000E + 00
0.00000E + 00
3.12500E - 01
1.50000E + 00
1.50000E -+ 00
0.00000E + 00
3.75000E - 01
1.50000E + 00
1.50000E + 00
0.00000E q- 00
4,37500E - 01
1.50000E + 00
1.50000E + 00
0.00000E + 00
5.00000E - 01
1.50000E -+-00
1.50000E + 00
0.00000E + 00
5,62500E - 01
1.50000E q- 00
1.50000E + 00
0.00000E q- 00
6.25000E - 01
1.50000E + 00
1.50000E + 00
0.00000E + 00
6.87500E - 01
1.50000E + 00
1.50000E + 00
0.00000E + 00
7.50000E - 01
1.50000E + 00
1.50000E + 00
0.00000E + 00
8.12500E - 01
1.50000E + 00
1.50000E + 00
0.00000E + 00
8.75000E - 01
1.50000E + 00
1.50000E + 00
0.00000E + 00
9.37500E - 01
1.50000E -b 00
1.50000E -b 00
0.00000E -b 00
1.00000E -t- 00
1.50000E q- 00
1.50000E + 00
0.00000E -k 00
* h = 1/16, c = 1.00000E-05, + ~u' ( 1, ¢3) - ¢2 = 0.00000E + 00
q~3 =
1.50000E+00,
u(1,¢3)
37
Chemical Flow Reactor
where
¢1
-'~
(tgb)
u(1) + 6u'(1) ----¢2,
=
0 and ¢2 = 3.
The initial guess ~b3 for u(0) is given by u(0) = ¢2 - 2 + ¢(¢1 the zeroth order asymptotic expansion U ( t , g ) = (~b2 -- 2 ) + t + t 2 4- g ( ¢ l
--
--
1) from
1)exp(-t/~).
The corresponding I V P is ~u" + u / = 1 + 2t,
0 < t < 1
u(0) = ¢3 ----¢2 - 2 + e ( ¢ 1 - 1),
(20a) u'(0) = - ¢ 1 ,
(20b)
whose solution is u ( t , ¢3) = ¢3 - ¢(¢z + 1 - 26)(1 - e x p ( - t / ¢ )
+ t + t 2 - 2~t
such that u(1, ¢3) + 6u'(1, ¢3) - ¢2 = - 6 , which is O(6). The IVP (20a, b) can be reduced into a system of the form ! U 1 - - U 2 = 0,
¢u~ + u2 = 1 + 2t,
Ul(0) = ¢3 u2(0) = - ¢ 1 ,
(21a)
(21b)
w h e r e u l = u a n d u2 = u t. Now the IVP (21a, b) can be solved by the schemes given in Section 2. It is observed that the scheme (4a-d) reduces to the scheme (6a-d) when applied to the IVP (21a, b). The numerical results are given in Table 4. A similar discussion can be carried out for the following test problem CU II 4- 'tt/ - - I t ~ O,
--u'(O) =
0 < t < t
u(1) + ¢u'(1) : ¢2,
(22a) (22b)
where ¢1 = 0 and ~b2 = 1. The numerical results are given in Table 5.
5.
CONCLUSION
From the above discussion, we conclude that the given chemical flow reactor problem can be solved by using shooting method. The given chemical reactor flow problem is transformed into an initial value problem using shooting method with an initial guess from zeroth order asymptotic
38
K. SELVAKUMAR
TABLE
t
*h =
u(t~)
u ( t D - ui
0 . 0 0 0 0 0 E 4- 0 0
9.37500E - 01
1.00000E4-
00
6.25000E
6.25000E
- 02
9.41406E - 01
1.02402E4-
00
8.26185E
- 02
1.25000E - 01
9.89758E -- 01
1 . 0 7 7 7 1 E 4- 0 0
8.79552E
- 02
- 02
1.87500E - 01
1.05940E + 00
1.14725E + 00
8.78540E - 02
2.50000E - 01
1.14181E + 00
1.22756E + 00
8.57525E - 02
3.12500E -- 01
1 . 2 3 3 8 6 E 4- 0 0
1 . 3 1 6 7 7 E 4- 0 0
8.29148E -- 02
3.75000E - 01
1 . 3 3 4 3 9 E 4- 0 0
1 . 4 1 4 2 0 E 4- 0 0
7.98066E - 02
4.37500E - 01
1.44298E + 00
1.51958E + 00
7.65988E - 02
5.00000E - 01
1 . 5 5 9 4 8 E 4- 0 0
1.63283E + 00
7.33542E - 02
5.62500E - 01
1 . 6 8 3 8 2 E 4- 0 0
L75391E4-00
6.25000E
- 01
1.81598E÷00
1.88282E+
6.87500E
- 01
1 . 9 5 5 9 6 E 4- 0 0
2 . 0 1 9 5 3 E 4- 0 0
6.35684E -- 02
7.50000E - 01
2 . 1 0 3 7 6 E 4- 0 0
2 . 1 6 4 0 6 E 4- 0 0
6.03030E - 02
8.12500E - 01
2.25937E + 00
2.31641E + 00
5.70374E - 02
8.75000E
- 01
2,42279E + 00
2.47656E4-
5.37715E - 02
9.37500E - 01
2.59403E + 00
2 . 6 4 4 5 3 E 4- 0 0
5.05056E - 02
1.00000E4-
2 . 7 7 3 0 7 E 4- 0 0
2 . 8 2 0 3 1 E 4- 0 0
4.72398E - 02
=
e
=
6.25000E
00 1/16,
¢3
=
9.37500E-01,
7.00963E - 02 00
u(1,~3)
00
6.68333E
4- e u ' ( 1 , ¢ 3 )
-
- 02
~2
-- 02
t
*h
4*
u~
~
~(t~)
~(t~) - ~i
0.00000E + 00
9.99990E - 01
1.00000E + 00
1.00136E - 05
6.25000E
- 02
1.00390E + 00
1.06640E + 00
6.24988E
1.25000E - 01
1 . 0 7 8 1 1 E 4- 0 0
1 . 1 4 0 6 1 E 4- 0 0
6.24975E - 02
- 02
1.87500E - 01
1.16015E + 00
1.22264E + 00
6.24963E - 02
2.50000E - 01
1.24999E + 00
1.31248E + 00
6.24950E - 02
3.12500E - 01
1 . 3 4 7 6 5 E 4- 0 0
1 . 4 1 0 1 4 E 4- 0 0
6.24938E - 02
3.75000E - 01
1.45311E + 00
1.51561E + 00
6.24925E - 02
4.37500E
- 01
1.56640E + 00
1.62889E + 00
6.24913E
- 02
5.00000E - 01
1 . 6 8 7 4 9 E 4- 0 0
1.74998E + 00
6.24900E
- 02
5.62500E - 01
1.81640E + 00
1.87889E + 00
6.24888E
- 02
6.25000E
- 01
1.95311E + 00
2.01560E + 00
6.24876E -- 02
6.87500E - 01
2.09765E + 00
2.16013E + 00
6.24864E -- 02
7.50000E - 01
2 . 2 4 9 9 9 E 4- 0 0
2.31247E + 00
6.24850E
8.12500E
- 01
2.41015E + 00
2.47263E + 00
6.24838E - 02
8.75000E - 01
2.57811E + 00
2.64060E + 00
6.24826E - 02
9.37500E - 01
2.75390E + 00
2 . 8 1 6 3 8 E 4- 0 0
6.24814E - 02
1.00000E + 00
2.93749E + 00
Z99997E
6.24800E
=
1/16,
÷ eu'(1,¢3)
e
=
1.00000E--05,
-~b2 = 1.00136E-05
¢3
=
4- 0 0
9.99990E+01,
u(1,¢3)
- 02
- 02
Chemical Flow Reactor
39
TABLE
t
5*
ul
u(tl)
u(tl) - ui
0.00000E
+00
3.90872E
-
01
3,87755E
- 01
3.11729E
- 03
6.25000E
- 02
3.91652E
- 01
3.96713E
- 01
5.06178E
- 03
1.25000E
-
01
4.08920E
- 01
4A5763E
-
02
6.84297E
- 03
1.87500E
-
01
4.32384E
- 01
4,39281E
-
01
6.89682E
- 03
2.50000E
-
01
4.58915E
- 01
4.65376E
-
01
6.46129E
- 03
3.12500E
-
01
4.87614E
- 01
4.93457E
-
01
5.84248E
- 03
3.75000E
- 01
5.18279E
- 01
5.23382E
-
01
5.10353E
- 03
4.37500E
- 01
5.50928E
- 01
5.55175E
- 01
4.24737E
- 03
5.00000E
- 01
5.85654E
-
01
5.88917E
- 01
3.26347E
- 03
5.62500E
- 01
6.22578E
-
01
6.24716E
- 01
2.13802E
- 03
6.25000E
- 01
6.61837E
- 01
6.62694E
- 01
8.56698E
- 04
6.87500E
- 01
7.03577E
-
01
7.02981E
- 01
5.96464E
- 04.
7.50000E
- 01
7.47955E
-
01
7.45717E
- 01
2.23833E
-- 03
8.12500E
- 01
7.95138E
-
01
7.91051E
- 01
4.08697E
- 03
8.75000E
- 01
8.45304E
- 01
8.39142E
- 01
6.16199E
- 03
9.37500E
-
01
8.98641E
- 01
8.90156E
- 01
8.48490E
- 03
1.00000E
4. 00
9.55351E
- 01
9.44272E
- 01
1.10788E
- 02
* h
=
=
•
e
8.03924E
1/16,
¢3
=
3.90872E+00,
u(1,¢3)
4- ~ u ' ( 1 , ¢ 3 )
-
¢2
- 03
t
u~
u(ti)
u(~i) - u~
0.00000E
+ 00
3.67916E
- 01
3.68038E
- 01
1.21415E
6,25000E
- 02
3.68650E
-
01
3.91724E
-
01
2.30743E
- 04 - 02
1.25000E
- 01
3.93909E
- 01
4.16977E
- 01
2.30682E
- 02
1.87500E
- 01
4.19360E
- 01
4.43858E
- 01
2.44982E
- 02
2.50000E
- 01
4.46545E
- 01
4.72472E
- 01
2.59270E
- 02
3.12500E
- 01
4.75486E
-
01
5.02930E
- 01
2.74442E
- 02
3.75000E
- 01
5.06302E
-
01
5.35352E
-
01
2.90496E
- 02
4.37500E
- 01
5.39115E
-
01
5.69864E
- 01
3.07490E
- 02
5.00000E
- 01
5.74053E
-
10
6.06600E
- 01
3.25474E
- 02
5.62500E
-
01
6.11254E
- 01
6.45705E
- 01
3.44509E
- 02
6.25000E
- 01
6.50866E
-
01
6.87331E
-
01
3.64656E
- 02
6.87500E
- 01
6.93043E
- 01
7.31641E
- 01
3.85979E
- 02
7.50000E
- 01
7.37925E
-
01
7.78807E
- 01
4.08546E
- 02
8.12500E
- 01
7.85770E
- 01
8.29013E
- 01
4.32430E
- 02
8.75600E
- O1
8.36685E
-
01
8.82456E
- 01
4.57709E
- 02
9.37500E
- O1
8.90898E
-
01
9.39344E
-
01
4.84463E
~ 02
1.00000E
4. 00
9.48622E
- 01
9.99900E
-
01
5.12775E
- 02
*h +
=
1/16,
eu'(1,¢a)
e -
=
1.00000E-04,
¢2 = 3,30031E-04
¢3
=
3.67916E+01,
u(1,¢a)
40
K. SELVAKUMAR
expansion for the solution. Then the resultant initial value problem is solved numerically using exponentially fitted finite difference schemes. The similarity between the given chemical reactor flow problem and the initial value problem is that the nonuniformity occurs only in the derivative but not in the solution. The advantage of this method is that there is no need for matrix inversion as in tridiagonal difference schemes and a three-point difference scheme is replaced by a two-point difference scheme (one-step method) using shooting method. When the exact solution of the initial value problem (Sa, b) is not known, using the numerical solution of (8a, b) one can check the condition (9). In Tables 2-5 the condition (9) is checked and given for the respective test problems to illustrate the procedure.
REFERENCES 1 J. V. Baxley, On singular perturbation of nonlinear two-point boundary value problems, Czechoslovak Math. J. 27:363-377 (1977). 2 E. P. Doolan, J. J. H. Miller, and W. H. A. Schilders, Uniform Numerical Methods for Problems with Initial and Boundary Layers, Boole Press, Dublin, 1980. 3 K. Selvakumar, Uniformly Convergent Difference Schemes for Differential Equations with a Parameter, Ph.D. Thesis, Bharathidasan University, India, 1992. 4 J.D. Lambert, Computational Methods in Ordinary Differential Equations, Wiley, 1973.
A Computational Procedure for Solving a Chemical Flow Reactor Problem Using Shooting Method K. Selvakumar
Department o] Mathematics Adhiparasakthi Engineering College Melmaruvathur-603 301 Chengai M.G.R. District Tamil Nadu, India Transmitted by John Cacti
ABSTRACT This paper presents a new computational procedure for solving a chemical flow reactor problem using shooting method. The convergence of the method is discussed and the numerical results are given to illustrate the computational procedure.
1.
INTRODUCTION Consider the following chemical flow reactor problem [1] L~(t) - ~ " + a(t)u' - b(t)u = / ( t ) ,
0 < t < 1, Bou(O) - -u'(0) = ¢I, Blu(1) - u(1) + eu'(1) = ¢~,
(la) (lb)
where e > 0 is a small parameter and a, b, and f are smooth functions satisfying a(t) >_ a > 0 and b(t) >_0 with respect to t. The solution of the problem (la, b) converges uniformly to the solution uo(t) of the reduced problem of (la, b) in the interval [0, 1] when ~ goes to zero: A P P L I E D M A T H E M A T I C S A N D C O M P U T A T I O N 68:27-40 (1995)
(~) Elsevier Science Inc., 1995 655 Avenue of the Americas, New York, NY 10010
0096-3003/95/$9.50 SSDI 0096-3003(94)00082-F
28
K. SELVAKUMAR
a(t)Uro(t)
-
b(t)uo(t) = f ( t ) ,
0
(2a) (2b)
Uo(1) = ¢2.
But their respective derivatives, in general, do not converge uniformly as goes to zero at t = 0. Because of this, classical schemes may not produce good approximations to the derivatives especially when e is small as in Table 1. T h e m a x i m u m principle and the stability result [2, 3] for the solution of the problem (la, b) are given in the following lemmas: LEMMA 1. Let v be a smooth function satisfying Lv(t) < 0, 0 <_ t < 1, Boy(O) >_ O, and BlV(1) > 0. Then v(t) >_ 0 for all 0 < t < 1. LEMMA 2. Let v be a smooth function. uniform stability estimate:
Then we have the following
{v(t)l <_ C(tBov(O)I + [Blv(1)] + max ILv(t)O t~[o,y] for all t in [0, 1]. The main aim of this paper is to solve the problem (la, b) using shooting method, t h a t is, the given problem (la, b) can be transformed into an initial value problem (IVP) then the resultant IVP is solved numerically using exponentially fitted finite difference schemes given in [3].
2.
SCHEMES F O R INITIAL VALUE P R O B L E M S
In this section, two finite difference schemes with variable and constant fitting factors for the IVP Ui - - U 2 = 0 ,
U 1(0) = ¢3
C ~ -~- a ( t ) u 2 -- b(~)Ul = f ( t ) ,
u2(0) = -(~1,
(3a) (3b)
where Ul = u and u2 ---- U I is given. The scheme with variable fitting factor is D+Ul,i - u2,~ = ((UlR(ti+l) - Uln(ti))/h) - u2R(ti) sa~(p)D+u% i + a(t~)u2, ~ - b(ti)ul, i = f h
(4a) (45)
Chemical Flow Reactor
29
TABLE
TEST
PROBLEM* --U'(0)
t
1
: Z U tr -F U r -- U ---- 0 , 0 < t < 1, ---- 0 , U ( 1 ) + ~ U ' ( 1 )
ui
---- 1
u(ti)
u ( t i ) - ui
0.00000E
+ 00
4.05826E
- 01
3.71469E
- 01
3.43570E
- 02
6.25000E
- 02
4.05826E
- 01
3.91354E
- 01
1.44721E
- 02
1.25000E
- 01
4.27691E
- 01
4.16332E
- 01
1.13598E
- 02
1.87500E
- 01
4.53751E
- 01
4.42911E
- 01
1.08399E
- 02
2.50000E
- 01
4.81793E
- 01
4.71188E
- 01
1.06058E
- 02
3.12500E
- 01
5.11620E
- 01
5.01269E
- 01
1.03506E
- 02
3.75000E
- 01
5.43300E
- 01
5.33271E
- 01
1.00282E
- 02
4.37500E
- 01
5.76942E
- 01
5.67317E
- 01
9.62520E
- 03
5.00000E
- 01
6.12668E
- 01
6.03536E
- 01
9,13203E
- 03
5.62500E
- 01
6.50605E
- 01
6.42067E
- 01
8.53872E
- 03
6.25000E
- 01
6.90892E
- 01
6.83058E
- 01
7.83479E
- 03
6.87500E
- 01
7.33674E
- 01
7.26666E
- 01
7.00849E
- 03
7.50000E
- 01
7.79105E-
01
7.73058E
- 01
6.04725E
- 03
8.12500E
- 01
8.27349E
- 01
8.22411E
- 01
4.93759E
- 03
8.75000E
- 01
8.78580E
- 01
8.74916E
- 01
3.66431E
- 03
9.37500E
- 01
9.32984E
- 01
9.30772E
- 01
2.21139E
- 03
1.00000E
+ 00
9.90756E
- 01
9.90195E
- 01
5.61357E
- 04
0.00000E
+ 00
0.00000E
+ 00
0,00000E
+ 00
0.00000E
+ 00
6.25000E
- 02
2.53641E
+ 00
3.86849E
- 01
2.14956E
+ 00
1.25000E
- 01
- 1.07797E
+ 01
4.12248E
- 01
1.11920E
+ 01
1.87500E
- 01
5.92668E
+ 01
4.38568E
- 01
5.88282E
+ 01
2.50000E
- 01
+ 02
4.66567E
- 01
3.08781E
+ 02
3.12500E
- 01
+ 03
4.96354E
- 01
1.62117E
+ 03
3.75000E
- 01
+ 03
5.28042E
- 01
8.51106E
+ 03
4.37500E
- 01
+ 04
5.61754E
- 01
4.46831E
+ 04
5.00000E
- 01
+ 05
5.97617E
- 01
2.34586E
+ 05
5.62500E
- 01
+ 06
6.35771E
- 01
1.23158E
+ 06
6.25000E
- 01
+ 06
6.76360E
- 01
6A6579E
+ 06
-3.08315E 1.62166E -8.51053E 4.46837E -2.34586E 1.23158E -6A6579E
6,87500E
- 01
3.39454E
+ 07
7.19540E
- 01
3.39454E
+ 07
7.50000E
- 01
- 1.78213E
+ 08
7.65477E
- 01
1.78213E
+ 08
8.12500E
- 01
9.35619E
+ 08
8.14347E
- 01
9.35619E
+ 08
8.75000E
-- 01
--4.91200E
+ 09
8.66337E
- 01
4.91200E
-4- 0 9
9.37500E
- 01
2.57880E
+ 10
9.21645E
- 01
2.57880E
+ 10
1.00000E
+ 00
9.24361E
- 01
9.80496E
- 01
5.613575E-02
*h =
1/16,
~ = 0.01
30
K. SELVAKUMAR Ul,0 : ¢3,
U2,0 ---- - - ¢ 1 ,
(4C)
where ci(P) -- pa(ti)/[1 - e x p ( - p a ( t i ) ) ] , f h = a( ti )u2R ( ti+ l ) -- b( ti)UlR ( ti )
p ---- h i e
(4d) (4e)
and Uln and u2a are the solution of the reduced problem of (3a, b) u~R - u2R = 0,
u l R ( 0 ) = ¢3
(ba)
a(t)u2R - b(t)ulR = f(t).
(55)
To reduce the computational time of the scheme (4a-e), a scheme with constant fitting factor is defined as follows: D+Ul# - u2# = ((UlR(~i+l) -- UlR(ti))/h) - u2R(ti) e a ( p ) D + u 2 # + a(~i)?L2,i+l - b(ti)l$1,i -- f h Ui,0 : ¢3,
U2,0 ~- - - ¢ 1 ,
(6a) (65) (6C)
where G(p) = pa(O)/[exp(pa(O))
- 1],
p = h/e
and UlR, u2R, and f h are defined as in (ba, b) and (4e). The above schemes are proposed in [3] for the IVP schemes are consistent, stable, and uniformly convergent reflect the asymptotic properties of the solution of (la, zero. The application of the schemes (4a-e) and (6a-d) next section.
3.
(6d)
(3a, b). These of order one and b) as e goes to are given in the
D E S C R I P T I O N OF T H E M E T H O D
In this section a new computational procedure to solve the problem (la, b) is presented using shooting method. For the problem (la, b), we make an initial guess ¢3 for u(0) from the zeroth order asymptotic expansion for the solution of the problem (la, b): u(t, c) -= Uo($) + c[(¢1 - Ulo(O))/a(O)]exp(-a(O)t/c),
(7)
where uo(t) is the solution of the problem (2a, b). Now, we denote u(t, ¢3) as the solution of the IVP ~ u " + a ( 0 u ' - b(t)u = f ( t ) ,
0 < t < 1
(8a)
31
Chemical Flow Reactor
~t(0) ~-~ (~3,
U'(0) ---~ --¢1
(8b)
such that [4] u(1, ¢3) + eu'(1, ¢3) - ¢2 = O(~).
(9)
As ~ goes to zero, u(1, ¢3) +eu'(1, ¢ 3 ) - ¢ 2 will tend to zero. The IVP (8a, b) can be reduced into a system of IVPs for first-order ordinary differential equations in which the first equation does not contain the parameter e and the second equation does contain a small parameter multiplied at the first derivative,
(10a)
u~ - u2 = 0,
eu~ + a(t)u2 - b(t)ul --- f(t), ~1(0) = ¢3,
where Ul : U and derivative but not given problem (la, solved numerically for (10a-c) is
0< t < 1
(10b)
(lOc)
us(0) = - ¢ 1 ,
u2 = u r. Here the nonuniformity occurs only at the in the solution. This is the similarity between the b) and the IVP (10a-c). Now the IVP (10a-c) can be using the scheme (4a~e). The finite difference scheme
D + u l # - u2,~ = ((UlR(t~+l)
-
-
ulR(t~))/h) - U2R(t,)
ea~(p)D+u2# + a(ti)u2,i - b(ti)ul,i = f h
ul,0 = ¢3,
i :> 0
u2,0 = --¢1,
(lla) (115) (Ilc)
where ul# = ui, cry(p), ulR, u2R, and fh are defined as in (4d), (5a-b), and (4e). Hence the given problem (la, b) can be solved numerically using shooting method. In the following the convergence of the solution u(t, ¢3) of the IVP (8a, b) to the solution u(t) of the problem (la, b) is given. THEOREM 1. Let u(t) and u(t,¢3) are the solutions of the problems (la, 5) and (8a, b), respectively. Then, for all t in [0, 1], we have
(12)
lu(t) - u(t, ¢3)t - c~.
PROOF. The proof of this theorem follows from the stability result (Lemma 2). From (la, b) and (8a, b) we have L[u(t) - u(t, ¢3)] : Lu(t) - Lu(t, 63) = f(t) - f(t) Bo[u(0)
= O,
- ~t(0, ¢ 3 ) ] ~-~ S o u ( O ) - S o u ( O , ¢ 3 ) m_ ¢1 - ¢ 1 = 0,
32
K. S E L V A K U M A R
and Bl[U(1) - u(1, ¢3)] - Blu(1) - Bltt(1, ¢3) = ¢2 - BlU(1, Ca) = ¢2 -- [U(1, ¢ 3 ) -t- eu'(1, ¢3)1.
From (9) we have IB1[u(1) - u(1, ~b3)]l _< 1¢2 - [u(1, ¢3) -~- eu'(1, ¢3)]I <_Ce. Now, using stability result, we have lu(t) - u(t, ¢3)1 --- Ce.
REMARK. W h e n e goes to zero, the solution u(t, ¢3) of (83, b) converges to the solution u(t) of (la, b). T h e following t h e o r e m gives the main result of this section, t h a t is, the convergence of the c o m p u t a t i o n a l m e t h o d . THEOREM 2. Let u(t) and Hi be the solutions of (la, b) and (lla~c), respectively. Then, we have, tu(td - u~l _< C(e + h).
PROOF.
Prom [3] we have the estimate for the solution of the scheme
(1 la-c) as
lu(t~, ¢3) - u~l <_ Ch. Prom T h e o r e m 1 we have, pu(t~) - u(t,, ~3)1 --- c e .
Therefore, I~(td - u~l < lu(td - u(t~, ¢3)1 + Ju(t~, ¢3) - ~ l <_ Ce + Ch = C(e + h).
33
Chemical Flow Reactor
REMARK. When the scheme (6a-d) is applied to the IVP (10a-c) the computational time will be reduced at each nodal points. From [3] we have lu(t~, ¢ 3 ) -
~1 < Ch,
where u(ti, ¢3) is the solution of the IVP (Sa, b) and us is the solution of the scheme (6a-d). Therefore, using the estimate (12), lu(t~) - u~l <_ l u ( t j
u(t~, ¢z)l
-
+
lu(t~, ¢3)
- u,I
< C(c + h).
4,
T E S T EXAMPLES AND NUMERICAL RESULTS EXAMPLE 1.
Consider the following homogeneous problem ~u" + u' = 0,
0
- u ' ( 0 ) = ¢1,
(13a)
u(1) + cu'(1) = ¢2,
(135)
where ¢1 = 0 and ¢2 = 0. The initial guess ¢3 for u(0) is given by u(0) = ¢2 +e¢1 from the zeroth order asymptotic expansion u(t,
~) = ¢2 + c ¢ l e x p ( - t / ~ ) .
The corresponding IVP is e u Ir + u ~ = 0,
(14a)
O< t < 1
u(0) = ¢3 = ¢2 + ~¢1,
~'(0) = - ¢ 1 ,
(14b)
whose solution is u(t,
¢3) = (¢3 - 6¢1) + c ¢ l e x p ( - t / ~ )
such that u(1, ¢3) + ~u'(1, ¢3) - ¢2 = 0.
The IVP (14a-b) can be reduced into a system of the form u~ - u2 = 0,
cu~ + u2 = 0,
ul(0)
= ¢3
u2(0) = -¢1,
(15~)
(15b)
34
K. SELVAKUMAR
where ul = u and u2 = u ~. The IVP (15a, b) can be solved numerically using the schemes given in Section 2. It is observed that the scheme (4a-d) reduces to the scheme (6a-d) when applied to the IVP (13a, b). Numerical results are given in Table 2. It is noted that this computational procedure solves exactly the sample problem (13a, b). EXAMPLE 2. able coefficient
Consider the following homogeneous problem with vari-
: u " + (1 + t)u' = 0,
0 < t < 1,
-ut(0) ~- (/)1,
(16a)
U(1) + eU'(1) = ¢2,
(16b)
where ¢1 = 0 and ¢2 = 1.5. The initial guess ¢3 for u(0) is given by u(0) = ¢2 + e e l , from the zeroth order asymptotic expansion
u(t, :) = ¢2 + : ¢ : e x p ( - t / e ) . The corresponding IVP is : u " + (1 + t)u' = 0,
0 < t < 1,
(17a)
u(0) = ¢3 = ¢2 + :¢1,
u'(0) = - ¢ 1 ,
(17b)
whose exact solution is
u(t, ¢3)
---- (¢3 -- : ¢ 1 ) q- : ¢ l e x p ( - - ( t
q-
t2/2)/C)
such that
u(1, Ca) + Eu'(1, Ca) - ¢2 = - : ¢ i e x p ( - 3 / 2 ~ ) = O(:). The IVP (17a, b) can be reduced into a system of the form u~ - u2 = 0, :U~ -~- (1 -[- t)~t2 ~-- 0,
u,(0) = ¢3
(18a)
~2(0) ---- --(~1,
(lSb)
where ul = u and u2 = u ~. Now the IVP (18a, b) can be solved numerically using the schemes (4a-d) and (6a-d). Table 3 gives the numerical result using the schemes (4a-d) and (6a-d), respectively, for the IVP (18a, b). EXAMPLE 3.
Consider the following nonhomogeneous problem
:u jt + u j = 1 + 2t,
0< t < 1
(19a)
Chemical Flow Reactor
35
TABLE 2* *
u~
u(t~)
u ( t i ) - u~
0 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + O0
0 . 0 0 0 0 0 E + 00
6 . 2 5 0 0 0 E - 02
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + O0
0 . 0 0 0 0 0 E + 00
1 . 2 5 0 0 0 E - 01
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00
0 . 0 0 0 0 0 E + 00
1 . 8 7 5 0 0 E - 01
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00
0 . 0 0 0 0 0 E + 00
2 . 5 0 0 0 0 E - O1
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00
0.0000OE + 0 0
3 . 1 2 5 0 0 E - 01
2 . 0 0 0 0 0 E + 00
2.O0000E + O0
0.0O000E + 00
3 . 7 5 0 0 0 E - 01
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00
0 . 0 0 0 0 0 E + 0O
4 . 3 7 5 0 0 E - 01
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00
0.0O000E + 0O
5.0000OE - 01
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00
0 . 0 0 0 0 0 E + 00
5 . 6 2 5 0 0 E - 01
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00
0 . 0 0 0 0 0 E + 00
6 . 2 5 0 0 0 E - 01 6 . 8 7 5 0 0 E - 01
2 . 0 0 0 0 0 E + 00 2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00 2 . 0 0 0 0 0 E + 00
0.0O000E + 00 0 . 0 0 0 0 0 E + 00
7 . 5 0 0 0 0 E - 01
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00
0 . 0 0 0 0 0 E + 00
8 . 1 2 5 0 0 E - Ol
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + O0
0 . 0 0 0 0 0 E + 00
8 . 7 5 0 0 0 E - 01
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00
0 . 0 0 0 0 0 E + 00
9 . 3 7 5 0 0 E - 01
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00
0 . 0 0 0 0 0 E + 00
1 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + O0
0 . 0 0 0 0 0 E + 00
*h = ¢ = 1/16, = O.O0000E + O0
¢3 = 2 . 0 0 0 0 0 E + 0 0 ,
u ( 1 , ¢ 3 ) + ¢ u ' ( 1 , ¢ 3 ) - ¢2
t
u~
u(t~)
u ( ~ ) - u~
0.0000OE + O0
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + O0
0 . 0 0 0 0 0 E + 00
6 . 2 5 0 0 0 E - 02
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00
0 . 0 0 0 0 0 E + 00
1 . 2 5 0 0 0 E - 01
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00
0 . 0 0 0 0 0 E + 00
1 . 8 7 5 0 0 E - 01
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00
0.00000E + 00
2 . 5 0 0 0 0 E - 01
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00
0 . 0 0 0 0 0 E + 00
3 . 1 2 5 0 0 E - 01
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00
0 . 0 0 0 0 0 E + 00
3 . 7 5 0 0 0 E - 01
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00
0 . 0 0 0 0 0 E + 00
4 . 3 7 5 0 0 E - 01
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00
0.00000E + 00
5 . 0 0 0 0 0 E - 01
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00
0 . 0 0 0 0 0 E + 00
5 . 6 2 5 0 0 E - 01
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00
0 . 0 0 0 0 0 E + 00
6 . 2 5 0 0 0 E - 01
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + O0
0 . 0 0 0 0 0 E + 00
6 . 8 7 5 0 0 E - 01
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00
0 . 0 0 0 0 0 E + 00
7 . 5 0 0 0 0 E - 01
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + O0
0.0000OE + 00
8 . 1 2 5 0 0 E - 01
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + O0
0 . 0 0 0 0 0 E + 00
8 . 7 5 0 0 0 E - 01
2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00
0 . 0 0 0 0 0 E + 00
9 . 3 7 5 0 0 E - 01 1 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00 2 . 0 0 0 0 0 E + 00
2 . 0 0 0 0 0 E + 00 2 . 0 0 0 0 0 E + 00
0 . 0 0 0 0 0 E + 00 0 . 0 0 0 0 0 E + 00
*h = 1/16, e = 1.00000E-05, + ~ u ' ( 1 , ¢ 3 ) - ¢2 = O.00000E + 00
¢3 =
2.0000OE+00,
u(1,¢a)
36
K. S E L V A K U M A R
T A B L E 3*
0.00000E + 00
1.50000E + 00
1. 50000E + 00
0.00000E + 00
6. 25000E - 02
1.50000E + 00
1.50000E + 00
0.00000E + 00
1. 25000E -- 01
1.50000E q- 00
1.50000E -I- 00
0.00000E q- 00
1.87500E - 01
1.50000E + 00
1.50080E + 00
0.00000E + 00 0.00000E + 00
2.50000E - 01
1.50000E + 00
1.50000E + 00
3.12500E - 01
1.50000E + 00
1.50000E + 00
0.00000E + 00
3. 75000E - 01
1.50000E + 00
1,50000E + 00
0.00000E + 00
4 . 3 7 5 0 0 E - 01
1.50000E + 00
1.50000E + 00
0.00000E + 00
5 , 0 0 0 0 0 E - 01
1.50000E + 00
1,50000E + 00
0.00000E + 00
5.62500E - 01
1.50000E + 00
1.50000E + 00
0.00000E + 00
6.25000E - 01
1.50000E + 00
1.50000E ÷ 00
0.00000E + 00
6.87500E - 01
1.50000E + 00
1,50000E + 00
0.00000E + 00
7.50000E - 01
1.50000E + 00
1.50000E + 00
0.00000E + 00
8.12500E - 01
1.50000E + 00
1.50080E + 80
0.00000E + 00
8.75000E - 01
1.50000E + 00
1.50000E + 00
0.00000E + 00
9.37500E - 01
1.50000E + 00
1.50000E + 00
0.00000E + 00
1.00000E -b 00
1.50000E + 00
1.50000E + 00
0.00000E -t- 00
*h -- ~ --= 1/18, = 0.00000E + O0
¢3
=
1.50000E+00,
u(1, q~3) + Eu'(1,~b3) -- ¢2
t
ui
u(t~)
u ( t i ) - u~
0.00000E + 00
1.50000E + 00
1.50000E + 00
0.00000E + 00
6.25000E - 02
1.50000E + 00
1,50000E + 00
0,00000E + 00
1.25000E - 01
1.50000E + 00
1.50000E + 00
0.00000E + 00
1.87500E - 01
1.50000E + 00
1.50000E + 00
0.00000E -b 00
2.50000E - 01
1.50000E + 00
1.50000E + 00
0.00000E + 00
3.12500E - 01
1.50000E + 00
1.50000E -+ 00
0.00000E + 00
3.75000E - 01
1.50000E + 00
1.50000E + 00
0.00000E q- 00
4,37500E - 01
1.50000E + 00
1.50000E + 00
0.00000E + 00
5.00000E - 01
1.50000E -+-00
1.50000E + 00
0.00000E + 00
5,62500E - 01
1.50000E q- 00
1.50000E + 00
0.00000E q- 00
6.25000E - 01
1.50000E + 00
1.50000E + 00
0.00000E + 00
6.87500E - 01
1.50000E + 00
1.50000E + 00
0.00000E + 00
7.50000E - 01
1.50000E + 00
1.50000E + 00
0.00000E + 00
8.12500E - 01
1.50000E + 00
1.50000E + 00
0.00000E + 00
8.75000E - 01
1.50000E + 00
1.50000E + 00
0.00000E + 00
9.37500E - 01
1.50000E -b 00
1.50000E -b 00
0.00000E -b 00
1.00000E -t- 00
1.50000E q- 00
1.50000E + 00
0.00000E -k 00
* h = 1/16, c = 1.00000E-05, + ~u' ( 1, ¢3) - ¢2 = 0.00000E + 00
q~3 =
1.50000E+00,
u(1,¢3)
37
Chemical Flow Reactor
where
¢1
-'~
(tgb)
u(1) + 6u'(1) ----¢2,
=
0 and ¢2 = 3.
The initial guess ~b3 for u(0) is given by u(0) = ¢2 - 2 + ¢(¢1 the zeroth order asymptotic expansion U ( t , g ) = (~b2 -- 2 ) + t + t 2 4- g ( ¢ l
--
--
1) from
1)exp(-t/~).
The corresponding I V P is ~u" + u / = 1 + 2t,
0 < t < 1
u(0) = ¢3 ----¢2 - 2 + e ( ¢ 1 - 1),
(20a) u'(0) = - ¢ 1 ,
(20b)
whose solution is u ( t , ¢3) = ¢3 - ¢(¢z + 1 - 26)(1 - e x p ( - t / ¢ )
+ t + t 2 - 2~t
such that u(1, ¢3) + 6u'(1, ¢3) - ¢2 = - 6 , which is O(6). The IVP (20a, b) can be reduced into a system of the form ! U 1 - - U 2 = 0,
¢u~ + u2 = 1 + 2t,
Ul(0) = ¢3 u2(0) = - ¢ 1 ,
(21a)
(21b)
w h e r e u l = u a n d u2 = u t. Now the IVP (21a, b) can be solved by the schemes given in Section 2. It is observed that the scheme (4a-d) reduces to the scheme (6a-d) when applied to the IVP (21a, b). The numerical results are given in Table 4. A similar discussion can be carried out for the following test problem CU II 4- 'tt/ - - I t ~ O,
--u'(O) =
0 < t < t
u(1) + ¢u'(1) : ¢2,
(22a) (22b)
where ¢1 = 0 and ~b2 = 1. The numerical results are given in Table 5.
5.
CONCLUSION
From the above discussion, we conclude that the given chemical flow reactor problem can be solved by using shooting method. The given chemical reactor flow problem is transformed into an initial value problem using shooting method with an initial guess from zeroth order asymptotic
38
K. SELVAKUMAR
TABLE
t
*h =
u(t~)
u ( t D - ui
0 . 0 0 0 0 0 E 4- 0 0
9.37500E - 01
1.00000E4-
00
6.25000E
6.25000E
- 02
9.41406E - 01
1.02402E4-
00
8.26185E
- 02
1.25000E - 01
9.89758E -- 01
1 . 0 7 7 7 1 E 4- 0 0
8.79552E
- 02
- 02
1.87500E - 01
1.05940E + 00
1.14725E + 00
8.78540E - 02
2.50000E - 01
1.14181E + 00
1.22756E + 00
8.57525E - 02
3.12500E -- 01
1 . 2 3 3 8 6 E 4- 0 0
1 . 3 1 6 7 7 E 4- 0 0
8.29148E -- 02
3.75000E - 01
1 . 3 3 4 3 9 E 4- 0 0
1 . 4 1 4 2 0 E 4- 0 0
7.98066E - 02
4.37500E - 01
1.44298E + 00
1.51958E + 00
7.65988E - 02
5.00000E - 01
1 . 5 5 9 4 8 E 4- 0 0
1.63283E + 00
7.33542E - 02
5.62500E - 01
1 . 6 8 3 8 2 E 4- 0 0
L75391E4-00
6.25000E
- 01
1.81598E÷00
1.88282E+
6.87500E
- 01
1 . 9 5 5 9 6 E 4- 0 0
2 . 0 1 9 5 3 E 4- 0 0
6.35684E -- 02
7.50000E - 01
2 . 1 0 3 7 6 E 4- 0 0
2 . 1 6 4 0 6 E 4- 0 0
6.03030E - 02
8.12500E - 01
2.25937E + 00
2.31641E + 00
5.70374E - 02
8.75000E
- 01
2,42279E + 00
2.47656E4-
5.37715E - 02
9.37500E - 01
2.59403E + 00
2 . 6 4 4 5 3 E 4- 0 0
5.05056E - 02
1.00000E4-
2 . 7 7 3 0 7 E 4- 0 0
2 . 8 2 0 3 1 E 4- 0 0
4.72398E - 02
=
e
=
6.25000E
00 1/16,
¢3
=
9.37500E-01,
7.00963E - 02 00
u(1,~3)
00
6.68333E
4- e u ' ( 1 , ¢ 3 )
-
- 02
~2
-- 02
t
*h
4*
u~
~
~(t~)
~(t~) - ~i
0.00000E + 00
9.99990E - 01
1.00000E + 00
1.00136E - 05
6.25000E
- 02
1.00390E + 00
1.06640E + 00
6.24988E
1.25000E - 01
1 . 0 7 8 1 1 E 4- 0 0
1 . 1 4 0 6 1 E 4- 0 0
6.24975E - 02
- 02
1.87500E - 01
1.16015E + 00
1.22264E + 00
6.24963E - 02
2.50000E - 01
1.24999E + 00
1.31248E + 00
6.24950E - 02
3.12500E - 01
1 . 3 4 7 6 5 E 4- 0 0
1 . 4 1 0 1 4 E 4- 0 0
6.24938E - 02
3.75000E - 01
1.45311E + 00
1.51561E + 00
6.24925E - 02
4.37500E
- 01
1.56640E + 00
1.62889E + 00
6.24913E
- 02
5.00000E - 01
1 . 6 8 7 4 9 E 4- 0 0
1.74998E + 00
6.24900E
- 02
5.62500E - 01
1.81640E + 00
1.87889E + 00
6.24888E
- 02
6.25000E
- 01
1.95311E + 00
2.01560E + 00
6.24876E -- 02
6.87500E - 01
2.09765E + 00
2.16013E + 00
6.24864E -- 02
7.50000E - 01
2 . 2 4 9 9 9 E 4- 0 0
2.31247E + 00
6.24850E
8.12500E
- 01
2.41015E + 00
2.47263E + 00
6.24838E - 02
8.75000E - 01
2.57811E + 00
2.64060E + 00
6.24826E - 02
9.37500E - 01
2.75390E + 00
2 . 8 1 6 3 8 E 4- 0 0
6.24814E - 02
1.00000E + 00
2.93749E + 00
Z99997E
6.24800E
=
1/16,
÷ eu'(1,¢3)
e
=
1.00000E--05,
-~b2 = 1.00136E-05
¢3
=
4- 0 0
9.99990E+01,
u(1,¢3)
- 02
- 02
Chemical Flow Reactor
39
TABLE
t
5*
ul
u(tl)
u(tl) - ui
0.00000E
+00
3.90872E
-
01
3,87755E
- 01
3.11729E
- 03
6.25000E
- 02
3.91652E
- 01
3.96713E
- 01
5.06178E
- 03
1.25000E
-
01
4.08920E
- 01
4A5763E
-
02
6.84297E
- 03
1.87500E
-
01
4.32384E
- 01
4,39281E
-
01
6.89682E
- 03
2.50000E
-
01
4.58915E
- 01
4.65376E
-
01
6.46129E
- 03
3.12500E
-
01
4.87614E
- 01
4.93457E
-
01
5.84248E
- 03
3.75000E
- 01
5.18279E
- 01
5.23382E
-
01
5.10353E
- 03
4.37500E
- 01
5.50928E
- 01
5.55175E
- 01
4.24737E
- 03
5.00000E
- 01
5.85654E
-
01
5.88917E
- 01
3.26347E
- 03
5.62500E
- 01
6.22578E
-
01
6.24716E
- 01
2.13802E
- 03
6.25000E
- 01
6.61837E
- 01
6.62694E
- 01
8.56698E
- 04
6.87500E
- 01
7.03577E
-
01
7.02981E
- 01
5.96464E
- 04.
7.50000E
- 01
7.47955E
-
01
7.45717E
- 01
2.23833E
-- 03
8.12500E
- 01
7.95138E
-
01
7.91051E
- 01
4.08697E
- 03
8.75000E
- 01
8.45304E
- 01
8.39142E
- 01
6.16199E
- 03
9.37500E
-
01
8.98641E
- 01
8.90156E
- 01
8.48490E
- 03
1.00000E
4. 00
9.55351E
- 01
9.44272E
- 01
1.10788E
- 02
* h
=
=
•
e
8.03924E
1/16,
¢3
=
3.90872E+00,
u(1,¢3)
4- ~ u ' ( 1 , ¢ 3 )
-
¢2
- 03
t
u~
u(ti)
u(~i) - u~
0.00000E
+ 00
3.67916E
- 01
3.68038E
- 01
1.21415E
6,25000E
- 02
3.68650E
-
01
3.91724E
-
01
2.30743E
- 04 - 02
1.25000E
- 01
3.93909E
- 01
4.16977E
- 01
2.30682E
- 02
1.87500E
- 01
4.19360E
- 01
4.43858E
- 01
2.44982E
- 02
2.50000E
- 01
4.46545E
- 01
4.72472E
- 01
2.59270E
- 02
3.12500E
- 01
4.75486E
-
01
5.02930E
- 01
2.74442E
- 02
3.75000E
- 01
5.06302E
-
01
5.35352E
-
01
2.90496E
- 02
4.37500E
- 01
5.39115E
-
01
5.69864E
- 01
3.07490E
- 02
5.00000E
- 01
5.74053E
-
10
6.06600E
- 01
3.25474E
- 02
5.62500E
-
01
6.11254E
- 01
6.45705E
- 01
3.44509E
- 02
6.25000E
- 01
6.50866E
-
01
6.87331E
-
01
3.64656E
- 02
6.87500E
- 01
6.93043E
- 01
7.31641E
- 01
3.85979E
- 02
7.50000E
- 01
7.37925E
-
01
7.78807E
- 01
4.08546E
- 02
8.12500E
- 01
7.85770E
- 01
8.29013E
- 01
4.32430E
- 02
8.75600E
- O1
8.36685E
-
01
8.82456E
- 01
4.57709E
- 02
9.37500E
- O1
8.90898E
-
01
9.39344E
-
01
4.84463E
~ 02
1.00000E
4. 00
9.48622E
- 01
9.99900E
-
01
5.12775E
- 02
*h +
=
1/16,
eu'(1,¢a)
e -
=
1.00000E-04,
¢2 = 3,30031E-04
¢3
=
3.67916E+01,
u(1,¢a)
40
K. SELVAKUMAR
expansion for the solution. Then the resultant initial value problem is solved numerically using exponentially fitted finite difference schemes. The similarity between the given chemical reactor flow problem and the initial value problem is that the nonuniformity occurs only in the derivative but not in the solution. The advantage of this method is that there is no need for matrix inversion as in tridiagonal difference schemes and a three-point difference scheme is replaced by a two-point difference scheme (one-step method) using shooting method. When the exact solution of the initial value problem (Sa, b) is not known, using the numerical solution of (8a, b) one can check the condition (9). In Tables 2-5 the condition (9) is checked and given for the respective test problems to illustrate the procedure.
REFERENCES 1 J. V. Baxley, On singular perturbation of nonlinear two-point boundary value problems, Czechoslovak Math. J. 27:363-377 (1977). 2 E. P. Doolan, J. J. H. Miller, and W. H. A. Schilders, Uniform Numerical Methods for Problems with Initial and Boundary Layers, Boole Press, Dublin, 1980. 3 K. Selvakumar, Uniformly Convergent Difference Schemes for Differential Equations with a Parameter, Ph.D. Thesis, Bharathidasan University, India, 1992. 4 J.D. Lambert, Computational Methods in Ordinary Differential Equations, Wiley, 1973.