On Newman's Numerical Technique for Solving Boundary Value

solution technique for solving two-point boundary value problems which consist of coupled, ordinary differential equations. His technique has not been...
0 downloads 0 Views 217KB Size
Ind. Eng. Chem. Fundam., Vol. 17, No. 4, 1978

367

On Newman’s Numerical Technique for Solving Boundary Value Problems

Ind. Eng. Chem. Fund. 1978.17:367-369. Downloaded from pubs.acs.org by UNIV OF WEST LONDON on 10/29/18. For personal use only.

Physical phenomena that interest chemical engineers can often be modeled mathematically as boundary value problems which consist of sets of nonlinear, coupled, second-order differential equations. Newman’s (1968) technique for solving such problems is reviewed and extended slightly.

Ten years ago, Newman (1968) presented a numerical solution technique for solving two-point boundary value problems which consist of coupled, ordinary differential equations. His technique has not been used widely perhaps because his derivation is not transparent. The purpose of this paper is to present an alternate derivation of the equations programmed by Newman (1973a) in an attempt to promote the use of his technique, to present an alternate method of treating boundary conditions, and to present (in Newman’s notation) an equation for linearizing nonlinear, finite difference equations. A linear set of coupled, second-order, ordinary differential equations can be written as

,()

d2ck{x)

+ bitk(x)

dx2

dcfc(x)

B(l) A(2) B(2) D(2) A (j) B (;) DO') B (NJ) '



"

C(l)

G(l)

C(2)

G(2)

CO')

GO)

=

+ di¡k(x)ck(x)

dx

G(NJ)_

C{NJ)_ or as

where i represents the equation number and varies from 1 to n for a complete set of equations. Using central difference expressions to approximate the derivatives, these equations become k=\

1) + Bl¡k(j)C,k(j) +

Ai,k(j)Ck(j

Dij,(j)Ck(j + where

Umax

=



=

Biik(j)

=

=

for all

Cfc(xy)

j

for

,*(*,)

-

2

=

for

-2aiik(Xj) + h2di k(xj)

j

j

(3)

DUkO')

LU T

=

Z

=

(13)

to NJ =

to NJ

2

-

aitk(xj) +

for

|h^(x;)

j

2

=

to NJ

-

1

b(l)

(6)

1

(14)

(15) L£ = Z The unknown T can be obtained by solving eq 14 once L, U, and £ are known. Using Newman’s notation, the matrix equations that result from decomposing M are

(4)

1

-

U T=£

EO)

=

B(l)

=

for

-b_10) DO) aO)

and

for

AO)

=

j

(16)

to NJ

;'

=

1

=

2

to NJ

-

(17)

1

(18)

and

G,0)

Let

(12)

becomes

(5) =

Z

which by letting

NJ) Ck(j)

, 0')

T

G,·(/') (2)

=

1)

MT=

X n matrices (B(l), etc.) where the elements of M are and the elements of T and Z are n-dimensional vectors. Equation 12 can be solved by decomposing M into lower X n and upper triangular matrices L and U (with identity matrices along the diagonal)

us

=

for j

h2gl(Xj)

=

2

to NJ

-

consider first simple boundary conditions

of the form

*=i

which for mesh point

euk{x)ck 1 can

be written

k

=

BO) + aO) EO

£0)

(9)

as

b_10) [GO')

=

GO')

T.BUk(NJ)Ck(NJ)

=

Gi(NJ)

These governing difference equations (eq 2, 9, and 10) be written conveniently in matrix form as

j

2

=

to NJ

(19)

=

£0) + E0)

-

=

b-HDGd)

a0)£0 C{NJ) GO'

1)] for

-

=

+1)

(20)

j

=

to NJ

2

(21) (22)

£(2VJ)

for

j

=

NJ

-

1

to

1

(23)

These (and the following) governing matrix equations can be solved with a high-speed digital computer using Newman’s (1973a) subroutine BAND(J). Now let us consider general boundary conditions of the

(10) can

0019-7874/78/1017-0367$01.00/0

for

1)

-

£(1)

as

1

and for mesh point NJ k=l

=

where aO) and bO) are the nonzero elements of L and -EO) are the off-diagonal, nonzero elements of U. Once the decomposition equations have been solved, forward and backward substitution can be used to solve eq 15 and 14, respectively

(8)

G.U)

=

bO)

on ck

/,(x)

=

EBiJk(l)Cfc(l)

(7)

1

form

©

1978 American Chemical Society

368

Ind. Eng. Chem. Fundam., Vol. 17, No. 4, 1978 n

dc*

0) ™ dx

fc=i

+ eiJk(x)ck

=

(24)

/¿(x)

To use central difference approximations for the derivatives in eq 24 Newman introduced the concept of an image point. That is, he let the points j = 1 and j = NJ be image points and the physical boundaries be at j = 2 and j NJ 1. Then, for j = 1 eq 24 becomes =

k=l

-

,*( ) ( )

+

a,*(DA(2)

+

xiikckm

=

G.d)

(25)

where Pi,k(x)

Bhk( 1)=

Xi,k

A>(D

=

(26)

=

2

(That is, eq 36-39 replace those given earlier by eq 17,18, 21, and 23.)

Let us now consider a special case. If the boundary conditions are such that ,· ( ) = 0 for some but not all i, then special treatment is required. One way to proceed is to use image points as needed and then eliminate explicit dependence of the system of equations on the values of the unknowns at the image points [see White (1977) for an example of this technique]. Another way is to use the following approximation for the derivatives (instead of image points and central differences) -ch(xj + 2h) + 4ck(xj + h)

dct dx

2h

Substituting eq 40 and

(27)

he¡_*(x)

Ad) and for

j

=

hfi(x)

(28)

2) +

Aiik(NJ)Ck(NJ

-

where Yij,

~Bi,k(NJ)

(30)

-

AUk(NJ)

=

heik(x)

(31)

Gi(NJ)

=

hfi(x)

(32)

and

Now the governing difference equations in matrix form are as shown in eq 33. Proceeding as before yields the

B(l) D(l) X

=

G;(l)

=

-n*

A>( 1)

=

DO)

Y

A (NJ)

ZAi:k(j)Ck(j

-

1) +

C(l)

G(l)

C(2)

G(2)

=

C(j)

GO)

(45)

2 pitk(x)

1)

=

G¡(j) (46)

dF;

_

_

t>Ck(j

;

Buk(j)

=

Di,k(j)

for l

G(NJ)m

(34)

Y

(35)

=

[D(2) + a(2) x]

—b_1(2) =

A (NJ) + yE(IVJ

-

(36) (37)

2)

=

yí(NJ

-

2)

-

a

(NJ)

t(NJ

-

1)] (38)

and

(1) + E(1)A2) + xA3)

(47)

=

-

~

1) C*(i)=C4°(i)

and

1,

j

+

1

(39)

=

9Ck{j)

j

-

(48) Ch(l)=Ck°(l)

and

1,

j

+

1

dF:

(33)

-B_1(1)X y

=

=

~Ai,k(NJ)

=

for l

=

=

(49)

dCk(j + 1) c*(0=c*°(i) j 1, j, and ; + 1 -

and

AO') x

AD

(44)

dF,

substitutions

-

hfi(x)

=

Biik(j)Ck(j) +

for l

B (NJl

previous matrix equations with the following additions and

£(NJ) b~1(NJ)[G(NJ)

(43)



Finally, to solve a set of nonlinear, finite difference equations each of which is expressed as F¡ = 0, one can expand F¡ about trial values and proceed as before. In such a case, each of the governing finite difference equations is of the form

0')

BO)

C(NJ)

(NJ)

(42)

Pi,k(x)

=

G¡(NJ)

'

a

+ hei:k(x)

|,()

where

AO)

=

(41)

Di,k(j)Ck(j +

A(2) B(2) D(2)

E(2)

+ heiM{x)

and

Pi,k(x)

=

a,

NJ into

with

1) +

Bhk(NJ)Ck(NJ)= Gi(NJ) (29)

=

=

Bi'k(NJ)

+ 0(h2) (40)

similar equation for j

3Pi,k(x)

=---—

Bi:k( 1)

NJ becomes

Yi,kCk(NJ-

a

eq 24 yields eq 25 and 29

and =

3ck(x,)

-

=

~F¡° + EAlM(j)Ck°(j

-

1) +

BiM(j)Ck°(j) +

) °0

+ 1) (50)

where the superscript o designates a trial value. If nearly correct trial values are selected, iteration will yield values accurate to within a specified tolerance limit. Newman (1973a) enumerates boundary value problems which can be solved using his technique and discusses briefly alternate solution techniques—see also White et al. (1975), Roberts and Shipman (1972), Liu (1967), Viliadsen and Michelsen (1978), and von Rosenberg (1969). In addition, eigenvalue problems can be solved using Newman’s technique—see, e.g., Newman (1973b), Fedkiw

Ind. Eng. Chem. Fundam., Vol. 17, No. 4, 1978

and Newman (1977), and Smyrl and Newman (1974). Recent applications of Newman’s technique have been reported by Trainham and Newman (1978, 1977), Alkire and Ng (1977), Alkire and Varjian (1977), White et al. (1977), Fedkiw and Newman (1977), Viswanathan and Chin (1977), and Choi et al. (1976).

Acknowledgment This work was supported by the Center for Energy and Mineral Resources at Texas A&M University. Professor John Newman’s comments and suggestions are acknowledged gratefully. Literature Cited Alkire, R., Ng, P. K„ J. Electrochem. Soc., 124, 1220 (1977). Alkire, R„ Varjian, R., J. Electrochem. Soc., 124, 388 (1977). Choi, K. W., Bennlon, D. N., Newman, J., J. Electrochem. Soc., 123, 1616 (1976). Fedkiw, P., Newman, J., AIChE J., 23, 255 (1977). Liu, S„ Chem. Eng. Sc/., 22, 871 (1967). Newman, J. S., “Electrochemical Systems", Prentice-Hall, Englewood Cliffs, N.J., 1973a.

369

Newman, J., Electroanal. Chem., 6, 187 (1973b). Newman, J., Ind. Eng. Chem. Fundam., 7, 514 (1968). Roberts, S. M., Shipman, J. S„ “Two-Point Boundary Value Problems: Shooting Methods", American Elsevier, New York, N.Y., 1972. Smyrl, W. H„ Newman, J„ J. Electrochem. Soc., 121, 1000 (1974). Trainham, J. A., Newman, J., J. Electrochem. Soc., 124, 1528 (1977). Trainham, J. A., Newman, J., J. Electrochem. Soc., 124, 1220 (1977). Vllladsen, J., Mlchelsen, M. L., "Solution of Differential Equation Models by Polynomial Approximation", Prentice-Hall, Englewood Cliffs, N.J., 1978. Viswanathan, K., Chin, D-T., J. Electrochem. Soc., 124, 709 (1977). von Rosenberg, D. U., “Methods for the Numerical Solution of Partial Differential Equations", American Elsevier, New York, N.Y., 1969. White, R. E., Ph.D. Thesis, Department of Chemical Engineering, University of California, Berkeley, Calif., 1977. White, R., Trainham, J. A., Newman, J„ Chapman, T, W„ J. Electrochem. Soc., 124, 669 (1977).

Department of Chemical Engineering Texas A&M University College Station, Texas 77843

Ralph E. White

Received for review January 23, 1978 Accepted June 2, 1978