J . Phys. Chem. 1988,92, 3163-3173 in terms of a simple coding scheme. For the anisotropic Kepler problem (AKP) in two dimensions, infinite sequences of binary numbers are the natural scheme. The characteristic features and the physical properties of such a dynamical system seem all to be contained in the code, not in the usual surface of section. This principle applies in particular when one tries to understand how a harshly chaotic system is connected to its quantum-mechanical ancestor. The classical limit of Feynman's path integral and, in particular, the trace formula appear to be the only general guides available in such an endeavor. This approach works beautifully for the free motion on a surface of constant negative curvature, where Selberg's trace formula (which states that the two sides of (1 7) are strictly equal) establishes a direct relation between the spectrum of the Laplacian and the collection of all periodic orbits (some useful introductions are quoted in ref 3). The mathematical proofs in this case of hard chaos, however, take advantage of the high degree of symmetry in the underlying space
3163
and cannot be generalized. Nor can anything as strong as Selberg's trace formula be valid for the many harshly chaotic systems. The two sides of (17) are nevertheless expected to be very close quantitatively, with the singularities with respect to E on the right-hand side very near the poles on the left-hand side. Any possible deviations have to be small compared to the distance between the poles. In particular, the imaginary parts of the singularities on the right-hand side have to be small compared to the distance between them in the complex plane. The example of a mathematical model for the AKP in terms of spin chains with pairwise interactions shows that this last requirement imposes restrictions on the dynamical system which are not understood. They have to do with some delicate statistical properties of the periodic orbits which do not seem to hold for some values of the parameters in the spin-chain model. This mysterious feature in a harshly chaotic system can be called its third entropy, which measures something beyond its toplogical and its metric entropy.
Dynamics for Nonconservative Systems: Ergodicity beyond the Microcanonical Ensemble Julius Jellinek Chemistry Division, Argonne National Laboratory, Argonne. Illinois 60439 (Received: August 1 7 , 1987; In Final Form: December 3, 1987)
A generalization of the procedure suggested recently by Nos6 for dynamical simulation of a canonical ensemble is presented. It is shown that there are infinitely many different dynamics capable of mimicking the canonical ensemble for a physical system provided the corresponding dynamics for an extended system consisting of the physical system and a "bath" is ergodic in a particular sense. An invariance property of ergodicity under scaling of the time is established. It is shown how the time scaling can be used to "tune the dynamics into ergodicity" with respect to any desired distribution function in a given phase space. The dynamics of the extended system implies an apparent trajectory for the physical system which is quite different from the one generated by the Hamiltonian dynamics for a closed physical system. A general unified scheme for dynamical simulation of any statistical mechanical ensemble is formulated.
I. Introduction: Background and Overview Use of computers as a tool of scientific inquiry resulted in a remarkable extension of our practical abilities to study detailed dynamical as well as statistical properties of physical systems. The systems treated vary in size from an aggregation of a few particles, such as clusters with free or constrained boundary conditions,' to systems with enough degrees of freedom to mimic bulk matter.2 The two most important computational techniques for this purpose are those of molecular dynamics3 (MD) and Monte Carlo4 (MC). The M D procedure follows the time evolution of a system while the M C method simulates its ensemble behavior. In an attempt to keep the results of these two approaches in harmony, the fundamental issue of ergodicity, the relation between the time-averaged and ensemble-averaged properties of a system, emerges as a practical problem. The fundamental question in the general theory of dynamical systems can be formulated as follows5 (we consider only classical systems): does the dynamical system of interest generate an invariant measure, Le., a distribution function in the phase space, and if it does, what is this measure? (1) See, e&; Jellinek, J.; Beck, T. L.; Berry, R. S. J . Chem. Phys. 1986, 84,2783. Davis, H. L.; Jellinek, J.; Berry, R. S. Ibid. 1987,86,6456. Beck, T. L.; Jellinek, J.; Berry, R. S. Ibid. 1987, 87, 545 and references therein. (2) See, e&; Nost, S.;Yonezawa, F. J . Chem. Phys. 1986,84, 1803 and references therein. (3) Kusick, J.; Berne, B. In Statistical Mechanics, Part E Berne, B., Ed.; Plenum: New York, 1977; Vol. 6, Chapter 2. (4) Valleau, J. P.;Torrie, G. M . In Staristical Mechanics, Part A ; Berne, B., Ed.; Plenum: New York, 1977, Vol. 5, Chapters 4 and 5. (5) Eckmann, J.-P.; Ruelle, D. Rev. Mod. Phys. 1985, 57, 617.
0022-3654/88/2092-3 163$01.50/0
It is of at least as much interest to flip the question and to ask what types of dynamics (more explicitly, what dynamical equations of motion, if any) are implied by a given statistical ensemble or distribution. The reason is that one would like, ultimately, to be able to select from all those dynamics yielding the desired equilibrium properties a particular dynamics that will represent in a simple way the nonequilibrium properties of the system in interaction with a specific, complicated environment. For example, if one replaces the real heat bath coupled to a system by a very simplified model, one would like to know how to select the most accurate representation of that bath, in order to compute not only the macroscopic equilibrium properties but also time-dependent properties such as time-correlation functions and response functions. A conceivable way to establish a correlation between a dynamics and a statistical ensemble is to capitalize on the ergodicity property, either proven or assumed. If a description is restricted to the conventional phase space I'(ij,g),every point of which represents a mechanical state of the system, only two of the traditional statistical mechanical ensembles, the microcanonical and the canonical, can be defined through their phase-space density distributions. Other ensembles can also be represented by distribution functions but in an augmented phase space obtained from the conventional "mechanical" phase space by adding additional quantities, e.g., ~ o l u m eas , ~dynamical variables. (6) McQuarrie, D. A. Statistical Mechanics; Harper & Row: New York, 1976. (7) Andersen, H. C. J . Chem. Phys. 1980, 72, 2384
0 1988 American Chemical Society
3164
The Journal of Physical Chemistry, Vol. 92, No. 11, 1988
Originally the ergodic p r ~ p e r t yand ~ , ~probleml0-I3 were formulated for the microcanonical ensemble, Le., conservative systems under the assumption that energy is the only uniform integral of motion, apart from the fixed number of particles and volume. A way to generate a continuous deterministic dynamics for this ensemble is to use Hamilton’s equations (or their analogues). The ergodic property is then, in accordance with the von NeumannBirkhoff theorem,’*13 a consequence of the metric indecomposability of the whole of an energy surface (“strong erg~dicity”’~) or of part of it.Is Encouraged by the fact that the canonical distribution function is also defined in the mechanical phase space, one may attempt to extend the formulation of the ergodic problem by trying to define equations of motion that would mimic a canonical ensemble, Le., that would generate a particular Boltzmann distribution of probability density exp(-Ho(?j$)/k7+) E’( in the phase space, where Z(T) is the canonical partition function. In other words, one can ask whether it is possible to generate trajectories that would sample the accessible part of the phase space in accordance with a particular canonical distribution function. [For realistic physical systems the accessible part of the phase space is frequently only part of the whole phase space.] Of course, the equations of motion of a phase point in this case, if they exist, cannot be those of Hamilton. A number of suggestions were put forward by different aut h o r s , ’ ~ ’ ~regarding -~~ dynamics for a canonical ensemble-more strictly, an ensemble for which the number of particles N , the volume V, and the temperature Tare all fixed. These suggestions can be grouped in three classes. The common idea behind the papers of the first class1“22,25is to alter Hamilton’s equations of motion so as to satisfy the constant-T condition in a certain sense, preserving, however, the deterministic character of the system. The alterations were effected by introducing additional force terms or by rescaling the velocities or into equations of motion1s-20~22*2s momenta of the particle^;^^^^^^^^^^^ this last procedure results in the loss of continuity of trajectories. The characteristic of the second class is its use of a stochastic Langevin type of dynamics.2’,26-zs Finally, the third class can be viewed as a hybrid approach utilizing ideas from both aforementioned classes. Thus, Andersen’s p r o c e d ~ r ee.g., , ~ consists of propagating a trajectory for a short time interval using Hamilton’s equations: assigning
(8) Boltzmann, L. Wein Eer. 1868, 58, 517. (9) Gibbs, J. W. Elementary Principles in Statistical Mechanics; Yale University Press: New Haven, CT, 1902. (Collected Works; Yale University Press: New Haven, CT, 1948; Vol. 11.) (10) von Neumann, J. Proc. Natl. Acad. Sci. U.S.A. 1932, 18, 70, 263. (1 1) Birkhoff, G. D. Proc. Natl. Acad. Sci. U.S.A. 1931, 17, 650, 656. Birkhoff, G. D.; Koopman, B. 0. Ibid. 1932, 18, 279. (1 2) Khinchin, A. I. Mathematical Foundations of Statistical Mechanics; Dover: New York, 1949. (1 3) Miinster, A. Statistical Thermodynamics; Academic: New York. 1969. (14) Tisza, L.; Quay, P. M. Ann. Phys. ( N . Y . ) 1963, 25, 48. (15) In accordance with Khinchin’s assertion,’* the whole energy surface can never be metrically indecomposable. (16) Woodcock, L. V. Chem. Phys. Lett. 1971, 10, 257. (17) Abraham, F. F.; Koch, S. W.; Desai, R. C. Phys. Ret. Lett. 1982, 49, 923. (18) Hoover, W. G.; Ladd, A. J. C.; Moran, 9. Phys. Reu. Lett. 1982, 48, 1818. Ladd, A. J. C.; Hoover, W. G. Phys. Reu. E : Condens. Matter. 1983, 28, 1756. (19) Evans, D. J.; J . Chem. Phys. 1983, 78, 3297. Evans, D. J.; Morris, G. P. Chem. Phys. 1983, 77, 63. Phys. Lett. A 1983, 98A, 433. (20) Evans, D. J.; Morris, G. P. Comp. Phys. Rep. 1984, 1, 297. Evans, D. J.; Hoover, W. G. Annu. Rev. Fluid Mech. 1986, 18, 243. (21) Heyes, D. M. Chem. Phys. 1983, 82, 285. (22) Haile, J. M.; Gupta, S . J . Chem. Phys. 1983, 79, 3067. (23) Evans, D. J.; Morris, G. P. J . Chem. Phys. 1984, 81, 3749. (24) Haile, J. M. J . Chem. Phys. 1984, 81, 3750. (25) Brown, D.; Clarke, J. H. R. Mol. Phys. 1984, 51, 1243. (26) Schneider, T.; Stoll, E. Phys. Rev. E : Condens. Matter 1978, 17, 1302. 1978, 188 6468. (27) Hiwatari, Y.; Stoll, E.; Schneider, T. J . Chem. Phys. 1978,68 3401. (28) Berendsen, H. J. C.; Postma, J. P. M.; van Gunsteren, W. F.; DiNola, A,; Haak, J. R. J . Chem. Phys. 1984, 81, 3684. (29) NosC, S . J . Chem. Phys. 1984, 81, 511. Mol. Phys. 1984. 52. 255.
Jellinek at the end of this interval a momentum randomly chosen from a Boltzmann distribution of momenta to one, also randomly selected, particle of the system; and sequentially repeating the deterministic propagation and stochastic momentum-assignment steps. The comparative advantages and drawbacks of the different M D schemes for a constant-T ensemble are extensively discussed in the so we shall not dwell on those here. In the context of the ergodic problem, the primary interest is, of course, in studies belonging to the first class. The most satisfactory dynamical scheme of this class is, indubitably, the one suggested recently by In fact, all the other continuous deterministic dynamics for constant temperature can be obtained from it by imposing additional constraints. The two essential points of Nose’s method are (a) expressing the physical momenta of a system as “virtual” momenta scaled by an additional degree of freedom s, which is introduced to represent a thermal bath, and (b) choosing a particular form of the Hamiltonian for the extended system consisting of the physical system and the ‘.bath“. The system with virtual dynamical variables will be called the “virtual system”. The remarkable feature of Nose’s Hamiltonian is that it produces a dynamics which generates a canonical distribution in the phase space of the physical system provided this dynamics is ergodic in the usual microcanonical sense in the extended (nonphysical) phase space of the virtual system and the “bath”.z9~30 This is in contrast to all other deterministic methods heretofore suggested. It seems to be generally accepted that only the particular Hamiltonian of Nose possesses the outstanding property just mentioned. This Hamiltonian defines a specific embedding of a surface of constant energy in the extended phase space of the virtual system and the “bath”. The Hamiltonian dynamics on this surface translates into a specific non-Hamiltonian evolution of the physical system. Both the Hamiltonian dynamics and the non-Hamiltonian evolution depend parametrically on the energy E of the extended system and the “mass” Q of the single degree of freedom s representing the “bath”. But to the extent these parameters are specified, the dynamics are defined uniquely. As a consequence, one might get the impression that there is a specific correlation between a canonical ensemble and a dynamics in the extended space, namely, that established by NosC’s Hamiltonian. One of the purposes of this paper is to analyze in detail the problem of dynamical simulation of a canonical ensemble and to introduce a general solution to this problem. It turns out that even if the parameters E and Q in Nose’s Hamiltonian are specified, that Hamiltonian still can lead to many different dynamics. Moreover, there exist infinitely many other Hamiltonians with all the properties of Nos&’sHamiltonian but which generate different dynamics. In fact, almost any Hamiltonian for the extended system can, in principle, lead to a dynamics capable of mimicking the canonical ensemble for the physical system. We shall also give a unified prescription for dynamical simulation of any statistical ensemble in which the number of particles is conserved. In the next section we rephrase tersely the essence of NosC’s approach. In section I11 we derive the different dynamics generated by Nosi’s Hamiltonian, consider the possible relations between these dynamics, and clarify the implications of the scalings and weightings involved. A general scheme for M D simulation of a canonical ensemble is formulated in section IV. In section V we consider the probtem of ergodicity and analyze the structure of a dynamics ergodic with respect to a canonical ensemble. A general unified scheme for dynamical simulation of any statistical ensemble is formulated in section VI. A brief summary is presented in section VII. 11. The Nose Scheme
Following NOSE,^^ we use primed coordinates ij’ 5 ( i j l q and momenta j ’ = for a physical system; unprimed 4‘ {glj and j = E,] are reserved for a virtual system (i enumerates the particles in the physical and virtual systems). We also consider two ex(30) Hoover, W. G. Phys. Rea. A 1985, 31, 1695.
The Journal of Physical Chemistry, Vol. 92, No. 11, 1988 3165
Dynamics for Nonconservative Systems
tended systems. One is that of a physical system and a “bath”; the other consists of a virtual system and a “bath”. The “bath” is represented by a single dimensionless coordinate s and its conjugate momentum ps. In this notation (qi’,fii$(iji,fii),(&’,fii’,s,ps), and (&,jji,s,ps)are the phase spaces of the physical system, the virtual system, and the two corresponding extended systems, respectively. The correlation between a physical and a virtual system, and consequently between the two extended systems, is established through the relations
&’ = qi,
3i’= @,/s
(1)
The dynamics is introduced into the extended phase space (q‘,,fii,s,psJ, in which qi and j j i and also s and p s are considered as canonically conjugate dynamical variables, by postulating Nost’s Hamiltonian:
Ps2 + 4({iji))+ - + gkT In s 2m,’s2
H=C- p’? I
2Q
(2)
where m,’are the masses of the particles in the physical system (mi = mi’s2can be viewed3*as a varying “virtual” mass); 4 is the potential energy for both the virtual and physical systems; Q is a “mass” associated with the s coordinate; k is Boltzmann’s constant; Tis an external parameter viewed as temperature: and, finally, g is a dimensionless parameter, the value of which depends in Nosi’s f o r m ~ l a t i o non~ ~the number of particles in the virtual (equivalently, physical) system and on the particular choice of time for trajectory propagation. [In section I11 we show that in fact any real finite nonzero value of g is admissible in ( 2 ) . Through appropriate time scaling the generated dynamics can, in principle, always be forced to simulate the canonical ensemble for the physical system.] Hamilton’s equations of motion implied by the form ( 2 ) , which generate trajectories in the space (qi,p’i,s,ps), are
conserve the value of H’, and thus the corresponding trajectories are confined to an energy surface in the extended space (lji’Ji’,s,ps). Nost’s original non-Hamiltonian equations generate trajectories in the phase space (&’&’,s$~’ = p,/s). The advantage of using p s (and not ps’= p , / s ) in conjunction with the scaling dt’ = dt/s is, as noticed also by Evans and H ~ l i a nthat , ~ ~(4a), (4b), and (4d) form a closed system of equations, and one may avoid solving (4c) while studying the time evolution of the physical system described by iji’ and fii’. Hoover30 approached Nosi’s equations of motion from a somewhat different point of view. It is easy to see, however, that Hoover’s form of these equations [(6) of ref 301 is identical with (4a), (4b), and (4d).33 The relevance of the dynamics discussed above for the canonical ensemble becomes clear from the following property of NOSE’S Hamiltonian ( 2 ) : there exists a value of g such that the (SIweighted, I finite and real; see section 111) microcanonical partition function Z,,calculated on any energy surface spanned in the extended phase space (&,fii,s,ps)by the equation ff(q‘,P’,s,p,)- E = 0
(6)
is proportional to the canonical partition function Zc
In (7) we omitted the constants unimportant in the present context and
is the energy of the physical system. Thus ZIgt$iJJ’sI = C(E,T,Q,g)@j’*@t1 PC
(9)
In Nosi’s original formulation the value of g which leads to equality 9 is 3N 1, if the unscaled time t is used, or 3N for the scaled time t’. The proportionality constant C depends on the and {iji’,p^i? parameters E , T , Q, and g. The superscripts (&,3i,s,ps) in (9) indicate the phase spaces in which the corresponding partition functions are calculated. It follows from (9) that for any phase function A(q’,p”) = A(ij,jj/s) defined on the physical phase space (&’,fiij
+
( A( i t , $ 3) p 3 , 1 = ( A(q‘,p’/s) ) p t
Using relations (1) and defining a scaled time t’-through the regular uniform time t as dt’ = dt/s, we can translate Hamilton’s equations (3) into non-Hamiltonian equations in the extended phase space (iji’,fii’,s,ps}: dqi’ - = - di’ dt’ mi’
(4a)
(10)
where (AT stands for the corresponding ensemble average of the quantity A ; the specifications of the ensembles and the phase spaces are given by the subscripts and superscripts, respectively. Thus, canonical ensemble averages of physical quantities can be calculated as their [weighted; see section 1111 microcanonical averages in the nonphysical extended phase space (C?,,fii,s,ps).The latter, in their turn, can be obtained, by invoking the ergodicity hypothesis, as dynamical time averages along the trajectories generated by Hamilton’s equations (3) or non-Hamiltonian equations of the type (4): 1
(A(ij,jj/~))f~@ = ~Tlim -~ m ~JTl’T~0Al ( g ( t ) , F ( t ) / s ( t ) )dt
(11)
[Nosi’s prescription is to use the value g = 3 N + 1 in (3d) and g = 3 N in (4d)l. This establishes the relation between the dynamics generated by Nose’s Hamilionian (2) and the canonical ensemble for the physical system. The Hamiltonian ( 2 ) has the following form in the space fi;,s,psl:
(q;,
111. Dynamics for NOS&’S Hamiltonian. Implications of the Scalings and Weightings In the previous section Hamilton’s equations of motion (3)in the space (&,p’i,s,ps)were translated into non-Hamiltonian equations (4) in the space (q‘i’,j5i’,s,ps).The scaled time t’was used in the
Although equations of motion (4) are not Hamiltonian, they (31) Ray, J. R.; Rahman, A. J . Chem. Phys. 1985, 82, 4243
(32) Evans, D. J.; Holian, B. L. J . Chem. Phys. 1985, 83, 4069. (33) In contrast to what is claimed by Evans and Holian,’* Hoover3o used in effect the variable p s and not ps’ = p J s ,
3166
The Journal of Physical Chemistry, Vol. 92, No. 11, 1988
space {ijll,jj,’,s,p,). Nos6 wrote the non-Hamiltonian equations (his employing equations in “real ~ a r i a b l e s ” in ~ ~the ) space {~l’,jjlr,s,ps~ an additional scaling p,’ = p s / s . As a motivation for transforming to the extended phase space of the physical system and the bath, one can cite the desire to propagate directly the physical coordinates &’and momenta PI’. In order to understand the implications of the scalings for the dynamics and the phase-space distribution, we shall consider separately general scalings of p s and of time. Consider first a general scaling of p , by a real nonvanishing differentiable function a ( s ) Pf’
= P,/.(S)
(12)
Jellinek [To avoid cumbersome notation, dynamical variables held constant in the partial differentiation in (16)-( 19) have been partially omitted. However, it is important to remember to hold primed variables constant when differentiating with respect to a primed variable and to hold unprimed variables constant when differentiating with respect to an unprimed variable.] Now taking into account (3), we conclude that the rhs of (20) is exactly the Liouville operator L in the space {q‘i,jji,s,ps}.Thus L’= L
(21)
It can readily be verified that equality 21 remains valid also for L’generated by equations of the type of_(L3)in the space (iji,$i,s,p,r = p,/a(s)) or in any general space {Qi,Pi,s,Ps), where
and write the non-Hamiltonian equations of motion in the space {ijl’,j?l’,s,~s~ using the unscaled time t. Equations for g,’,P,’, and s are obtained from (4a), (4b), and (4c), respectively, through dividing the latter by s and using the relations (12) and s dt’= dt: dq‘l’ -=dt
PI’
(13a)
m,’s
andfix(s), hix(s),and u(s) are real nonvanishing differentiable functions of s; i = 1, ..., N, X = x, y , z. It follows from (21) that eiL’f = eiLI (23) The remaining equation for ps’ takes the form dp,‘ - d(.P,/a(s)) 1 dp, P, d4s) ds _ --- dt
dt
a ( s ) dt
.2($)
ds
dt
(14)
which, taking into account ( l ) , (3c), (3d), and (12), can be rewritten as
The Liouville operator L’ generated by (1 3) is
Using (1) and (1 2) and relations
($),
=
(
s) ;
we can rewrite L’defined by (15) in terms of q, P, s, and p,. The result is
Le., the time evolution of any phase function A ( i j r ( t ) g f ( t ) )= A ( @ ( t ) g ( t ) / s ( t )is) the same in all the different phase spaces considered and is uniquely defined by Hamilton’s equations (3). This result is, of course, a trivial consequence of the fact that the non-Hamiltonian dynamics in the different phase spaces are by construction simply images of the Hamiltonian dynamics in the space (iji,$i,s,ps].We gave an explicit proof of equality 23 in order to stress two points: first, that there is nothing unique about spaces {iji’,&’,.~,p,’= p S / s )and {&’,@i’,s,ps] considered by N o ~ i ;Hoover,3o ~ and Evans and H ~ l i a nand ; ~ ~second, that the possible scalings of p,, (12), are not correlated with the transformations (1) at all. It is important, however, to realize that although the time evolution of any phase function is indifferent to the choice of the phase space, the different phase spaces-are not equivalent. The space {~,,ji,s,ps) is special because, by our choice, only in this space are the coordinates and momenta of the extended system canonically conjugate d p m i c a l variables. Selection of a Hamiltonian H(&j,s,p,) then defines the dynamics in this space uniquely. By contrast, knowledge of the Hamiltonian in terms of the (primed) variables of any other phase space is insufficient to determine the equations of motion; because the primed variables are not canonically conjugate, one must also know the transformation from the primed to the unprimed variables. We stress again that equality 23 is obtained assuming that the same unscaled time t is used in all the phase spaces. In earlier work,29,30*32 however, the non-Hamiltonian equations of motion were always written employing ihe scaled time t ’defined through the differential relation dt’ = d$s. An attempt was made29to relate this particular scaling of time to the transformation (1) of the momenta. This might leave the impression that the time scaling is defined by the transformation of the phase space or vice versa. In what follows we shall show that the time scaling is a totally independent procedure and is not related to or defined by the correlation between the physical and virtual systems, (l), or the transformations of the phase space (iji,Pi,s,ps), (22). That the time scaling has a nontrivial effect on the phase-space density and the corresponding partition function was noticed already by NosiZ9and discussed by Evans and H 0 1 i a n . ~Here ~ we shall show that the heretofore unrecognized but important implication of time scaling in the phase space (iji,&,s,p,)is the weighting of the microcanonical distribution generated by the Hamiltonian unscaled-time dynamics; that the latter generates a microeanonical distribution, at least in subregions of an energy shell in the space (q‘i,~i,s,ps), is guaranteed by the von Neumann-Birkhoff t h e ~ r e m . ’ ~A’ quite ~ remarkable consequence
The Journal of Physical Chemistry, Vol. 92, No. 11, 1988 3167
Dynamics for Nonconservative Systems of this weighting is an invariance feature of the ergodicity property of a dynamics under the time-scaling transformation. To formulate and prove this invariance, we consider a general time scaling dt’ = dt/@(s)
furnish canonical averages of physical quantities. Indeed, following Nosb’s procedure,29 we can write
(24)
where @(s) is a real continuous nonvanishing function. The t’ time average A of a phase function A(ij’(t’),$’(t’)) = A(ij(t’),jj(t’)/s(t’)) can be represented as 1 T’ A E lim A(ij(t’),j?(t’)/s(t’)) dt’= T’-- T’ O
=
1
s3N+1
gkT
@(SI
- S d 8 ’ I d a ’ s d p s I d s -X
-s
In obtaining the rhs of (29) we used ( l ) , (2), and (8) and the relation Invoking ergodicity in the usual [microcanonical] sense in the space (iji,j5i,~,p5} we can rewrite (25) as
A=
(A(q‘,jj/S)/@(S))lj~~i,5”.)
( 1 /@(s))ljJ?”,.”P.J
where so is the solution of the equationf(s) = 0. Choosing @(s) as a power of s,
--
I real
(31)
g = 3N + 1 - 1
(32)
@(s) = SI,
1
- JdP’ ZFC
Jdq‘ JdP,
JdS
6(H@’?7,p5,s) - E ) @(s) ( A (q‘,jj/s) )i$JJ’tl
-
and defining g as we deduce from (29)
(26)
The subscript wfic in the rhs of (26) stands for the ”weighted microcanonical”. Note that although (25) and (26) are written for a special dependence of A on the phase variables, they are valid for any function defined on (iji,jji,sLps). The equality of the time average A and of the ensemble aveFage ( A ) w Finc (26) can be viewed as a definition of ergodicity in the extended, weighted microcanonical sense. The wfic ergodicity implies equality between the @(s)-scaled-timeaverage and the phasespace average of any function defined on the space {iji,jji,s,p,} in which the [unnormalized] density distribution [d.d.] is given by a @(s)-weightedmicrocanonical density:
Z$$tJPxl = C(E,T,Q,N,I) s d j j ’ f d i j f e-[HoU’,a‘)/kTl =
c(E,T,Q,N,[)z~B~’,I?I~(33) where the proportionality constant (34) Thus, for any real 1 the weighted microcanonical density distribution d.d. =
+ 6((iji})+ Ps2 - + (3N + 1 - 1)kT In s - E 2Q
S‘
It follows from (25) and (26) that if an unscaled-time dynamics is ergodic in the usual microcanonical sense, the corresponding scaled-time dynamics is ergodic in the weighted microcanonical sense. Repeating the steps taken in deriving (25) and (26), one can easily prove that the converse statement is also true; Le., ergodicity in the weighted microcanonical sense when the scaled time is used implies ergodicity in the microcanonical sense for the unscaled time dynamics. We shall refer to this invariance property of ergodicity as its scaling invariance. It should be noted, however, that (A(q,p’/S))~;$Ss.J # (A(q‘,jj/s))ljJ?Pn5.P.l
(28)
unless @(s) = 1. The last inequality reflects the fact that if g in (2) is chosen so that the microcanonical distribution function 6(H(ij,jj,s,ps) - E ) , defined in the space (ijl,jJi,s,ps),generates a canonical distribution exp(-Ho(ij’,p”)/kT) in the physical phase space (ijif,Bij,Le., g = 3 N + 1 , then the weighted microcanonical distribution function, (27), with the same value of g, will not generate a canonical distribution in the space (ij,’,P,j. It follows from the scaling invariance of ergodicity and from inequality 28 that if a dynamics is known to be ergodic in eithei the unscaled-time or scaled-time sense, the unscaled-time average of any physical quantity cannot be equal to its corresponding scaled-time average [the word corresponding is used here to stress that the same value of g is used in both of the dynamics]. Nevertheless, inequality 28 does not preclude the possibility that a weighted microcanonical distribution, (27), might generate a canonical distribution in the physical phase space. Suitable weightings correlated with appropriate choice of the value of parameter g will guarantee that the distribution function (27) will
(35) generates a canonical distribution in the physical phase space. The average value of any physical quantity A(q‘$/s) = A(q‘’,p”)calculated for distribution (35) with any real l will have the same value and will be equal to the canonical ensemble average of this quantity: (A(ij,P’/s))~;{yP.l = (A(q‘’,P”? )p!”
(36)
The average on the Ihs of (36) is calculated for the distribution function (35), which depends parametrically on I. This distribution function confines the extended nonphysical system to an energy surface in the space (iji,jjlrs,ps],but the principle of equal a priori probability is not satisfied by it. This, however, does not cause any conceptual difficulties since the points on-this surface do not represent states of a closed physical system. We stress now that one may select-any real value of g in Hamiltonian (2). Representing this value through (32) we obtain the unscaled-time Hamilton’s equations of motion:
3168 The Journal of Physical Chemistry, Vol. 92, No. 11, 1988 If the dynamics defined by (37) is ergodic, then [almost] every trajectory will generate in the space {q‘&,s,pS]the microcanonical density distribution d.d. =
+ t$({;J)
SZ + P+ (31%’ + 1 - I)kT In s - E 2Q
(38) which does not lead, unless 1 = 0, to a canonical distribution in the physical phase space (;i’,$ij. In other words, the time averages of physical quantities calculated along trajectories generated by (37) will not be equal to their canonical ensemble averages. Introduce now the scaled time t ’ as d t ’ = dt/s‘
(39)
and rewrite (37) in terms oft’:
ds = dt’
P2 Q
If (37) generate the microcanonical distribution (38), then in accord with the scaling invariance of ergodicity the dynamics defined by (40) will generate the w ~ distribution c (35). But the latter, as proved above, produces canonical ensemble averages of physical quantities A(q“,j’)= A(;,$/s). Thus time averages of A(q(t’),fi(t’)/s(t’)) calculated along trajectories generated by (40) will be equal to their canonical ensemble averages. We have shown here that any real value of g is admissible in the Hamiltonian (2); for any such value parametrized through (32) there exists a time scaling, (39), such that the non-Hamiltonian scaled-time dynamics, (40), mimics the canonical ensemble behavior of the physical system, provided the Hamiltonian dynamics of the extended system, (37), is ergodic in the usual microcanonical sense. Although (40) are n~n-Hamiltonian,’~ they generate trajectories on the same energy surface as the corresponding Hamilton’s equations (37); the correspondence is established by the same value of I. Indeed dH(q‘,j,s,p,) = dt ’
dff(&p’,s,ps) =0 dt
Jellinek generating the wwc distributions, (35). Since, as shown above, the latter lead for all values of 1 to the same canonical distribution in the phase space of the physical system, the different I-labeled dynamics, (40), which generate the corresponding I-dependent density distributions ( 3 9 , all furnish the same time-averaged value of any physical quantity A(ij(t?,j?(t’)/s(t?).This time average of A(;’,*’) is equal to its canonical ensemble average. Thus many different dynamics may simulate the same canonical ensemble behavior of a physical system. This characteristic suggests an indicative test for ergodicity of dynamics with respect to a canonical ensemble. Occurrence of two or more different dynamics, (40), corresponding to different values of I , which produce the same time-averaged value of a physical quantity A is an indication of ergodicity of each of these dynamics; the likelihood that the time averages are equal by accident is very small. If different dynamics are ergodic with respect to a specified canonical ensemble in the physical phase space, the corresponding time averages for any physical quantity will be equal. Thus,. occurrence of different dynamics which produce the same time-averaged values for more than one physical quantity may be interpreted as a further support for ergodicity of each of the dynamics. The practical application of this test for ergodicity may be complicated by questions such as how long to run the trajectories and how close the values of the averages should be to be considered equal. The usual wisdom is to run the trajectories until convergence is achieved within a predefined accuracy [say, O.l%] and to consider the calculated averages equal if they coincide within this accuracy. Defining more precise criteria of convergence is a nontrivial important task. In order to establish a general correlation between the corresponding unscaled-time ( t ) and scaled-time (t’) dynamics for a given value of g in Hamiltonian (2) we consider again a general time sc;iling, (24). Equations 3 rewritten in terms of the scaled time t’ take the form dGi _ -dt ’
$i
@(SI
--
mi’s2
It follows immediately from (3) and (42) that
L = P(S)L (41)
It is trivial to see that this result remains valid for any general time scaling, (24). O_ne_can rewrite (40) in terms of any other extended phase space {Ql,C,s,P,}, (22). This, as shown above, does not alter the dynamics. Different values of g in Hamiltonian (2) imply different dynamics. More explicitly, (37) generate different trajectories for different values of I; see (32). Some, all, or none of the dynamics corresponding to different I are ergodic;35those that are generate the distribution functions (38). Due to the scaling invariance of ergodicity those values of 1 for which (37) define evolutions ergodic in the microcanonical sense will define scaled-time dynamics, (40), (34) By including t and -H into phase space as conjugate ‘coordinate” and *momentum”, one can recast the equations into Hamiltonian form. See, e.g., Szebehely,. V. Theory of Orbirs; Academic: New York, 1967. (35) We call a Hamiltonian dynamics nonergodic if the accessible part of an energy shell cannot be covered completely with a single trajectory. In accord with the von Neumann-Birkhoff theorem, however, even in the worst case a complete coverage can be achieved by selecting a large enough number of different initial conditions. ( 3 h ) Posch. H 4.: Hoover. W G : Vesely, F. J . Phys. Rei.. A 1986. 23. 4’5 >
where
is the scaled-time Liouville operator and L is the unscaled-time Liouvillian defined by the rhs of (20). We represent now the scaled-time propagator exp(ilt9 of a phase point (q‘,,@,,s,p,)in the form32 exp(ilt’) = expR( L‘iLd#)
(45)
where expR is the time-ordered exponential, and use (43) and the general relation (24) between the unscaled time 7 and scaled time 7‘ to obtain: exp(iLt’) = exp(ilt)
(46)
Thus for any scaling of the time, the corresponding [to the same value of g] unscaled-time and scaled-time dynamics not only evolve on the same energy surface but they generate the same trajectories. The only effect of time scaling on dynamics is that a fixed unscaled-time interval At becomes a varying scaled-time interval At’:
The Journal of Physical Chemistry, Vol. 92, No. 11. 1988 3169
Dynamics for Nonconservative Systems (47)
A consequence of this is that the corresponding unscaled-time and scaled-time averages of any phase function cannot, in general, be equal. It is in this sense that the scaling of the time produces a new dynamics. This last statement about the inequality of the corresponding unscaled-time and scaled-time averages is a general one. It is in fact an extension of a similar but weaker conclusion derived above for ergodic dynamics by using scaling invariance of ergodicity and inequality (28). This general conclusion does not mean that either the unscaled-time or the scaled-time dynamics cannot simulate a canonical ensemble for the physical system. It means only that the different-time corresponding dynamics cannot generate the canonical distribution simultaneously. At most one of a set of corresponding dynamics can actually generate the canonical ensemble for a given temperature in the physical phase space. We reemphasize now that equality 46 is valid for any time scalings and any phase spaces, including the cases when the t-time and r’-time corresponding dynamics are written in different phase spaces [cf. (21)]. Evans and H 0 1 i a n ~(EH) ~ obtained (46) for the case dt’= dt/s [i.e., P(s) = s in (43)J and considering the phase for the time r and (g,’ = q‘,Jf’ = jjf/s,s,ps)for the spaces (q‘l,jjl,s,ps] time t’. From it they concluded that NOSE’Sdynamics, (3), and Hoover’s dynamics, in our notations (4), generate trajectories which follow the same paths [meaning, of course, after transformation to the same phase space; note that (4) take in the space (q‘,,p’,,s,pJthe form of (40) with I = 11. This is, however, not so. In order for the E H conclusion to be true, the two dynamics would have, in our terminology, to correspond; Le., they would have to have the same value of g . But then at least one of them would not generate a canonical distribution for the physical system. Nosb’s p r e ~ c r i p t i o nis~to ~ use different values of g in (3d) [ g = 3 N + 11 and in (4d) [g = 3N, the same value g = 3Nis obtained using (40d) with I = 11. With these different values of g (3) and (4) generate trajectories that are not related by (21), (43), and (46), and thus these trajectories do not follow the same paths. Each of the two dynamics is, however, capable [at least in principle] of producing can6nical averages for physical quantities. The main conclusion following from the considerations presented above is that Nose’s Hamiltonian (2) with g parametrized as in (32), the distribution functions (35), and the time scalings (39) suggest a wide variety of I-labeled dynamics, (40), for simulating a canonical ensemble. These dy_narflicscan be expressed in terms of any generalized pfiase space (Q,,P,,s,P,),(22). NOSE’Sdynamics, (3), and Hoover’s dynamics, (4a), (4b), and (4d), are just particular cases with I = 0 and I = 1, respectively. The ergodic properties of the different dynamics depend on the parameters Q and E . Earlier numerical tests29as well as our own show that the static quantities may remain insensitive even to large variations in the value of Q. Nosez9observed only weak sensitivity to the value of Q over a wide range even for a dynamical quantity-the velocity autocorrelation function. [Extremely small or large values of Q do have noticeable effects, however.] The requirement of ergodicity of a dynamics with respect to the corresponding distribution function (35) in the space (q‘f,jjl,s,psp,) is a sufficient condition for this dynamics to produce time averages of physical quantities equal to their canonical ensemble averages. In the next section we show that there are infinitely many other distribution functions leading to equalities of the type (9) and (10). Each of these distribution functions implies a different extended phase space dynamics furnishing the same canonical ensemble averages for the physical system. The distribution function generated by these dynamics in the physical phase space is in fact somewhat different from that for a true canonical ensemble.29 This is because the dynamics conserve the total virtual momentum Xp’, for both the unscaled [cf. (3b)l and scaled [cf. (42b)l times; for free boundary conditions also the total virtual angular momentum E(;, X jj,) will be conserved. The net effect of these conservation laws on the physical system is the imposition of additional constraints and the reduction of the number of degrees of freedom.
This should be taken into account in expressing the total kinetic energy of the physical system in terms of k T / 2 .
IV. Other Hamiltonians In the previous section we have shown that equalities of the type (9) and (10) can be obtained for Hamiltonian (2) with any nonzero real values of g. This follows immediately from the distribution function (35) and (32). Different values of g define different dynamics, (40), [cf. (32)], each of which is capable of mimicking the canonical ensemble for the physical system. Here we show that even the requirement for the s-dependent “potential energy” term to have the functional form In s can be removed. In fact, by appropriately choosing an s-weighting of the microcanonical distribution in the extended phase space (q‘,,jjl,s,p,},one can obtain equalities of the type of (9) and (10) for any smooth Hamiltonian of the extended system lthe requirement of smoothness is specified below]. Conversely, for any fixed sweighting of the microcanonical distribution function there are infinitely many different Hamiltonians which lead to equalities of the type (9) and (10). To prove these statements, consider a general weighted microcanonical distribution function
where w(s) is a nonvanishing continuous function. Consider also a general Hamiltonian pixz Piy2 ++
hi:($)
-> 2
hj,2(~) hjZ2(s)
+
wherefiA(s),hih(s),u(s), and u(s) are nonvanishing differentiable functions; i = I , ..., N, X = x, y , z. The weighted microcanonical ,, in the space (q‘,,jj,,s,p,p,) is defined by partition function ,Z
where H(q‘,p’,s,p,)is given by (49). Introducing primed variables in accordance with
we can rewrite (50) in the space {q‘i’,jjif,s,psl:
H’ in the last equation is
Assuming that the equation H’(“’,”’,s,p,?
-E
=0
(54)
has a single real solution so with respect to s $0
= u-’
(
iji’2
[E
- Ei-,2mi - $J((&?)-
where u-’ is the function inverse to L’, and using property (30) of the 6 function,. we obtain from (52)
3170
The Journal of Physical Chemistry, Vol. 92, No. 11, 1988
Jellinek
For illustration consider the Hamiltonian - 2
used by Nos&29as a specific example of the case in which the equalities 59 and 60 allegedly cannot be obtained. The wwc partition function in this case is
where
and u’(s) is the derivative of the function u ( s ) [we assume that u’(so) # 01. The requirement
where F(E,Q,T;p,’) is an integrable function of p:, together with (50)-(52) and ( 5 6 ) then immediately leads to the equality ZIB,B,.s.P,I wllc = C(E,Q,T)ZLqi’*@l’’
(59)
6(s
- [ E - Ho(P”Z? - PSZ/~QI / & T )
(62) w(s) The rhs of (62) is arrived at by using (l), (8), and (30). Choosing w(s) as
Interpreting q’ and p” as the physical coordinates and momenta, respectively, we obtain from (59) that for any physical quantity w(s) = (63) A(G’29 = A(&,(s)qi,fiy(s)qiy~z(s)qi:~, ~ i ~ / h i ~ ( s ) , p i ? / h i ~ ( s ) , p i : / we immediately obtain from (62) the equality 59 and thus the hil(s))), an equality holds: equality 60 for any value of g . The scaled-time non-Hamiltonian (/g)yi‘AI = (A)IB,*bi,S,P,i WPC (60) equations of motion for simulating a canonical ensemble for the physical system take in the space (iji,jji,s,ps)the form Thus the special property possessed by NOSE‘SHamiltonian ( 2 ) dqi s3N-2e-g~ gi is also valid for any Hamiltonian (49) [and distribution function - =-(48)] for which the requirement (58) and the assumption about (64a) dt ’ &? mir singleness of the real solution of (54) are satisfied. It is obvious that one can make an infinite number of different selections of the set of functions (hi,(s),hiy(s),hi,(s)), &,(s), fy(s)fiz(s)), u(s), and u(s) such that each set will be a solution of the functional equation (58) with a predefined form of the function w(s). This means that for any selected weighting in distribution function (48), defined by fixing the function w ( s ) , there is an infinite number of different Hamiltonians (49) such that each of them leads to the equalities 59 and 60. Each of these Hamiltonians defines a different dynamics [more precisely, many different dynamics for different time scalings-see the previous These equations can easily be rewritten in terms of physical cose_ction]in the space (iji,j3iis,ps)or in any other generalized space ordinates and momenta. (Qi,Pi,s,Ps) defined by (22); thef(s), h(s), and u(s) functions in Note, finally, that in Hamiltonian (49) not only the physical (22) are, in general, unrelated to those in (51). The different momenta but also the physical coordinates are expressed as their Hamiltonians define different embeddings of a surface of the same virtual counterparts scaled by the corresponding functions of s. constant energy in the phase space. The dynamics defined by those This has no effect whatsoever on the computational procedure different Hamiltonians (49) which satisfy the requirement (58) if free boundary conditions are used. In the case of periodic can be written in terms of the scaled time t’given by (24) with boundary conditions the only p i n t to be taken into account is that to a fixed computational cell in the physical configuration space P(s) = w(s). Each of the thus obtained t’-time dynamics will furnish the same time averages of physical quantities equal to their 1there corresponds a varying computational cell in the space canonical ensemble averages, provided the corresponding uns{GI] of the virtual coordinates. caled-time dynamics are ergodic in the usual microcanonical sense. V. Ergodicity and Structure of the Dynamics in the Phase These different dynamics can be used for indicative test of erof the Physical System Space godicity. Consider now the case of an arbitrary fixed Hamiltonian (49). Ideally one would like to find necessary and sufficient conditions Since for any givenfx(s), hix(s),u(s), and u ( s ) [i = 1, ..., N, X for a dynamics to be ergodic with respect to a canonical ensemble. = x, y , z ] one can always find a function w(s) satisfying ( 5 8 ) , In other words, one would like to have an analogue of the von we conclude that for any Hamiltonian (49) there is a distribution Neumann-Birkhoff theorem’&l3 which gives these conditions on function (48) that generates the equalities 59 and 60. A consean energy surface in the phase space. This problem is, however, quence of this is that for any Hamiltonian (59) there is a time beyond the scope of the present study. Instead, we shall analyze briefly the implications of the dynamics considered in the previous scaling defined by (24) with p(s) = w(s), where w(s) is the solution of the functional equation (58) for givenAh(s), hiA(s),u(s), and sections for the phase space of the physical system. u(s), such that the scaled-time non-Hamiltonian dynamics genThe ergodicity of a dynamics with respect to the corresponding erates the canonical distribution for the physical system if the weighted microcanonical ensemble in the extended phase space corresponding unscaled-time Hamiltonian dynamics of the ex{~l,jjl,s,,ps), which establishes its ergodicity with respect to the tended system is ergodic in the usual microcanonical sense. Thus canonical ensemble for the physical system, neither follows from nor is made more plausible by equalities 59 and 60. This property scaling of the time can be used to “tune” in principle any dynamics should be imposed as an additional requirement. To illustrate into ergodicity with respect to the canonical ensemble. In fact, as. we show in section VI, the time scaling can be used to tune this, consider the Hamiltonian a dynamics into ergodicity with respect to any ensemble, i.e., FZ PSZ distribution function, defined on the accessible part of the physical Wq‘$,s,p,) = XL, + qfJ((s‘ll) + - + k T In s (65) phase space. I 2m, 2Q
The Journal of Physical Chemistry, Vol. 92, No. 11, 1988 3171
Dynamics for Nonconservative Systems and the microcanonical partition function
Ps2
4((&}) + - + kT In s - E 2Q Assuming that the physical and the virtual systems are the same, i.e.
gi’ = 4i
-
(67a)
ti’= $i
(67b)
and applying to (66) the procedure used in the previous sections, we-arrive at equality 59. The Hamilton’s equations generated however, completely decouple into those by the Hamiltonian (a), for the virtual [equivalently, physical] system:
-dGi = - P’i dt mi’ d$i - = - -a+ dt a& and those for the “bath”: (69a j
Each of the two decoupled systems of equations (68) and (69) separately conserves the total energy of the virtual [physical] system and the “bath”, respectively, and thus the underlying dynamics on an energy surface in the extended phase space {&$,p,ps]is not ergodic. This example underscores the importance of the indicative test for ergodicity discussed in section 111. One could, of course, decompose the phase space of the physical system into layers of constant energy and run trajectories on each layer separately using (68)- The microcanonical average (A),(,!?’) of any physical quantity A could then be obtained on each constant-E’ layer as the time average of this quantity along the trajectory, provided the Hamiltonian dynamics (68) on each layer is ergodic in the microcanonical sense. The canonical average ( A ) , (7‘) would, finally, be calculated from the microcanonical averages (A),,(E’) in accordance with
Both partition functions Z,,(E’) and Z,( 13 in the last expression correspond to the physical system. The problem with this route to canonical averages is that the ratio Z,(E’)/Z, is, in most cases, not known. In addition this route can say nothing about dynamical quantities such as time-correlation functions. In order to get an idea about the possible structure of the manifold M’ traced out by the physical phase space image of a trajectory generating the canonical ensemble for the physical system, we restrict ourselves to the case when the underlying dynamics in the extended space {gi,p’,,s,pS]is ergodic with respect to a true microcanonical distribution function [i.e., w ( s ) = 1 in (48)]. This restriction allows us to invoke our knowledge of the structure of the manifold M which covers part or all [“strong ergodicity”14J5] of an energy surface in the space (~,,p’,,s,p,}and over which the dynamics is actually ergodic in the microcanonical sense. In accordance with the von Neumann-Birkhoff theorem’b13 M is metrically indecomposable, which means that it differs at most by a manifold of zero measure from the manifold M traced out by almost all trajectories started in M . These latter are the ergodic trajectories. The physical phase space image M’of the manifold M is the part of the physical phase space {&’,C,j over which the dynamics is ergodic with respect to the canonical en-
semble. Clearly M’is also metrically indecomposable in the sense that it is different from M’at most by a measure zero manifold of points. The projection of the space (gl,$I,~,ps) onto the space {Ljl’,3,1 is carried out by using (51). This projection is not isomorphic. Each point in the manifold M’is an image of infinitely many different points from the manifold M . The “physical trajectories” tracing out the manifold M’ are continuous since they are images of continuous trajectories in the space {Zjl,jji,s,ps}.But as a result of the lack of isomorphism between the physical and the extended phase spaces, the trajectories in these two spaces are different in an important respect: while the extended phase space trajectory does not cross itself, its physical phase space projection does. In fact every point of the manifold M‘is a bundle of incoming and outgoing trajectoriei. It is this recrossing that accomplfshes the Boltzmann weighting of the different points in the physical phase space. [Recall that the weighting produced by the trajectories can be effectively changed by scaling of the-time; see sections I11 and IV.]
VI. Other Ensembles An important heuristic merit of the general scheme formulated in section IV is in that it allows, in principle, the construction of dynamics for any ensemble of physical systems [for a possible exception see the end of this section]. We restrict ourselves first to the case in which the states of a pbysical system are represented by points of its mechanical phase space. Consider, for example, an ensemble of physical systems in contact with a heat bath but not yet in equilibrium. If the relaxation to equilibrium is slow, one can, for a finite time interval, attribute to this ensemble a nonequilibriumsteady-statedistribution function p&’,$’), where the subscripts t and T refer to the time interval and the temperature of the bath, respectively. In order to construct a dynamics for the ensemble characterized one should properly specify the Hamiltonian (49) by ~,,~(ij’,p’’), and the distribution function (48). In other words, the functions f;,(s), h,,(s), u ( s ) ,and U(S) in (49) and the function w(s) in (48) should be chosen in such a way as to guarantee the equality ZIgi&.Pd wllc
= C(E,Q, qZF‘#,? TT
(71)
where Zw,, is defined by (50) and Zp,,Tis the partition function corresponding to the distribution function p,,,(g’,$’): (72) A direct consequence of (71) is, of course, the equality
where ( )pl;T is the p,..ensemble average and A(G’,p’’)is an arbitrary physical quantity. Then any dynamics generated in the space (~&,s,p,}by the Hamiltonian (49) with properly selected functions J,(s), hih(s),u(s), and U(S) and with the appropriately w(s)-scaled time, which is ergodic with respect to the weighted microcanonical distribution funchon (48), produces time averages of physical quantities equal to their ptirensemble averages. The equation to be satisfied by those sets of functions f,,(s), hlh(S),u ( s ) , u(s), and w(s) which guarantee the validity of the equalities 71 and 73 can be obtained by repeating the steps which led to (58). The result iswhere F(E,Q,T;pi) is an integrable function of p i , and G(s) and so are defined by (57) and ( 5 9 , respectively. The feature that the ps‘dependence can be factored out, as in the rhs of (74) [or (58)], is not a necessary requirement. A less restrictive condition to be satisfied by the functionsJ,(s), hiA(s),u(s), u(s), and W ( S ) ,E
I“
SG(”(E,Q,T;ij’,p’’,p,?)
dp,’ = C(E,Q,T,)p,;T(q“,p”)
(75)
where C is a constant parametrically dependent on E , Q, and T . There are infinitely many sets of functionsl;,(s), hi,($), u(s), u(s),
3172 The Journal of Physical Chemistry, Vol. 92, No. 11, 1988
and w(s) which satisfy (74) or (75), and thus there are infinitely many different dynamics which, in principle, can mimic the P , , ensemble. Finally, we turn to the problem of dynamical simulation of any ensemble of physical systems. In general, the description of states of these systems may require variables additional to their mechanical coordinates and momenta; these, for example, may be thermodynamic variables. We define here three steps that comprise a unified scheme for dynamical simulation of any statistical ensemble. The first step, which can be called the augmentation step, consists of considering the additional variables of interest, e.g., the relevant extensive thermodynamic quantities, as additional “coordinates” and introducing “momenta” conjugate to these ”coordinates”. The new “coordinates” and “momenta’! are then added to the mechanical coordinates and momenta to form the augmented phase space. The latter thus can be denoted as {g,’,jj,’,V’,p+,...}, where V’ is the volume “coordinate”, p + is its conjugate “momentum”, etc. and the prime is used to show that these are physical quantities [as opposed to virtual quantities]. The idea of augmentation was first implemented by Andersen’ for simulating an isoenthalpic-isobaric ensemble. The second step can be called the extension step. In this step a single dimensionless coordinate s and its conjugate momentum p,’representing a “bath’! are added to mechanical or augmented phase space to form the extended phase space {q,’,jj,’,s,p,j or the extended and augmented phase space (4‘,’,jjp’,’,V’,p;, ...,s,ps’). The idea of extension was originally suggested by The third step is to consider an extended and augmented Hamiltonian
“ 2
... + -+ kTv(s) (76) PS
2QU2W ( i = 1, ..., W, X = x , y , z ) and a wpc distribution function Pwrc
=
Wmd,V,P+,...,S,PS) - E ) w(s)
(77)
In (76) W is the “mass” associated with the volume “coordinate” and P is the pressure parameter. Equations 76 and 77 are written in terms of virtual (unprimed) variables related to physical (primed) variables as
wherefix(s), hiA(s),j ( s ) , r ( s ) , ..., u(s) are nonvanishing differentiable functions. In a manner similar to that used in the previous sections we calculate now the wpc partition function
where G(s)
r(s)u(s)... j ( s ) ...w(s)u’(s)
hih(s) n-f;A(s) i.X
Jellinek and
~
so =
PV’-
... - .-],kT) 2Q
Imposing the requirement
C(sd = U E ,Q,W,., .*T;Pv’~...,PS
3~T,P,...(G ’J’,V’,..
(82 )
where L is an integrable function ofp;, ...,p,’and prg ,,,,(ij’,c’,V’,...) is the (T,P,...)-labeled physical phase space distribufion to be generated dynamically, we obtain:
Thus for any physical quantity A
and therefore each dynamics, generated by Hamiltonian (76) and a w(s)-scaled time t’, which satisfies the requirement (82) and produces in phase space (~i,p’i,V,pv,...s,ps} the distribution (77), furnishes time averages of physical quantities equal to their pT,p, ,-ensemble averages:
There are infinitely many different Hamiltonians(76) and defining the scaling of time functions w(s) [(24) with p(s) = w(s)] which satisfy (82), and thus there are infinitely many different dynamics for which (85) may hold. N O S Cformulated ~~ a dynamics for simulating a constant temperature-constant pressure ensemble by augmenting the specific Hamiltonian ( 2 ) . Many other dynamics can be introduced for this ensemble following the prescription given above. A note regarding the applicability of this unified scheme to the grand canonical ensemble should be made here. Implicit in the augmentation-step’is the assumption that the variables added to the mechanical phase space are continuous functions of time. The number of particles N in the system, a variable in grand canonical ensemble, does not satisfy this assumption. One could possibly resolve this problem by time propagating N ( t ) as if it were a continuous function and counting the number of particles by rounding off N ( t ) to its nearest integer. There is, however, a more severe complication preventing a straightforward formulation of dynamics for a grand canonical ensemble. This complication stems from the fact that a varying N ( t ) implies a constantly [and continuously] varying number of degrees of f r e e h m .
VII. Summary In this study we addressed the classical problem of-the relation between mechanics and statistical physics usually referred to as the ergodic problem. More explicitly, we explored-the question of how to generate dynamics for any given statistical mechanical ensemble. The answer to this question is furnished by a unified scheme which is a generalization of the approach suggested recently by N o d for dynamical simulation of a canonical ensemble. This generalization includes clarification of some important points related to implications of the scalings and weightings involved in the procedure. The main results can be summarized as follows. There are three basic elements in the scheme: (a) the extended and appropriately augmented Hamiltonian of the type of (76), (b) the weighted microcanonical distribution function of the type (77) defined in the extended and appropriately augmented phase space of virtual variables, and (c) the distribution function plR,ldefined in the appropriately augmented physical phase space and to be simulated dycamically [(RJ is the set of the control parameters R, char-
3173
J. Phys. Chem. 1988, 92. 3173-3181 acterizing the physical system]. The s-dependent scaling functions in the Hamiltonian, the weighting function in the weighted microcanonical distribution, and the distribution function in the physical phase space should satisfy a well-defined equation of the type (82) [cf. also (80)] in order to assure that the dynamics produced by the Hamiltonian and the appropriately scaled time iscapable, at least in principle, of mimicking the desired ensemble for the physical system. For any fixed pIRIlthis equation has an infinite number of different solutions with respect to the scaling functions and the weighting function. Thus for any distribution in the physical phase space there are infinitely many different dynamics which can, in principle, generate it. The dynamics can be written in terms of any generalized extended and augmented phase space obtained from the corresponding extended and augmented phase space of virtual variables by using transformations of the type (22). Rewriting a dynamics in a different phase space does not have any effect on it. The scaling of the time, however, does. The latter in fact is a way to generate new dynamics even when the Hamiltonian is fixed: although a phase point in its scaled-time evolution traces out the same trajectory as in its unscaled-time evolution, it does it with a different rate. It should clearly be understood that the scaling of the time, (24), on the one hand, and the scalings (78) which relate the physical and virtual variables and which establish an effective coupling between the different variables describing the system of interest and the “bath”, on the other hand, are totally independent. Although a Hamiltonian and a time scaling, Le., a dynamics, can satisfy (82) with a particular distribution plRII,there is no guarantee that this dynamics is ergodic, Le., that it indeed gen. erates in the physical phase space the distribution p I R , ~Because of this, use of the indicative test of ergodicity is recommended. The effect of time scaling on a dynamics ergodic with respect to a specified density distribution in a certain part of the extended and augmented phase space is to produce another dynamics which is also ergodic in the same part of the phase space albeit with
respect to a different, weighted density distribution. We call this invariance feature the scaling invariance of ergodicity. As long as only the averaged quantities are of interest, there are no fundamental reasons to prefer some of the dynamics to the others (we have in mind different dynamics for generating the same ensemble). This, together with the fact that there are no a priori indications for ergodicity, puts the practical considerations of computational simplicity and efficiency as guidelines for selecting particular dynamics in actual numerical studies. It is desirable to write the equations of motion for a specific time scaling in that phase space in which these equations acquire their simplest form. For example, they may partially decouple as in the case of (4). Computational studies of simple systems such as the harmonic o s ~ i l l a t o as r ~well ~ ~ as ~ ~of more complex systems will reveal the ergodic properties of different dynamics and the dependence of these properties on the parameters such as E , Q, and W. These studies may give also an indication of the comparative computational efficiency of different classes of dynamics. The advantage of the M D technique over the M C procedure is that the former can produce not only averaged values of quantities of interest but is able also to follow the detailed time evolution of the system studied and to furnish dynamical information on time-dependent correlations and relaxation processes. Now that we know how to generate many different dynamics for the same ensemble of choice, the task of selecting the particular dynamics which adequately mimics not only the averaged but also the detailed behavior of the system can be addressed. This is particularly important for bringing the simulations closer to realistic situations and for trying to mimic actual experimental conditions. The results of numerical studies using the scheme developed in this paper will be reported in future publications.
Acknowledgment. I thank Steve Berry for many useful discussions. This work was supported by the U S . Department of Energy, Office of Basic Energy Siences, Division of Chemical Sciences, under Contract No. W-31-109-Eng-38.
Photodetachment of Electrons from Dipolar Anions David C. Clary University Chemical Laboratory. Lensfield Road, Cambridge CB2 1E W, UK (Received: July 8, 1987;
In Final Form: Septem’ber 27, 1987)
We present a theoretical study on the spectroscopy of the photodetachment of electrons from dipolar anions. Systems examined include CH2CHO- and GH2CN- for which rotationally resolved experimental spectra have been reported previously. The teehnique involves the calculation of rotationally adiabatic potentials by diagonalizing a Hamiltonian containing a simple potential constructed from an electron-dipole interaction plus a damped polarization term. It is assumed that the electron moves on one of these rotationally adiabatic potentials. Calculations of the energies and widths of the shape and Feshbach Tesonances for; these potentials are identified with the positions and widths of the lines in the photodetachment spectra. A novel analysis of the Hamiltonian for the electron-dipole interaction leads to simple formulas relating to the resonance energies and widths. This analysis also provides a clear explanation of several experimental results including the regular rotational structure observed for the line positions in the photodetachment spectrum, the independenceof the line widths on the symmetric top projection quantum number k for CH,CN-, and the sharp increase of the line width with increasing k for CH2CHO-.
I. Introduction Recent experiments on the photodetachment of electrons from certain dipolar anions have revealed a well-resolved resonance
CH,CN-
L,CH,CN
+ e-
(1) Zimmerman, A. H.; Brauman, J. I. J . Chem. Phys. 1977, 66, 5823. (2) Jackson, R. L.; Zimmerman, A. H.; Brauman, J . I. J . Chem. Phys.
1979, 72, 2088.
0022-3654/88/2092-3173$01.50/0
and the acetaldehyde enolate ani0n~9~ CH2CHO- -k C H 2 C H 0
+ e-
(3) Jackson, R. L.; Hiberty, P. C.; Brauman, J. I. J . Chem. Phys. 1981, 74, 3705. (4) Lykke, K. R.; Mead, R. D.; Lineberger, W. C. Phys. Reu. Lett. 1984,
52, 2221. ( 5 ) Mead, R. D.; Lykke, K . R.; Lineberger, W. C.; Marks, J.; Brauman, J. I. J . Chem. Phys. 1984, 81, 4883.
0 1988 American Chemical Society
