J . Phys. Chem. 1992, 96, 3958-3962
3958
1/(N4)
+ 1/(1
- 4)
- 2x = 0
with the identification x = 2c/kBTfor the mean field interaction parameter, (where the numerical factor is obtained as (2- 2)/2, where Z is the coordination number of the lattice, Le., 2 = 6 in our case, ignoring the effects due to terminal units). It is observed how the simulation curve shows a considerably broadening around the critical point, in comparison with the theoretical curve. (This effect hinders a very precise determination of the values of aCand it may have affected to the determination of the exponent for eq 13.) The same broadening has been reported in comparisons between phase separation curves of experimental systems and their mean field prediction^.','^^^^ Then, our simulation curves reproduce the most important qualitative aspects of experimental data, when simulation and experimental results are both compared with the mean field theory. As is observed in Figure 7, the broadening in the critical region increases with the chain length, as the curves become simultaneously more asymmetric. A quantitative study of the curve characteristics in the critical region is somehow difficult because of the previously discussed broadening. From our numerical results for the spinodal curve, we have attempted to obtain the critical exponent (3, involved in the scaling law which describes the difference between the lower (15) Koningsveld, R.; Kleintjens, L. A. Macromolecules 1971, 4 , 637. (16) Chiu, G.; Mandelkern, L. Macromolecules 1990, 23, 5356.
and higher phase separation concentrations at a given temperature, @ff
- *f,2
@" - @'
-
(C/kBT- e/kgTc)BN-1/2+8/2
(15)
where renormalization group theory predicts /3 = 0.327,2*3J4 yielding the exponent -0.34 for the dependence on N . We have not been able to study the dependence with c/kBT, since most of our choices for this parameter do not range close to the critical point. However, using c/kBT = 0.35 and the three longest chains (for which we obtain phase separation at this temperature), we have obtained
in good agreement with the theory. However, only a formidable computational effort could enable one to reproduce accurately the different exponents involved in the critical scaling laws for polymer-solvent systems at the present stage. Nevertheless, we believe that the present work shows how many interesting features of the phase separation curves can be computed from simulated scattering experiments, following relatively simple procedures. Acknowledgment. This work has been partially supported by Grants PB0093-89 and PB0198-89 from the DGICYT (Spain). A.L.R. also acknowledges a fellowship from the PFPI (Becas Predoctorales).
Nonequilibrium Concentration Fluctuations in Nonldeal Reacting Systems Robert M. Mazo* and John Guslandert Institute of Theoretical Science and Department of Chemistry, University of Oregon, Eugene, Oregon 97403-5203 (Received: September 1 1 , 1991)
The birth-and-death formalism is applied to chemical reactions in nonideal (in the thermodynamic sense) systems. The combinatorial transition probabilities conventionally used in the ideal theory are replaced by activity-dependent transition probabilities. Exchange of reactants and products with reservoirs is explicitly taken into account. The master equation is expanded in inverse powers of the system volume, and a Fokker-Planck equation is derived in the large system limit. Anomalous fluctuations near instability points are discussed.
I. Introduction Ever since the early work of Wilhelmy,' the results of investigations in chemical kinetics have been described in terms of rate equations. These equations relate the rate of a reaction to the concentrations of the several reacting species. The concentrations involved are mean concentrations, averaged over volumes of the reacting mixture large on a molecular scale. Thus, the classical description of chemical kinetics is a mean-field approximation to the actual phenomena. It has, of course, long been recognized that chemical reactions take place through a series of more or less independent events, or collisions. Thus, fluctuations about the mean-field description, Le., about the solutions of conventional rate equations, are to be expected. Such fluctuations are expected to be quite small under normal circumstances, just as fluctuations about equilibrium states are very small under normal circumstances. Hence, although their existence has been recognized for a long time, they have been ignored as being of no experimental significance. Some theorists have studied fluctuations in chemically evolving svstems. The literature through 1967 has been reviewed bv McQuarrie.2 At approximately ihis time, the subject of chemical instabilities and chemical oscillations began to assume some Present address: Department of Chemistry, University of Montana, Missoula. MT 59812.
0022-365419212096-3958$03.00/0
pr~minence.~Since fluctuations are generally enhanced at instability or transition points, the study of chemical fluctuations took on new life. The theoretical treatment of fluctuations in reacting systems has been a mesoscopic one, intermediate between the macroscopic, phenomenological description and a purely microscopic, molecular scale treatment. The mathematical tool commonly used is the theory of stochastic birth-and-death processes. A master equation is assumed to govern the time evolution of the probability distribution of the system. The transition rates of this master equation depend on the populations of the several reacting species in a way suggested by mass action kinetics. Thus, there is an assumption in the theory that the reacting system is thermodynamically ideal. This can clearly be seen from the circumstance that mass action kinetics necessarily leads to an equilibrium state in which the equilibrium quotient is expressed in terms of concentrations rather than activities. With an appropriate master equation in hand, one is faced with the problem of solving the equation and drawing chemical con(1) Wilhemy, L. Pogg. Annalen 1850, 81, 413. (2) McQuarrie, D. J . Appl. Prob. 1967, 4 , 570. (3) Prigogine, I.; Nicolis-,G. J . Chem. Phys. 1967, 46, 3542. Lefever, R.; Nicolis. G.; Prigogine, I. J . Chem. Phys. 1967, 47, 1045. Turing, A. M. Philos. Trans. R.SOC.London 1952, 8237, 37.
0 1992 American Chemical Society
Nonequilibrium Concentration Fluctuations clusions. Several methods of solution have been used. One of these is to form the generating function of the probability distribution function of the molecular population^.^^^ One obtains a partial differential equation for this generating function. Another method is to carry out an asymptotic expansion of WKB type on the master e q ~ a t i o n . ~Still another is to carry out a direct asymptotic expansion on the master equation with the fluctuating variables properly scaled.6 If we turn our attention to real, rather than ideal, systems, then further interesting effects occur. This problem has been studied in the macroscopic description by Li, Nicolis, and Frische7 These authors show that new steady states and bifurcations can arise when thermodynamic nonideality is taken into account. The question arises of how to take fluctuations into account in the treatment of nonideal systems. It is not merely a matter of replacing particle numbers or concentrations by activities, because one will want to hold the concentrations of some of the chemical species, e.g. the major species, constant, just letting the intermediates fluctuate. Since activities depend, in general, on the concentrations of all of the species present, holding activities constant and holding concentrations constant are not equivalent. So one must reformulate the problem from the beginning. Furthermore, however one reformulates the problem, one must be sure that the solution reduces to the known expression for equilibrium fluctuations when the constraints are such that the system is in fact a t equilibrium. In this paper, we present a formulation of fluctuations in a chemically reacting open system which obeys all of the necessary subsidiary conditions of which we are aware. For the implementation of general ideas we have used the van Kampen system size expansion.6 The choice between this and the WKB type expansion of Kubo et aLS was largely a matter of taste. The generating function method, while applicable in principle, turns out to be convenient only if the transition probabilities are polynomials in the particle numbers. In nonideal systems, polynomial transition rates are not the case. The paper is organized as follows. Section I1 describes the model we propose and formulates the master equation for this model. Section I11 contains the system size expansion for this master equation and its solution, and section IV contains discussion and concluding remarks. 11. The Model In macroscopic chemical kinetics, one conventionally writes the kinetic equation as a
B
where ci is the concentration of species i, and u,,~(uB,J is the rate of formation (destruction) of i by the ath (8th) elementary step in the reaction. This rate is conventionally written as
The Journal of Physical Chemistry, Vol. 96, No. 10, 1992 3959 rate laws lead to an approximate, rather than exact, equilibrium expression, they cannot themselves be exact. Let us emphasize that we are not here discussing the quality of the approximation, which indeed has shown itself to be adequate for many applications. Rather, we are drawing attention to a matter of principle. From a formal point of view, there are several ways in which one could write down *corrected” rate laws. One is to let the rate “constants”, k,, be dependent on concentration. Unless one knows a priori what this concentration dependence is, such an artifice is useless for prediction or interpretation. Another approach is to retain the structure of (2), with constant k’s, but replace all concentrations (e,) by activities (aj). This approach has the advantage that the rigorous equilibrium expression is obtained in the long time limit. Furthermore, the degree of closeness of the revised rate laws to the original ones (2) can be easily assessed in any given case. We are not aware of any rigorous demonstration that merely substituting activities for concentrations in (2) is the correct thing to do. On the other hand, it does have much to recommend it, and it is supported by transition-state theory, which posits that reactants are in equilibrium with the activated complex.s Therefore, it has been adopted by many author^.^ We shall also adopt it, but we do not want to do so without emphasizing its provisional character. Now, let us shift our attention to the mesoscopic, or stochastic, treatment of kinetic rate laws. In stochastic theory, the hypothesis is made that the chemical kinetic processes can be modeled by a birth-and-death process obeying a kinetic equation written schematically as
where P(X,t) denotes the probability that the system has composition X at time t , X denotes the particle numbers of the various species present, and r denotes the changes of these numbers in the several elementary steps. W is the transition rate for the change indicated by its arguments. The sum is over all elementary reactions contributing to a mechanism, including flow in and out of external reservoirs. Here and in the following, upper case letters, e.g. X, will stand either for the name of species X or for the number of particles of species X present; it will always be clear from context which meaning is meant. In the treatment of ideal systems, the forward transition rate for a given elementary reaction vAA
for the transition ( A , B, ..., M ,N , ...)
(4) Nicolis, G . ; Babloyantz, A. J . Chem. Phys. 1969, 51, 2632. (5) Kubo, R.;Matsuo, K.; Kitahara, K. J . Stat. Phys. 1973, 9, 5 1 . (6) van Kampen, N . G . Stochastic Processes in Physics and Chemistry; North-Holland: Amsterdam, 198 1. (7) Li, R.-S.; Nicolis, G.; Frisch, H. J . Phys. Chem. 1981, 85, 1907. Li, R.-S.; Nicolis, G . Ibid. 1981, 85, 1912. Li, R.-S.; Nicolis, G . J . Chem. Phys. 1982, 77, 2365.
-
( A - V A , B - VB, ..., M
+
VM,
(4)
N
+ V N , ...)
is taken to be kA(A - 1)
where k, is the rate constant, and u,, is the number of moles of j occurring on the left-hand side of the a t h elementary reaction. It is well recognized that the exponent of cj is only ua, for an elementary step. It is clear that this formulation of the phenomenological rate equations, valuable as it has proved during 130 years of use, cannot be completely correct. For as the reaction tends toward equilibrium, the rate equations become equivalent to the well-known equilibrium constant expression written in terms of concentrations. Rigorously, however, the equilibrium constant expression must be written in terms of activities, not concentrations. Since the
+ v ~ +B ... * vMM + vNN + ...
... ( A - V A + 1) B ( B - 1) ... ( E - v g v*!
Vg!
+ 1) ...
(5)
where k is a rate constant for the given reaction. The backward rate has a similar structure. k for the forward reaction is proportional to Y(”A+”B+...)+I, where Vis the volume of the system. It is the assumption embodied in eq 5 which ensures that the system described is thermodynamically ideal. We adopt here the hypothesis that a nonideal system is described by a master equation of the same form (3) but that the transition rates are to be constructed by multiplying all particle numbers, A, B, ..., by activity coefficients. For example, for the elementary reaction 2X anything, one would have, according to (3,kX(X- 1)/2 as the transition rate multiplying P(X,t). We 1)/2. The activity coefficients (y) are use instead ky,y,,X(Xgenerally functions of the mole fractions of the various components;
-
(8) Moreau, M. Physica 1981, 102A, 389. ( 9 ) Keizer, J. Statistical Thermodynamics of Nonequilibrium Processes; Springer-Verlag: New York, 1987.
Mazo and Guslander
3960 The Journal of Physical Chemistry, Vol. 96, No. 10, 1992
Reactant reservoirs
Reaction chamber
Product reservoirs
Figure 1. Schematic depiction of the reaction chamber with reservoirs of products and reactants. The concentrations in the reservoirs are fixed and do not fluctuate.
they are intensive quantities. By yx we mean that X is used in computing all mole fractions. By yx-, we mean that X - 1 is used in computing all mole fractions, etc. Note that, in the ideal case, the transition rates are polynomials in the particle numbers. In the nonideal case, this is no longer true, in general. It is usually the case that the logarithms of the activity coefficients depend relatively simply on the mole fractions, not the activity coefficients themselves. This means that the generating function method, which has been so useful in solving the master equation for ideal systems, is not particularly helpful in the present case. There is one more point to be introduced to complete the model. We wish to be able to study steady states. This implies that we must have a source of reactants, and a sink for products. Otherwise, the only steady state possible is equilibrium. We therefore introduce a set of reservoirs for reactants (Ri)and a set for products {Pi), The concentrations in the reservoir are not supposed to fluctuate; no reaction takes place in the reservoirs. The transition rate of a reactant (product) species into the system from its reservoir is assumed to be k?p,(res) (klpai(res)). Here lower case letters refer to concentrations, just as upper case letters refer to particle numbers. The rate constants k:,orp are volume independent. The transition rates out of the reaction volume into the reservoirs are assumed to be k,& (k;+,Pi) for reactants (products). The general scheme is indicated schematically in Figure 1. The master equation has the form
a large system with co fmed! It is conventional, however, and more convenient in writing formulas, to pretend that Y ’itself is our small parameter, and we shall follow this convention. The system size expansion for systems described by a master equation has been developed in detail by van Kampen.6 Dekker’O showed how the van Kampen expansion could be made valid in the neighborhood of an instability. However, Dekker only considered the case of one independent variable. He also assumed certain symmetry properties of the jump moments (defined below). Both of these restrictions need to be removed for applications to chemical reactions. Van den Broeck et al.I1 have treated the case of the FokkerPlanck equation of a two-variable system near an instability. Although they start from a Fokker-Planck equation, their problem is formally very similar to ours, and we have benefited by their methods. We write the master equation in the form of eq 3 and assume that the transition probabilities, W, are expressible in the form (7)
where x = X/ V, the set of concentrations, and r is the set of jumps. Let us define jump moments, a(“), by a(”)(x)= Xr“w(x,r) r
where r“ = rr...n factors ...r. Furthermore, let us write x =%
+4
P(X,t) = p(4,t)
(8)
That is, j z is some reference set of concentrations and 4 represents fluctuations about this reference set. p depends on X parametrically. There is a great deal of freedom to choose the set j t to simplify the resulting equations. The master equation can then be written a,p(x,t) = t-iCw(x-cr,r) p(x-tr,t) - w(x,-r) p(x,t) (9) r
Introducing the shift operators Ei = exp(-e a/aXi) d,p(x,t) =
c-’C(nEri- 1) w(x,-r) p ( x , f ) r
where R denotes reactants, P products, and X intermediate species. Thus, all species in the reaction volume may fluctuate according to the stochastic dynamics of the mechanism. The reservoirs are supposed to be controllable and do not fluctuate. For simplicity, we have not allowed reactants to flow into the reactant reservoir nor products out of the product reservoir. Detailed balance, therefore, does not hold for these processes. To summarize, we generalize the master equation which has been used for ideal systems in two respects. We associate with each particle number factor in the transition rate an appropriate activity coefficients, and we introduce product and reactant reservoirs, by means of which we can keep the system in a nonequilibrium steady state. Neither of these generalizations has been derived, in any sense, from a more fundamental theory, but then neither has the ideal master equation itself been so derived. We believe that our proposed generalizations are reasonable, and they are the simplest reasonable ones which we could think of. It is well worth working out their consequences. 111. System Size Expansion
In this section, we take the master equation, eq 6, and carry out a systematic expansion in the small parameter c = VI,where V is the size of the system. Of course, in picking a small parameter, one should in principle always choose a dimensionless parameter in order to keep the ‘smallness” invariant to the units one is using. Although VI does not obey this criterion, expanding in VI is equivalent to expanding in XO-l where Xois some typical where number of particles in our system. That is, VI = Xo-lc0, co is a typical concentration, and one is interested in the limit of
i
nJJ,
+ +
Here the notation is as follows: j! = hl = j , j , ...; Vj:d) stands for the complete contraction 13(a/aXlyl(a/aX2Y2...,Ci). The usual procedure is to scale the fluctuation quantities 4 according to = t1/2qin the limit of small t. As is well-known,6 this works well away from an instability but not at or near an instability. A different scaling is called for. An instability can occur when one (or more) eigenvalue of Vu(’) vanishes. In this paper, we shall assume that only one eigenvalue vanishes. Another possibility is that the real parts of two eigenvalues simultaneously vanish. This is a Hopf bifurcation and leads to a periodic or oscillating reaction. We shall not discuss this possibility here, nor shall we treat the case of a degenerate zero eigenvalue. In order to isolate, as far as possible, the effect of the vanishing eigenvalue, we transform to a coordinate system in which Vcu is diagonal: 4 l. These { coordinates are then assumed to scale as follows. Let subscript 1 refer to the coordinate associated with For all other rs, we the small eigenvalue, A,. Then = set {,= taqj. The 7’s are supposed to remain fixed as e 0. Furthermore, we set AI = ccA’, where A’ also remains fixed as t
-
-
-
0.
Inserting these scaled variables in (lo), and expanding the jump about %, one finds moments, d), (10) Dekker, H. Physica 1980, 103A, 5 5 , 80. ( 1 1) Van den Broeck, C.; Malek Mansour, M.; Baras, F. J . Star. Phys. 1982, 28, 557.
The Journal of Physical Chemistry, Vol. 96, No. 10, 1992 3961
Nonequilibrium Concentration Fluctuations
fl is known, @ is determined. However, in some cases, it may turn out that this relation gives an e exponent in the dominant streaming term which is negatiue. This is unacceptable, so the procedure must be changed to choose a to make this exponent zero. Then the exponent of the diffusion term is fixed, and a may not be Or it may turn out that a = but that there are other streaming terms of order eo besides the Ai(a/hi)qi terms. This would be the case, for example, if a$!)1 = 0, while a$;:), # 0 for some i. Many more possibilities exist. It is best to consider each case separately, as it arises, using the principles enumerated above. As previously stated, the macroscopic kinetic law is given by R = a ( l ) ( x ) . For chemical systems, therefore, some of the terms a$!i,which represent the second-order terms in an expansion of a($ about 8, will be nonzero; these will certainly arise from second-order terms in the kinetics. Therefore, except in very extraordinary circumstances, there is no need to go beyond the terms retained in eq 11, even in the case of nonideality. a and
Note that the terms in e-] vanish, since we have chosen j2 to satisfy in the usual way.6 We have the macroscopic equations k = only displayed terms containing second derivatives or less. The question is now, how should one choose the scaling exponents a and @? Since q is the critical variable and, as we shall see, the slow variable, let us integrate over all vi, j > 1. One then finds an equation for p ( q l , t ) :
Here ( qj) is the conditional average of qj given ql, etc. Now take @= and c = This choice will be discussed below. With this choice, (12) becomes, to lowest order in e
a=
Thus, the probability distribution of q1 relaxes on a slow time scale of to its steady value. Of course, q1 does not evolve completely independently of the other variables, since eq 13 contains the conditional average ( q l q j ) ,i, j > 1. With the choice of a and 0 made above, eq 11 becomes, to leading order in e
Since the right-hand side of eq 13 is independent of ql, the solution has the form p(q,t) = 4(q2....,t ) X +(ql). On the time scale of the evolution of q2, ..., the distribution of ql does not change, its evolution being O(e113) slower. Thus, we may identify $(ql) with the p ( q l , t ) previously obtained (eq 12). Thus, we have found that the noncritical variables, q2, ...,relax to their Gaussian steady-state values on a time scale of the order Acl, where A, is the smallest of X2, .... The critical variable, q I , relaxes to a non-Gaussian steady-state value in a time of order longer. The conclusions stated above about the magnitudes of the exponents and the time scale of the evolution are not universal and depend in an essential way on the whether or not the coefficients a& vanish. As given, the conclusions are valid for the case a# # 0 for some i. If @ < a as we expect on physical grounds, the lowest power of e which occurs in the a dhl, streaming terms (i > 1) is 2 a - j3. It is also required that afll # 0, so that the smallest power o f t occurring in the a/aql streaming terms is @. If these conditions on the coefficients are violated, the exponents will be different. The general procedure for determining the exponents is as follows. First, choose a and @ so that the diffusion term with e order 1 - 28 in the dg(ql,t) equation balances the dominant streaming term. This gives a relation between a and 8; it does not determine them uniquely, of course. Then balance the diffusion terms in eq 10 (with e order 1 - 2a) against the since the streaming terms. Normally, this will lead to a = leading streaming terms are of order eo. Since the relation between
I
IV. Discussion This paper actually has two relatively independent themes. The first is the introduction of nonideality into the master equation description of homogeneous chemically reacting systems. The second is the reduction of such a master equation to a FokkerPlanck equation in the limit of a large system. In fact, the second theme was implemented in a rather general way, independent of the detailed form of the master equation being considered. We proposed to replace particle numbers (or concentrations) in the transition rates of the master equation by activities. Although this idea is an old one in the context of rate equations, we have emphasized its provisional nature. Furthermore, because the activity of a given species depends, in general, on the particle numbers of all species, we proposed to replace the constraint of fied particle numbers of the major species by a constraint of fiied feed rates of reactants and products from external reservoirs. The internal fluctuations of all species are then determined by the assumed stochastic dynamics. Of course, this could be done in the ideal case, too, but does not seem to be necessary there. Our treatment of the reduction to a Fokker-Planck equation is deeply influenced by the work of Dekker'O and Van den Broeck et al." Dekker made a thorough study of the expansion of a master equation near a critical point. H e considered only a onedimensional system but treated correction terms to the lowest order results in some detail. Van den Broeck et al. started from a Fokker-Planck equation whose diffusion term is proportional to a small parameter, e, and studied how the solutions scale with e for cases of two variables. We have considered here the expansion of the master equation for more than one independent variable. The formulation leads to equations of the same form as those of Van den Broeck et al. The extension from two to many variables causes only notational complications. We claimed that one will usually find for the scaling exponents a and 8, which determine how the fluctuations scale with system size, a = 1/2 and @ = At first sight, this seems to contradict the examples of Van den Broeck et al. In fact, however, our conclusions are consonant with theirs. The apparent contradiction arises from the fact that the examples emphasized special cases with dramatic effects on the fluctuations, whereas our discussion centered on the case we expect to be most common in chemical systems. For example, their example B, which has a = and j3 = I/.,, is precisely the case a\!)lI= 0 and a$:,)l# 0 mentioned in the previous section. For this case, we also get and I/., for the exponents. In interpreting all of these results, one must bear in mind that the variables Vi, or their unscaled versions ci, are not the chemical variables of interest but are linear combinations of them. Recall that these variables were chosen to diagonalize the linear part of the streaming terms and so isolate, in part, the unstable mode. In general, a given chemical species will have contributions from all (or at least many) of the {variables and so will contain terms scaling both as e" and as e@. Since @ < a,in general, the eo terms will dominate. Only if, for some reason depending on the reaction
3962
The Journal of Physical Chemistry, Vol. 96, No. 10, 1992
mechanism, the coefficients of all the €0 terms vanish will a chemical species concentration scale with ea. In order to apply our results, one must (a) choose a chemical mechanism and decide how the activity coefficients are to be modeled, (b) compute the jump moments a(]) and d2), and (c) use the results of section 111 or of ref 11 to compute the fluctuations. Far away from an instability (all eigenvalues of Va(') # 0), the scaling X1 = h'ec is not relevant. One should scale all fluctuations by as in the ideal case far from an instability. In that case, the fluctuations of the original chemical species become Gaussian with variance u given by = -2Va(1).(a(2))-l
Mazo and Guslander Appendix In this appendix, we show that our result for the probability distribution in a steady state agrees with that given by statistical mechanics when the steady state is, in fact, the equilibrium state. For this, it suffices to consider the system far from an instability, since an equilibrium state is stable. Therefore, it is not necessary to make the linear transformation to the 5 variables. Furthermore, a = 0 = ' I 2(and c = 0). The master equation (1 1) becomes, in the steady state
-va('yn):v(€p)
Acknowledgment. This work was supported, in part, by NSF Grant CHE-8609377.R.M.M. thanks Professor G. Nicolis for very useful discussion about this problem. (12) de Groot, S.; Mazur, P. Non-Equilibrium Thermodynumics; NorthHolland: Amsterdam, 1962; p 101. (13) Elson, E. L.; Magde, D. Biopolymers 1974, 13, 1. (14) Lemarchand, H.; Nicolis, G. Physic0 1975,824 521. (15) Gardiner, C. W.; McNeil, K. J.; Walls, D. F.; Matheson, I. S.J . Sfat. Phys. 1976, 14, 309. (16) Kuramoto, Y . Prog. Theor. Phys. 1974, 52, 718.
('41)
The solution to this equation which vanishes at infinity is
(15)
as follows from the usual form of a Fokker-Planck equation.12 Even though fluctuations far from an instability are small, they are measurable and have yielded, in some cases, useful information on both reaction rates and diffusion ~0efficients.l~ Near an instability, the fluctuations will be enhanced and will also give useful information about the instability. To our knowledge, measurements of fluctuations in a reacting system near an instability have not yet been made. Perhaps the greatest shortcoming of this development is that it focuses on global fluctuations of the reacting system. What is more interesting, and what is more closely related to what is measured, are local fluctuations, fluctuations in small volume elements (still on a mesoscopic scale). To study such matters, one must take diffusion into account, as has been done in ideal systems already.1"16 We hope to return to this problem in a future publication. In order to keep the fluctuations in the system from being artificially damped by constraints of constancy of certain concentrations, we had to couple the system to external reservoirs of products and reactants. But such a coupling introduces its own fluctuations. This is true even in equilibrium. Average values of thermodynamic quantities do not depend on the ensemble used to compute them, but fluctuations depend on the ensemble. In the nonequilibrium case, the fluctuations due to mass flow will be sensitive to the mechanism, while in equilibrium this is not the case. It seems that there are two possible scenarios. Either the fluctuations due to mass flow coupling are small compared to fluctuations due to chemistry, or they are not. In the former case we have the chemical effects more or less isolated for study. In the latter, we have the opportunity of studying how the system interacts with its surroundings, a problem of considerable current interest. To study this, one needs a more realistic model of the mass transport than that used here. Let us emphasize one more time that a basic element of the theory here presented, the expression of the transition probabilities in terms of factorial activities instead of factorial concentrations, is not unique. It gives the correct equilibrium limit, and is suggested by transition-state theory, but is by no means definitely established as correct. Nevertheless, these transition probabilities are reasonable and are in the same spirit as those used in the ideal case. Their consequences are well worth pursuing.
+ '/2a(2)(R):VVp = 0
p = ( 2 ~ ) - ~ / ~ ( u)-ll2 d e t exp(-'/2Q.u-1.€)
(A21
where u-I is given by eq 15. At equilibrium, the distribution function has the same structure as (A2), but u-' at equilibrium is given by
('43) where xi is X j / V , not mole fraction. So, to show agreement between the mesoscopic and microscopic formulas, we need only show that (A2) reduces to (A3) at equilibrium. We shall assume, for the purpose of simplifying the notation, that we have only one elementary reaction. In this case cy(]) = ri(wr- wb) ('44) where wf and wb are the forward and backward transition rates, respectively; ri is the change in the number of species i in the reaction. Similarly
+ wb)
ah*)= rirj(wf
We are working in the limit where we can consider the 5 variables as continuous; this is implied by even writing (Al). Hence we may take
Therefore
Here u~ and ukb are the number of molecules of species k taking part in the forward or back reaction, respectively. We might call them pseudo stoichiometric coefficients. For example, in the hypothetical elementary reaction 2A A + B, UfA = 2 and UbA = 1. Since this reaction is stoichiometrically equivalent to A B, the true stoichiometric coefficient is one. At equilibrium, W f = W b and U k f - Ukb = rk. Thus
-
-
By (AS) this is
or, in matrix form, V a ( l )= (j3/2)a(').Vp, which is (A4). In any normal case, there will be, of course, more than one elementary reaction. To adapt the foregoing development to that case, one need only add indices to the pseudo stoichiometric coefficients to indicate the reaction to which they refer and then sum over all elementary reactions. Thus, we have shown that our formulation of the problem is consistent with the equilibrium limit, a necessary condition for its validity.
