Modeling of complex chemical reactions in a continuous-flow reactor: A Markov chain approach

Chemical Engnwcrtng Science, Printed in Great Britain.

Vol. 43, No. IO, pp. 2807- 2815, 1988.

0009 2509/88 IF3.00+0.00 C 1988 Pergamon Press plc

MODELING OF COMPLEX CHEMICAL REACTIONS IN A CONTINUOUS-FLOW REACTOR: A MARKOV CHAIN APPROACH S. T. CHOU,’ Department

of Chemical

L. T. FAN*

Engineering, Durland

and R. NASSAR’ Hall, Kansas State University, U.S.A.

Manhattan,

KS 66506,

(Received 18 February 1987; accepted 17 March 1988) Abstract-To simulate the dynamics of a chemicaliy reacting system as a Markov chain, the states of the chain need to be properly identified. In the present approach, a molecule is viewed as an “object”, “entity” or system. The transformation of the molecule from one species to another is visualized as the transition of this “entity” from one form to another. Furthermore, appropriately selected chemical species or entities in the mixture serve as the states of the chain. The selection of such species is subject to the stoichiometric constraint based on the atomic balance. The approach is illustrated with several chemically reacting systems. The results are in good agreement with the known results obtained from the deterministic approach.

lNTRODUCTiON

approaches have been widely employed for simulating the dynamics of chemical reactions in a batch or flow reactor because of their probabilistic nature [e.g. Singer (1953), Bartholomay (1958a, b), Montroll and Shuler (1958), Isida (1960), McQuarrie (1963, 1964, 1967), Darvey and Staff (1966), Darvey et al. (1966), Fredrickson (1966), Oppenheim et al. (1968, 1977), Gillespie (1976, 1977a, b, 1981), Nassar et al. (198 1) and van Kampen (1981)]. In such approaches, a chemical reaction is considered as a Markov process with a continuous time parameter and countable states. Isothermal first-order reactions give rise to linear first-order differential equations which can be solved in a straightforward manner. Nevertheless, more often than not, higher-order or nonisothermal nonlinear governing differential reactions yield equations. Analytically solving such differential equations is extremely difficult, if not impossible. For this situation, it may be advantageous to resort to the use of a Markov chain model [e.g. Too ef al. (1983)]. It has been shown that the Markov chain is an effective tool for simulating dynamics of complex chemical reactions, particularly in a flow reactor [e.g. Miller (1972), Pippel and Philipp (1977a, b), Formosinho and Miguel (1979) and Too et al. (1983)]. To simulate the dynamics of a chemically reacting system as a Markov chain, the numbers of molecules of chemical species participating in the reactions may be identified as the states of the chain. Alternatively, appropriately selected chemical species participating in the reactions may be identified as the states of the chain. It is demonstrated here that the selection of these species needs to be subject to the stoichiometric constraint based on the atomic balance Stochastic

fDepartment of Statistics. tTo whom correspondence should be addressed. 2807

for nonlinear reacting systems. The present work represents, in essence, an extension of the approach proposed by Fredrickson (1966) to such systems. MODEL

FORMULATION

Consider a chemically reacting system containing a mixture of chemical species, Aj, j= 1, 2, . . . , J. A chemical reaction among Aj induces the change in the state of the mixture of the system. A typical reaction, characterized by a set of stoichiometric coefficients, sj and rj, can be written as [e.g. van Kampen (198l)j S,A, +s,J4,+

. . +s,/4,

&A,

+r,l4, + . . . +r,A,.

(I)

In what follows a molecule is considered as an “object”, “entity” or system [e.g. Fredrickson (1966) and van Kampen (1981)]; a molecule of type Aj is regarded as the “entity” at level j, and any reaction involving it is viewed as a transition of this “entity” from this level to another, or vice versa. Thus, the collection of all transitions involving the “entity” forms a Markov chain. Note that the levels serving as the states of the resultant Markov chain must be specified through judicious identification of the chemical species such that the system or “entity” can occupy one and only one state after each transition. Let the level of the “entity” or the state of the system at time M At be the random variable Y(m) and pij(m, mf 1) be the one-step transition probability that the ‘*entity” at level i at time m At will be at level j at time (m + 1) At inside the reactor, i.e.

pij(m,m+l)=Pr[

Y(m+l)=jl

Y(m)=i].

(2)

Furthermore, let Pia(m, m + 1) denote the transition probability that the “entity” at level i at time m At will be discharged to the surrounding environment (dead or absorbing state) from the reactor at time (m + 1) At.

2808

S. T. CHOU et al.

The one-step transition probabilities of this Markov chain are best exhibited collectively in the form of a matrix as shown below:

P(m, m-t

i

Ptr(m, m+ 1)

PIz(m, m+ 1)

.

.

P,,(m,

m+

1)

Pdm,

Pzl(m,

Pt2(m,

.

.

Pdm,

m+

1)

PaAm.

l)=

m+

1)

L

m+

1)

P,z(m,m+l)

.. . .. .

0

0

Pdm,m+l)

p ij( m,m+l)+pi,(m,m+l)=l, . . . , I;

m=O,

1, 2,. . .

(4)

Suppose that the absolute (or state) probabitity distribution of the “entity” among the selected I levels and the dead state after m-step transitions, pi(m), j= 1, 2, . . , 1, d, is denoted as a row vector, P*(m), i.e. P*(m)=

CpItm)pz(m)

Its initial distribution P*(O)=

. . . M4pdWl.

is denoted

CPLOP2(0)

Then the distribution

1)

1)

P,Am,m+l) 1

0

where I(!< J) is the total number of the independent levels serving as the states of the Markov chain. Note that all elements in this matrix are nonnegative, and the elements of each row sum to unity, i.e.

i= 1,2,

i1 m+

m+

(3)

I PII(m,m+l)

$,

which have come from the n,(m) “entities” at level i at time m At, due to the transition of these “entities” from level i to level j, and let the random variable N&n, m

(5)

as P*(O), i.e.

+ 1) represent the number of these q(m) “entities” discharged into the surrounding environment (dead state) from the reactor during the time interval [m At, (m + 1) At]. If these n,(m) “entities” transform among the selected I levels or exit from the reactor independently during the time interva1 [m At, (m + l)At], then each of these “entities” must be at one of the 1 levels, e.g. level j, inside the reactor with a probability of pij(m, m+ 1) or outside the reactor with a probability of pid(m, m + 1) at time (m + I) At. Therefore, these ni(m) “entities” will be multinomialiy distributed among the states 1, 2, . . . , 1 and d, with parameters q(m), m+l), Piztm, m+l), . . . , p&n, m+l), 1) at time (m+ 1)At [e.g. Rohatgi (1976)],

Pilb m+

. . . P,WP,Wl.

after the first transition is given

by

_I

Pr[Ni,(m,

m+ l)=ni,,

l)=~,

N,,(m,m+

..., q(m)!

Ni,(m,m+l)=nirlNi(m)=ni(m)]= Pj(l)=

5

i=l

Pi(O)Pij(O,

1X

j=

1, 2, .

. 1k d

(6)

I

ij=l

or, in matrix form, by P*(l)=P*(o)P(o,

1

X

1).

(7)

Similarly,

-i

n j=l

CPij(4

m+

111””

If

“id

PiAm,

m+

1)

1

.

(11)

Then, we have P*(2)=P*(l)P(l,

2)

1

(8)

ni(m)

=

C

P*(m+

l)=P*(m)P(m,

m+ 1)

E[Nij(m,

(9)

i= 1,2,

or i

i=l

nij+

p, .( m ) p,,..( m,m+l),

j=l,

2,. .

.I, d. (10)

Assume that n,(O) molecules of type Ai, i.e. n,(O) “entities” at level i, are initially in the reactor. Moreover, let the random variable X,(m) be the number of molecules of type Ai (or “entities” at level i) entering the reactor during the time interval [m At, (m + 1) At], and let the random variable N:(m) denote the number of molecules of type Ai in the reactor at time m At; -u,(m) and n,(m) are the realizations of X,(m) and N;(m), respectively. Furthermore, let the random variable Nij(m, m + 1) represent the number of “entities” at level j (or molecules of type Ai) inside the reactor at time (m + 1) At,

m+

2,

i=l,

flfd,

i= L

and, in general:

pj(m+l)=

p,Jm,

i.e.

. . ,I

l)]=n,(m)pij(m,

. . . , 1; j=l,

m+

(12)

l),

(13)

2,. . . , I, d

and var[Nij(m,

m + I)] = ni(m)pij(m, X Cl -PPij(W i= 1, 2, . . . , I; j=

Note

m + 1) m+

111,

(14)

1, 2, . . . , 1, d.

that

I Nj(m+l)=

c

i= I

Nii( m,m+l)+Xj(m),

j=l,

2, _. _, 1 (15)

and the random variables Nij(m, m + l), Nkj(m, m + 1) and Xi(m), i#k, are independent. Thus, the conditional mean and variance of Nj(m+ l), given that N,(m)= q(m) and X,(m)= x,(m), i= 1, 2, _ _ , J, are

Modeling of complex chemical reactions expressed, respectively, E[Nj(WI+

= E

l)[Ni(m)=ni(m),

i=l,

2,. . . , ;I]

1

N,j(m, m +I)

5

as [e.g. Rohatgi

=

5

E[Nij(?&

at level i and any molecule being outside the reactor (i.e. in the surrounding environment or dead state) is an “entity” at level d. Note that levels 1, 2, 3 and d of the “entity” or equivalently chemical species A,, and A,, and the surrounding environment d are the states of the Markov chain. Thus, the transition probability matrix for this Markov chain with n;(m) “entities” at level i, i.e. n,.(m) molecules of type Ai existing in the reactor at time m At is expressed as [see Too et al. (1983); also see the Appendix]

(1976)]

A,

xi(m)=xi(m),

1

+E[Xj(m)]

i=l

WI+ l)]

+ Xj(m)

i=l

1 P(m, m+ I)=

iz,ni(m)Pij(

=

m,m+

r

1 -k,n,(m)Ar-_~1At 0

0

d

0

0

j=

1, 2,. _ . , 1

l)lNi(m)=ni(m),

i= 1, 2,. . = Var

=

i

i-

$I

E[iV,(m+

=

1, 2,.

m+

-Pij(%

., 1

1)J (17)

for the case where E[Xj(m)] =o.

= xj(m)

and Var[xj(m)]

A chemically reacting system considered here is either a continuous stirred-tank reactor, CSTR, or a batch reactor, initially containing n,(O) molecules of type Ai. Let the number of molecules of this type entering the CSTR during the time interval [m At, (m + 1) At] be xi(m), and the intensity of each molecule of any type exiting from it at any time be p [e.g. Fan et al. (1979)]; naturally, both xi(m) and ,u are zero for the batch reactor. Example 1 competitive, Consider the following reactions carried out in a CSTR: +A,-

consecutive

fl A2

+A5 (18)

d,+d,-

At 1 1

i)~N,(m)=n,(m),

Xi(m)=xi(m),

t

n,(m)pij(m,m+

I)+xj(m),

j=l,

2, 3.

(20)

ECN,(m+

l)lNj(m)=ni(m),

X,(m) = x,(m),

i= 1, 2, 3, 4, 51

=n4(m)C1-~Atl--n,(m)p,,(m,m+1)

EXAMPLES

A,

fi

0

_ (19)

Provided that q(m) molecules of type Ai (i= 1, 2, 3, 4, 5) exist inside the reactor at time m At and xi(m) molecules of type Ai enter the reactor during the time interval [m At, (m + 1) At], the conditional mean number of chemical species A, inside the reactor at time (m + 1) At is obtained through the stoichiometric constraint as

m + l)] m + l)Cl

ni(m)Pij(m,

I-/IAt

1

i=l

Nij(m, m + I)]+varCXj(m)l

Var[Nij(m,

FAt

i= 1, 2, 3, 4, 5]

I

jil j=

[

PAt

k,n,(m)At

From eqs (16) and (18), we have

X,(m)=x,(m).

_,.J]

d

0

l-k,n,(m)At--fAf

(16)

var[Nj(m+

3

k, n4(m) At

0

and

=

2

1

21 3

l)+xJm),

2809

-n2(m)pZ,(m,

m+ l)+x&)

(21)

where pij(m, m + 1) is the 0th element of P(m, m + 1) given in eq. (19). Similarly, we have for chemical species A, ECN,(m+

l)lNJm)=n,(m),

Xi(m)=x,(m),

i= 1, 2, 3, 4, 51

+n,(m)p,,(m,

m+ 1)+x,(m).

(22)

Since the reactor initially contains n,(O) molecules of type Ai (i = 1, 2, 3,4, 5), the mean number of molecules of any type Ai in the CSTR at time (mt 1) At is obtained by iterating the sequence in eqs (20)-(22).

Ir2 A,+A,.

As mentioned in the preceding section, we regard’a molecule of type Ai inside the reactor to be an “entity”

Example 2 This example is identical to example 1 except that no reacting species exists initially in the reactor.

2810

S.

T. CHOU et al.

Example 3 This example considers the same reactions as those in examples 1 and 2. However, the reactions are carried out in a batch reactor instead of a CSTR. Example 4 Suppose that the following batch reactor:

reactions take place in a

Levels I, 2 and 3 of the *‘entity”, or equivalently chemical species A,, A, and A,, serve as the states of the Markov chain. The transition probability matrix for this Markov chain is expressed as 2

1

r

0

I -kk,At-kg,(m)At

k,n,(m) At 1

0

AS in example

I, we have i= I, 2, 3, 43

1)1N,(m)=n,(m),

E[Nj(m+

= $ ni(m)pij(m, i= I

m-t I),

j=

I, 2, 3

(25)

and ECN,(m+

l)INi(m)=ni(m),

m+

have the same value for the first-order reaction; both have the same dimension of l/time, i.e. the frequency. k is obtained by dividing k* by the product of the volume of the reacting mixture, V, and the Avogadro number, N,, for the second-order reaction. The solution and calculation schemes for example 1 are detailed in the Appendix.

DISCUSSION

I).

(26)

The mean number of molecules of any type in the reactor at time (mt l)At is obtained by iterating the sequence in eqs (25) and (26). Note that the governing equations, eqs (20), (21), (22), (25) and (26), are based on the numbers of molecules of the participating chemical species. The mean, E [Cj(m)], and the variance, Var [Cj(m)]. of the molar concentration of an individual chemical species, Aj, in the reactor are obtained, respectively, as [e.g. Rohatgi (1976)]

ECC,(m)l = ECNj(mVVNoI

= EtNj(m)l/vN,

(27)

and Var[Cj(m)]

1

(24)

i= 1, 2, 3, 41

=~4(m~+~2(m)~z,(m,m+l)--n,(m)p,,(m,m+l) -n2(m)pz3(m,

SIMULATION

The examples are numerically simulated in terms of temporal variation of concentrations of the reactants and products. The results obtained with the present approach (with At = I .O min) and the deterministic approach are compared in Figs f-4 for examples 1,2,3 and 4, respectively. To facilitate comparison, the reaction rate constants, k*s, shown in Figs l-4 are expressed in terms of the molar concentrations of the reactants. The values of the ks in eqs (19) and (24) are obtained from the k*s. It is worth noting that k and k*

3

k, n4(m) At

P(m, m+ 1)=2

the molar concentration of an individual chemical species is extremely small; it is approximately of the order of magnitude of 1O-12 [e.g. Fox and Fan (1987); also see the Appendix]. This indicates that the fluctuations of the molar concentration are negligible.

= Var [ Nj(m)/ VNo] = Var [ Nj(m)]/(

VN,)’

(28)

where V is the volume of the reacting mixture, and N,, the Avogadro number. It is worth noting that the standard deviation (the square root of the variance) of

The results of simulation from the present approach are in good agreement with those from the deterministic approach. This seems to indicate that the present approach is valid. The examples also indicate that the present approach is computationally effective. In applying a Markov chain model to simulate the dynamics of a chemically reacting system, the states of the chain should be identified in such a way that the system occupies one and only one state after each transition with a certain probability: this is the essential foundation of the Markov chain theory [e.g. Bharucha-Reid (1960)]. The states of a Markov chain may be specified, as in a conventional Markov process (i.e. Markov chain with At-O) simulation of such a system, by row vectors, the elements of which are the numbers of molecules of the chemical species partlclpating in the reactions [e.g. McQuarrie (1963, 1964, 1967), Oppenheim et al. (1968, 1977), Gillespie (1976, 1977a, b, 1981) and van Kampen (1981)]. As such, the states of this chain are countable. Nevertheless, the size of the chain is essentially infinite, rendering it cumbersome in obtaining the mean and/or higher moments of the number of molecules of each chemical species.

Modeling 1.0

.

A ,,+A,

.

A,+A,

l

0.8

z

of complex chemical reactions

\

nC

k,

C,, = C,. = 0.5 M

k; ---A,+4

c,,

= C3, = C,, = 0.0

C,,

= Czo = C,

iA2+As

C,r=C,,=l.OM

Determ‘ ~tnlsw-. moael .

-

b .*nvr -0,

.

k; = 0.025/M-min.

_ MS&r.., ...l...v~ chain model

:‘-.-._._._,

2811

_

Time

= 0.01 S/M.min

k; v=

r_rkq=:Ol/min._._~.

M

= 0.0 M

1OOOI

c,

(min)

Fig. 1. Temporal variation of concentrations of the reactants and products as functions of time for the competitiveeconsecutive reactions, eq. (18), in a CSTR with the initial condition Ci,, = Cc,-, i = i, 2, 3, 4, 5.

*

._

0.8

I

A,+A,

-

5

Y

l A0

C,,=0.5M,C,,=l.OM

kt

A,+A,---A,+&

Deterministic x+

C,< = C,, = C,, = 0.0 M

. kz -A,+A,

Markov

C 10 = C,,

= C,,

= C,,

k; = O.OZWM~min,

model

q = 10

chain model

160 Time

Urnin,

= C,,

= 0.0 M

k; = 0.015IM.min v

=

1000

I

200 (min)

Fig. 2. Temporal variation of concentrations of the reactants and products as functions of time for the competitive-consecutive reactions, eq. (18). in a CSTR with initial condition Ci, =O, i= 1. 2, 3, 4, 5.

Chemical species participating in the reactions have also been designated as the states in Markov chain simulation of a chemically reacting system [e.g. Miller (1972), Formosinho and Miguel (1979), Too et al. (1983) and Antia and Lee (l985)]. However, the chemical species to be identified as the independent states of the Markov chain must be judiciously selected. For example, if all the chemical species in the reactor are identified as the states of a Markov chain, the chain may occupy simultaneously two or more states after one transition, thereby causing the sum of the elements in a row of the transition probability matrix to exceed unity [e.g. Formoshinho and Miguel (1979) and Antia and Lee (1985)]. This violates the probability axiom [e.g. Bharucha-Reid (1960)].

CONCLUDING

REMARKS

The present paper describes an approach for circumventing the difficulties involved in simulating a chemically reacting system where the reactions are nonlinear. In this approach a molecule is regarded as an *‘object” or “entity”. This “entity” manifests itself as different chemical species. In other words, all the chemical species in the reactor are viewed as different forms or levels of this “entity”. The transformation of the molecule from one species to another is visualized as the transition of this “entity” from one level to another. Furthermore, the levels of this “entity” serving as the independent states of the Markov chain are identified subject to the stoichiometric constraint based on the atomic balance. The present approach

S. T. CHOU et al.

2812

leterministic l AOX*

Markov

-X-X

-x-x

/X-x

model

chain

model _X-X-X

C,,

= 0.5 M,

c,,

= c,,

C,, = 1.0M

k;

= O.OX/M~min

k;

= 0.01 S/M.min

= CSO = 0.0 M

l_a-a-*-*

Time

Fig. 3. Temporal

variation of concentrations competitive-consecutive

of the reactants and products as functions reactions, eq. (IS), in a batch reactor

I 1.8

-

Deterministic

. b xo

Markovchain A,+%

t

A,

gg 6 .5 E E 5 u

-*Cx

(min)

2.0

t=

C,

C ,D = C 2O = C, = k; =

1.2 *\_c

+

of time for the

model model

f

&A2

3

&

A3

1.00 M Cao = 0.00 M 2.00 M k; =O.OOl/M.mir 0.015/min

1.0 -

l

--._

.c4

C, 0.4 -

0.0

20

0

I

40

I 60

I 80 Time

Fig. 4. Temporal

can

probably

tion of more probabilistic

variation

be generalized. could

I 120

I 140

Simultaneous to violation

(min)

be prevented

constraints based on various physical laws such as atomic, mass and energy

as functions

by imposing and chemical balances.

Acknowledgements-The authors wish to acknowledge the Engineering Experiment Station of Kansas State University for partial financial support of this work and the constructive criticism of one of the reviewers.

of time for the

NOTATION

occupaof the

I 180

I 160

of concentrations of the reactants and products reactions given by eq. (23) in a batch reactor.

than one state, leading axiom,

I 100

type

Ai Cj(m)

c JI

molar species

concentration Aj in the reactor

molar species

concentration Aj in the feed

of

3

expected variable

chemical

at time mAt chemical of

stream

initial molar concentration cal species Aj in the reactor

cjO

EC

of molecule

value of the given or random vector

of chemirandom

2813

Modeling of complex chemical reactions k N(m) m + I)

Nfj(mv

Nkj(m, m+

N* nj(m)

p(m)

1)

reaction rate constant time-sequence of random vectors random variable representing the number of molecules of type A j inside the reactor at time (m+ l)At which have come from ni(m) molecules of type Ai due to reaction random variable representing the number of molecules of type A, exiting from the reactor during the time interval [mAt, (m+ l)At] which have come from ni(m) molecules of type Ai Avogadro number number of molecules of type Aj in the reactor at time m At transition probability matrix at time mAt one-step

Pij(m)

level

transition

i to level

probability

j at time

volumetric flow rate volume of the reacting

4

V

from

m At mixture

in the

reactor Var[

of the given random variable or random vector random variable representing the number of molecules of type Aj entering the reactor during the time interval [m Ar, (m + 1) At] flow rate in terms of the number of molecules of type Aj intensity of a molecule exiting from the CSTR during a time interval At, variance

]

Xj(m)

xj(m)

P

/c=q/V REFERENCES Antia, F. D. and Lee, S., 1985, The effect of stoichiometry on Markov chain models for chemica1 reaction kinetics. Chem. Engny Sci. 40, 1969-t 97 1. Bartholomav, A. F.. 1958a. Stochastic models for chemical reactions:-I. Theory of the unimolecular reaction process. Bull. math. Biophys. 20, 1755190. Bartholomay, A.‘F.: l958b. Stochastic models for chemical reactions: II. The unimolecular rate constant. Bull. math. Biophys. 21, 363-373. Bharucha-Reid, A. T., 1960, Elements of the Theory oj Markov Processes and Their Applica;ions, pp. IO-1 3. McGraw-Hill, New York. Darvey, 1. G. and Staff, P. J., 1966, Stochastic approach to first-order chemical reaction kinetics. J. them. Phys. 44, 990-997. Darvev, I. G., Ninham, B. W. and Staff, P. J.. 1966. Stochastic models for second-order chemical reaction kinetics. The equilibrium state. J. them. Phvs. 45. 214552155. Fan; L. T., Fan, L. S. and Nassar,k.. 1979, A stochastic model of the unsteady stage age distribution in a flow system. Chem. Enyng Sci. 34, I i 7221174. Formosinho, S. J. and Miguel, M. M., f979, Markov chains for plotting the course of complex reactions. J. them. Educ. 56, 582-585. Fox, R. 0. and Fan, L. T., 1987, Comments on a stochastic approach to the analysis of chemically reacting systems. Chem. Engng Sci. 42, 1861-1862. Fredrickson, A. G., 1966, Stochastic triangular reactions. Chem. Engng Sci. 21, 687-691. Gillespie, D. T., 1976, A general method for numerically simulating the stochastic time evolution of coupled chemical reactions. J. camp. Phys. 22, 403434. CES

43:

10-Q

Gillespie, D. T., 1977a, Concerning the validity of the stochastic approach to chemical kinetics. J. statist. Phys. 16. 31 l-318. Gillespie, D. T., 1977b, Exact stochastic simulation of coupled chemical reactions. J. phys. Chem. 81,234&2361. Gillespie, D. T., 1981, Conditioned averages in chemical kinetics. J. them. Phys. 75, 704709. fsida, K., 1960, The stochastic model for unimolecular gas reaction. Bull. them. Sot. Japan 33, 1030-1036. McQuarrie, D. A., 1963, Kinetics of small system I. J. them. Phys. 38, 33-436. McQuarrie, D. A., 1964, Kinetics of small system II. J. them. Phys. 40, 29 14-292 I.

McQuarrie, D. A., 1967, Stochastic approach kinetics. J. appl. Prohl. 4, 413478.

to chemical

Miller, P. J., 1972, Markov chains and chemical processes. J. them. Educ. 9, 222-224. Montroll, E. W. and Shuler, K. E., 1958, The application of the theory of stochastic processes to chemical imetics. Adv. them. Phys. 1, 361-399. Nassar, R.,- Fan, L. T., Too, J. R. and Fan. L. S., 1981, A stochastic treatment of unimolecular reactions in an unsteady state continuous flow system. Chem. Engng Sci. 36, 1307-1317. Oppcnheim, I., Shuler, K. E. and Weiss,G. H., 1968, Stochastic and deterministic formulation of chemical rate equation. J. them. Phys. 50, 460466. Oppenheim. 1.. Shuler, K. E. and Weiss, G. H.. 1977. Stochas-tic Processes in Chemical Physics: the Master Equation, pp. 431492. MIT Press, Cambridge, MA. Pippel, W. and Philipp, G., 1977a, Utilization of Markov chains for simulation of dynamics of chemical systems. Chem. Engng Sri. 32, 543-549. Pipped, W. and Philipp. G., 1977b. An improved aDDrOaCh of simulating chemical reactions by a ic?arkov-chain cell model. Chem. Enana Sci. 32, 1535 1536. Rohatgi, V. K., 1976, An Introduction to Prohahility Theory and Mathematical Statistics, pp. 197-199. John Wiley, New York. Singer. K.. 1953. Application of theory of stochastic processes to the study of irreproducible chemical reactions and nucleation processes. J. R. statist. Sot. B15, 922106. Too, J. R., Fan, L. T. and Nassar, R., 1983, Markov chain models of complex chemical reactions in continuous flow reactors. Corn; them. Engng 7, I-12. van Kampen, N. G., I98 I, Stochastic Processes in Physics and Chemistry. pp. 180-208. Elsevier. New York. APPENDIX:

ELABORATION

OF EXAMPLE

This example considers the following secutive reactions carried out in a CSTR:

1

competitive.

con-

r: A,+A,---+

t:

A,+A,The rates of these reactions

A,+A,

(Al)

A,+A,.

(A2)

are expressed.

respectively,

as

;c,

=lc:c,c,

(A3)

$c,=c:c,c,.

(A41

and

The operating C,,

conditions =C,,

C,,=C,,=

are such that

=0.5 M. 1.0 M.

C,~=C,,=C,,=C,,=C,,=C,,=OM,

2814

S. T. CHOU et al. k: =0.025/M kr = 0.0 q=‘lO

. min,

Dividing both and N, yields

by the product

of v

n,(m) C,(m)-C,(m+l)=Ckl(VN v~ At 01C,(m)

1S/M min, I/min,

1

and

V= 10001.

=Ck,( VNo)ICl(m)C4(m)At.

Ident$cation qf the states of the Markov chain As mentioned in the text, we regard a molecule of type Ai inside the reactor as an “entity” at level i, and any molecule outside the reactor (i.e. in the surrounding or dead state) as an “entity” at level d. Note that an “entity” can occupy one and only one level after each transition. Levels 1,2, 3 and d of the “entity or equivalently chemical species A,, A, and A,, and the surrounding environment d are the states of the Markov chain. Derivation of the transition probabilities Consider two consecutive times, m At and (m + 1) At, where At is sufficiently small so that what follows is valid. An “entity” at Ievel 1 at time m At may move to level 2 but not to level 5 which is not a state of the Markov chain due to reaction (Al); it may move to level d by exiting from the reactor or it may remain at level I at time (m+ 1)At with probabilities governed by certain probabilistic laws. Let k, and k, be the intensities (or frequencies) of occurrcncc of reactions (Al) and (A2), respectively. k, is interpreted as the intensity of the effective collision of molecules of types A, and A, leading to the occurrence of reaction (Al), and k, as that of the effective collision of molecules of types A, and A, leading to the occurrence of reaction (AZ). Note that, in general, k, and k, correspond but are not identical to ky and k:, respectively. Suppose that ni(m) molecules of type Ai, i = 1,2,3,4,5: exist in the reactor at time m At. The probability of a molecule of type A, to react with any of the n.,(m) molecules of type A, during the time interval [m Af,(m+ l)At] is n,(m)k, Ar. In other words, the transition probability of an “entity” at level 1 at time m At being at level 2 at time (mt l)At through reaction (Al) can be expressed as plz(m,m+

I)=k,n,(m)Ar.

(A5)

The probability of a molecule of type A, in the reactor (i.e. an “entity” at level 1) at time mAt exiting from the reactor (i.e. moving to level d) at time (mt 1) Ar is p Ar [e.g. Fan et al. (1979)], or p,d(m,m+l)=~Af. Consequently,

(Ah)

we have pll(m,

m+ l)=

1 -kk,n,(m)At-pAt

(A7)

By the same argument, we eventually obtain the transition probability matrix for this Markov chain with n,(m)“entities” at level i in the reactor at time m At; it is 1

P(m,m+l)=

which

sides of this expression

2

1

I -k,n,(m)At

2

0

3

0

0

d [

0

0

-pAr

k,n,(m)

At Ar - ~1At

I - k,n,(m)

Recovery of ks from k’s Since n1 (m) molecules of type A, exist in the reactor at time m Ar and the probability that each of these molecules react with any of the n*(m) molecules of type A., during the time interval [m At, (m + 1) At] is k, n4(m) At: I)=nl(m)k,n,(m)At =k.n.(mln.(m)At.

(A91 I

(AlO)

eqs (A3) and (AIO) gives k:=k,(VN,,)

(All) Similarly

(A121 In general k* k= for a yth-order

(A13)

1

( VN,)‘_ reaction.

Fluctuations of molar concentrarion Expanding eq. (17) in the text system yields VarCN,(m+

l)IN,(m)=n,(m),

for

the present

reacting

X,(m)=x,(m),

i= 1, 2, 3, 4, 51 m+l)~l-pll(m,m+l)]

=n,(m)p,,(m, Var[N,(m+

l)INi(m)=n,(m).

(Al4) X,(m)=x,(m),

i= 1, .2, 3, 4, 51 =n~(m)P~,(m,

1) [l -pZ2(m,

mf

m-t

+n,(m)p,,(m,m+!)Cl-_p,,(m,

I)]

(A151

m+l)l

and VarCN,(m+

l)IN,(m)=n,(m),

Xt(m)=x,(m).

i= 1, 2, 3, 4, 51

= n3(m)p33(m, m + 1) L1 -pa3(m,

m+

l)]

+n,(m)p,,(m.m+1)Cf-_p,,(m,m+l)l. The conditional variances obtained, respectively, as var[N4(m+

of N4 (m+

l)IN,(m)=n,(m),

1) and

(J416) N,(m+

1) are

X,(m)=x,(m),

i= 1, 2, 3, 4, 51

(A81 11

=n,(m)p,,(m,

mf

3

d

0

PAt

k,n,(m)At

1) Cl -pP44(m,

m+

111

(Al71

,uAt

0 -pAt

is eq. (19) in the text.

n,(m)-nl(m+

Comparing

1 PAt

and VarCN,(m+ i=

l)lN,(m)=n,(m),

X,(m)=x,(m),

I, 2, 3, 4, 51

=n5Cm)pss(m.

m-t

tJL1 -pss(m, l)[l

+n,(m)p,,(m,

m+

-pls(m,

+n2Wp2,(m.

mf 1)Cl --P2sh

m+

I)] m+

l)]

m+ 111.

(AffJ)

Modeling Note

that

for the present

example

P,s(‘n,

m+l)=p,,(m,

m+

1)

Pzs(W

m-t

m+

1)

l)=p,,(m,

P44(~~.n~+l)=1-p,*(m,m+1)-pp,,(m,m+1) - P,,(W

m +

1)

pss(m,m+l)=l--p&m, P_,d(m,

m+

l)=p,,(m,

m+l) mt

I)=pAr.

of complex

chemical

reactions

2815

The simulation results (with Ar = 1.0 min) show that the 1. I), i= 1. 2, 3, 4, 5, m=O. of N,(m+ conditional variances 2..... ranee from 0.002 x (6.02 x IO*“) to 0.022 x (6.02 i 10z3), or equivalently the donditional ‘variances of C,(m +1),i=l.2,3,4,5.m=0.1,2 .,._, rangefrom3.0x10-“Jto 3.7 x 1O-3z M*. The standard deviation of the molar concentration range from 5.48 x 10-l’ to 1.92 x lo- lh M, thus indicating that the fluctuations of the molar concentration are negligible.