The Chemistry Maths Book, Second Edition

20.10 Systems of differential equations 591

general problem is that of Nfirst-order initial value problems for Nfunctions,y

1

(x),

y

2

(x), =, y

N

(x):

(20.59)

with given initial valuesy

1

(x

0

),y

2

(x

0

), =

Such systems of first-order initial value problems occur in the chemical kinetics

of multistep chemical reactions, and are therefore of considerable interest to the

chemist (boundary value problems are generally much more difficult). They are

solved numerically by applying one of the methods described in Section 20.9 to each

equation in turn at each step. For example, for the pair of equations

y′(x) 1 = 1 f

1

(x, y, z), z′(x) 1 = 1 f

2

(x, y, z) (20.60)

with initial conditionsy(x

0

) 1 = 1 y

0

andz(x

0

) 1 = 1 z

0

, Euler’s method leads to the sequence

of recursion steps

(1) y

1

1 = 1 y

0

1 + 1 hf

1

(x

0

, y

0

, z

0

), z

1

1 = 1 z

0

1 + 1 hf

2

(x

0

, y

0

, z

0

)

(2) y

2

1 = 1 y

1

1 + 1 hf

1

(x

1

, y

1

, z

1

), z

2

1 = 1 z

1

1 + 1 hf

2

(x

1

, y

1

, z

1

)

(3) y

3

1 = 1 y

2

1 + 1 hf

1

(x

2

, y

2

, z

2

), and so on

The procedure is essentially the same as that for the single equation.

Stiff equations

A stiff differential equation is one whose solution contains terms that differ greatly in

their dependence on the independent variable. For example, the second-order initial

value problem

(20.61)

has the general solutiony 1 = 1 ae

− 10 x

1 + 1 be

+ 10 x

and particular solutiony 1 = 1 e

− 10 x

. A numerical

solution of the equation by one of the methods discussed so far gives a solution that

behaves correctly for values of xclose to x 1 = 10 but ‘explodes’ like e

+ 10 x

as xincreases

because of rounding and truncation errors that inevitably lead to a small admixture of

the unwanted term; that is, the numerical solution has the form

y 1 ≈ 1 e

− 10 x

1 + 1 εe

+ 10 x

(20.62)

where εis not zero.

Problems involving stiff equations occur in kinetics when several elementary

processes have very different rate constants. Thus, problem (20.61) is equivalent to

the pair of first-order problems

(20.63)

dy

dx

z

dz

dx

=, = 100 y

dy

dx

yy y

2

2

=, 100 ()0 1 0 10=,′=−()

yx′ == ,,,,, =,,,

dy

dx

fxyy y i N

i

i

iN

() ( )

12

...... 12