6468
J . Phys. Chem. 1990, 94, 6468-6472
Tests of Certain Approximations in Computations Modeling the Formation of Clusters in Supersonic Flow Lawrence S. Bartell* and Richard A. Machonkin Department of Chemistry, University of Michigan, Ann Arbor, Michigan 48109 (Received: March 5. 1990)
Two approximations made in previous treatments modeling the nucleation and growth of molecular clusters generated in nonequilibrium nozzle flow are examined along with the sensitivity of results to uncertainties in input parameters. The first approximation, namely, the assignment of identical temperatures to all clusters at a given point in the flow, was found to be not seriously in error. The kinetics of mass and heat transfer must be suitably taken into account, however, to obtain the correct difference between the cluster temperature and that of the surrounding carrier gas. The second approximation, according to which the Kelvin equation is applied to critical nuclei but not subsequently to growing clusters, has an appreciable effect. When small clusters are given their enhanced volatilitiesat all times during growth, the slower growth and heat evolution into the medium delay termination of nucleation and result in somewhat cooler clusters that are smaller when fully grown. Although the thermal accommodation coefficient is crucial in the kinetics of heat transfer between clusters and gas, the calculated thermal history of clusters does not depend strongly upon the exact value of this poorly known coefficient. Likewise, final results are relatively insensitive to asssumptions about other uncertain parameters. Overall, results are robust enough to make computer modeling a valuable aid in the interpretation of the types of cluster phases encountered in experimental studies.
Introduction Large molecular clusters produced by homogeneous nucleation from the vapor have provided a novel and effective approach to the study of condensed matter.’q2 They have made it possible to investigate the structure of greatly supercooled 1iq~ids.j.~ They have yielded new crystalline phases, some of which have not been able to be reproduced in the b ~ l k . ~ They - ~ * ~have also offered the opportunity of studying homogeneous nucleation rates for phase changes in liquids and solids on a microsecond time scale.7 In electron diffraction studies of large clusters, nucleation has been induced by the cooling accompanying the adiabatic expansion of a seeded carrier gas flowing supersonically through a miniature Lava1 nozzle. Of particular interest in such studies has been the control the experimentalist can exercise over the condensed phase generated by adjusting conditions of the f l o ~ . ~ *It~ is* *reasonable to suppose that the basis of this control over the cluster phase lies in the range of temperatures experienced by the molecular aggregates at some critical stage of growth. Unfortunately, the miniature size of the nozzle makes it impractical to attempt to follow experimentally the nonequilibrium gas dynamics during the period of nucleation and growth of the clusters. Therefore, a computer modeling of the flow was initiated in an effort to shed light on the processes involved that we would otherwise be ignorant of. As has been described e l s e ~ h e r e ,this ~ . ~line of attack proved to be fruitful even though certain approximations of unknown severity were introduced to render the complex sequence of events amenable to analysis. Examples of the utility of the computer modeling include some unexpected implications. Computations of the developing temperatures of condensing particles yielded persuasive evidence that clusters of certain materials grew in one phase but changed phase spontaneously before they were examined by electron diffractiom2 This prediction was subsequently confirmed by direct ob~ervation.~ The purpose of the present paper is to examine the effects of some of the simplifications incorporated in the previous analyses (1) Bartell, L. S. Chem. Reo. 1986, 86, 49. (2) Bartell, L. S.; Harsanyi, L.; Valente, E. J . J . Phys. Chem. 1989, 93, 6201. (3) Valente, E. J.; Bartell, L. S.J . Am. Chem. SOC.1984, 80, 1451. (4) Bartell, L. S.; Sharkey, L. R.: Shi, X . J . Am. Chem. SOC.1988, l/O, 7006. (5) Harsanyi, L.; Bartell, L. S.; Valente, E. J . J . Phys. Chem. 1988, 92, 4511. (6) Bartell, L. S.; Valente, E. J.; Caillat, J . C. J . Phys. Chem. 1987, 91, 2498. Bartell, L. S.; French, R. J. Reu. Sci. Insrrum. 1989, 60,1223. (7) Bartell, L. S . ; Dibble, T.S . J . Am. Chem. Soc. 1990, / f 2 , 890. ( 8 ) Bartell, L. S . ; Harsanyi, L.; Valente. E. J. NATO AS1 Ser., Ser. B (Clusters) 1987, 158, 37. (9) Bartell, L. S. J . Phys. Chem. 1990. 94, 5102.
0022-3654/90/2094-6468$02.50/0
of cluster formation. These simplifications are discussed in the next section in the context of the nonequilibrium processes governing the formation of clusters.
Aspects of Cluster Formation To Be Investigated A detailed account of the expressions characterizing the nucleation and growth of clusters in nozzle flow has been presented el~ewhere.~-l’Also documented was the implementation of the theoretical relations into a practical program for modeling the course of events taking place in the formation of cluster^.^ This material will not be repeated in detail here. Computer modeling of cluster formation begins by applying a variant of the classical “liquid drop” model of nucleation theoryI2-l4 to a rapidly expanding, cooling element in the flow. It turns out that the equivocal and controversial basis of nucleation theory, itself, does not seriously compromise the inference of the quantities of greatest significance in the present research. Of far greater importance are the relations adopted to describe the decidedly nonequilibrium growth of clusters in the expanding, supersaturated vapor. Of greatest concern for the purposes of interpreting the type of phase (or cluster structure) produced is the thermal history of clusters in the flow. The simplified laws of cluster growth commonly adopted in simulations of cluster formation lead to conclusions in serious disagreement with experimental observations in this laboratory, as related in the preceding papers in this Therefore, the focus in the present paper is upon what one can hope to learn about the growth of clusters by modeling supersonic expansions. Only those relations are introduced in the following text that are needed to define the relevant terms. Clusters are considered to be born when fluctuating aggregates, or germs, forming and evaporating by chance collisions, happen to reach a critical size such that addition of another molecule would lead to unlimited spontaneous growth. In particular, in a supersaturated vapor with partial pressure p corresponding to a degree of supersaturation s =P/Pm
(1)
a critical cluster is one whose vapor pressure matches the ambient partial pressure p of the subject vapor. The augmentation of the critical cluster’s vapor pressure p(r*) above the bulk vapor pressure p m is a consequence of its small radius r* and the resultant Laplace pressure 2 u / r * exerted upon its interior. The augmentation is given by Kelvin’s equation (IO) Hill, P. G.J . Fluid Mech. 1966, 25, 593. ( I I ) Wegener, P. P.; Wu, B. J . C. Ado Colloid InterfaceSci. 1977, 7, 325. (12) Volmer, M.: Weber, A. 2. Phys. Chem. 1925, 119, 277. (13) Becker, R.; Doring, W. Ann. Phys. 1935, 24, 719. (14) Ostwatitsch. K . Z . Angew. Math. Mech. 1942, 23, 1
0 I990 American Chemical Society
The Journal of Physical Chemistry, Vol. 94, No. 16, 1990 6469
Formation of Clusters in Supersonic Flow p ( r * ) / p = = exp(2uum/r*kT) = s
(2)
where u is the surface tension and u, is the volume per molecule in the cluster. This application of bulk thermodynamics to embryonic droplets calculated to contain fewer than a dozen molecules is, of course, a great extrapolation whose justification is only as strong as its success in accounting for observations. Classical nucleation theory, modified by a multiplicative constant r to confer empirical adjustment,ISJ6relates J , the rate of formation of clusters per unit volume per unit time, to temperature and degree of supersaturation by the e x p r e s s i ~ n ' ~ ~ ' ~
J = r K exp(-AG*/kT)
(3)
where the free energy barrier AG* is given by AG* = 4 ~ ( r * ) ~ u / 3
(4)
AG* = '~3?ru3u,2/(kTIn s ) ~
(5)
which becomes
in view of eq 2. The preexponential factor K of eq 3 adopted in the present treatment is that proposed by V ~ l m e r . ' ~ *Defeating l~ any hope of obtaining absolute accuracy in an ab initio treatment of nucleation rate is, among other things, the presence of the cube of the surface tension in the exponent of eq 3. Any small error in u has an extraordinarily large effect upon J. It is widely believed that u depends upon drop size in a manner more or less as given by Tolman's r e l a t i ~ n ' ~ . ~ ~ u(r,T) = udT')/(l
+ 26/4
(6)
Because the parameter 6 has a very poorly known magnitude, the surface tension of very small drops is not well-known. In fact, even the form of the variation of surface tension with drop size is still a matter of dispute.lO*zO-zzWe circumvent the problem to some extent by adopting a plausible magnitude for 6 and then adjusting the factor r in eq 3 to give approximate agreement with experiment under some standard condition. As shown in ref 9, such a compensation is reasonable. Once the critical clusters appear, they begin to grow by accretion of matter in accord with the gas kinetic collision rate, subject to several considerations discussed in detail in ref 9. Heat of condensation is evolved into the medium, affecting the nucleation rate markedly and, thereby, strongly influencing the number of other clusters produced, which in turn evolve heat as they grow. All of this, of course, has a complicated effect upon the temperature of the gas surrounding the clusters. Even more important to a newly born cluster, however, is its own heating caused by its growth. The only available mechanisms for a cluster to get rid of the latent heat of condensation on the time scale of events in the flow are partial evaporation or heat removal by thermal accommodation with the surrounding gas. Explicit coupled differential equations describing the kinetics of heat and material transfer are given in ref 9. Two parameters of concern in the present test of approximations are the thermal accommodation coefficient a of the carrier gas c in collisions with clusters of substance a and the heat capacity C, of the growing clusters. Because there is little information in the literature concerning thermal accommodation coefficients for gases colliding with nonmetals, we have taken Bade's formulaz3 (7)
based on relative molecular masses, as a reference for comparison. (15) Abraham, 0.;Kim, S. S.;Stein, G. D. J . Chem. Phys. 1981,75,402. (16) Reiss, H.In Advances in Chemical Reaction Dynamics; Rentzepis, P. Cohen, E. R. J . Chem. Phys. 1968, 48, 5553. (17) Zeldovich, J. Zh. Eksp. Teor. Fiz. 1942, 12, 525. (18) Tolman. R. C. J . Chem. Phys. 1949, 17, 333. (19) Tolman, R. C. J . Chem. Phys. 1949, 17, 118. (20) Wakeshima, H. J . Phys. Soc. Jpn. 1961, 16, 6. (21) Sinanoglu, 0. J . Chem. Phys. 1981, 75,463. (22) Hooper, M. A.; Nordholm, S. J . Chem. Phys. 1984.81, 2432. (23) Baule, B. Ann. Phys. 1914, 44, 145.
It gives plausible values, but the sensitivity of results to its choice needs to be investigated. Similarly testable is the selection of parameter C, (a dimensionless constant, which added to the vapor heat capacity per mole yields the cluster molar heat capacity in units of R ) . This heat capacity enters directly in the explicit computation of changes of cluster temperature TDthrough the relation equating the time rate of change of energy retained in the cluster to the difference between the rate of energy incident upon it (by condensation) and energy removed (by evaporation and thermal accommodation), or For convenience, calculations are based on the rate of change per unit cluster area. Alternatively, the cluster temperature can be calculated implicitly by the more common approximation of assuming that the left-hand side of eq 8 is so much smaller than either of the two right-hand terms that it is reasonable to adjust the cluster temperature to equalize Ei, and Eout.9*11,24,z5 That is, the dependence of evaporation rate and thermal accommodation on cluster temperature is such that a cluster temperature can be taken to correspond closely to a net energy flow of zero. Results of the two methods are compared in the present work. After clusters are fully grown and comparatively warm, appreciably warmer than the surrounding gas, they cool quite rapidly in the continuing adiabatic expansion. Cooling is caused partly by thermal accommodation with the remaining gas and, partly, by e v a p o r a t i ~ n . ~In ~ ,many ~ ~ cases the cooling is sufficiently large that the clusters find themselves in a supercooled phase. What form they are observed to be in when probed downstream by a laser or electron beam depends, then, upon the kinetics of phase transformation.
Specific Approximations To Be Tested As discussed above, the derived information of greatest interest for the present purposes is the time evolution of the temperature of clusters, particularly during their growth. As shown elsewhere? the present approach, when incorporating certain simplifications, yields fairly robust results, results that do not depend unduly upon the substantial uncertainties of nucleation theory itself. Of uncertain quantitative influence are complications associated both with the adoption of a gas dynamic model of "one-dimensional flow" that neglects the partial mass fractionation known to occur somewhere in the flow6 and with the adoption of a simplified correction for the boundary layer at the nozzle wall. These problems, discussed in ref 9, will not be tested in the present paper. Assumptions To Be Considered. Factors to be examined will include the sensitivity of results to reasonable variations in (1) the poorly known thermal accommodation coefficient (cf. eq 7), (2) the heat capacity constant C, used in eq 8, and ( 3 ) the magnitude of Tolman's parameter 6. In addition, cluster temperatures TD will be calculated by direct updating via eq 8 and also, in some cases, by setting eq 8 to 0 and finding Thlr the temperature required to establish the implied energy balance. Of more fundamental concern, however, are two other approximations introduced in the earlier calculations? First, it was assumed that all clusters produced in a given reference element flowing through the nozzle have a common temperature, irrespective of where they happened to be nucleated in the flow. The present treatment has a provision for keeping separate thermal histories of clusters nucleated at different positions in the nozzle. Second, prior computations handled vapor pressures of clusters discontinously. As discussed when introducing eq 2, critical clusters were considered to have vapor pressures in excess of the bulk vapor pressure owing to their Laplace pressure. This assumption is at the heart of the classical nucleation theory, of course, But as soon as clusters started to grow, their vapor pressures were abruptly switched to those of the bulk at the same temperature. The present treatment (24) Knauer, W. J . Appl. Phys. 1987, 62, 841. (25) Merritt, G. E.;Weatherston, R. C. AIAA J . 1967, 5, 721. (26) Gspann, J. In Physics of Electronic and Atomic Collisions; Datz, S., Ed.; North Holland Press: Amsterdam, 1982; pp 79-96. (27) Klotz, C. E.Nature 1987, 327, 222.
6470 The Journal of Physical Chemistry, Vol. 94, No. 16, I990
Bartell and Machonkin
TABLE I: Comparison of Models 1-3" for Supersonic Expansions of Benzene Vapor in Argon Carrie# at Cbaracteristic Points in the Flow' frac cond = 0.05 frac cond = 0.50 model J" at x X TD T X TD T X ( R ) ? nm 25 I 158 0.16 232 173 0.28 11.2 I 0.76 x 1025 0.145 25 I 158 0.16 232 I72 0.28 2 0.76 x 1025 0.145 11.9 0.20 223 I58 0.34 9.3 0.157 239 144 3 1.05 X IO2' "See text for explanation of models. bInitial temperature, pressure, and benzene mole fraction 298 K, 2.7 bar, and 0.04, with nozzle shape of ref 28 and r = 2 X 6 = 0.1 nm. CAt positions x (cm) along the nozzle from its throat, J,,, is the maximum nucleation rate in m-3 s-', TDis the mean cluster temperature ( K ) and T, the gas temperature at fractions 0.05 and 0.50 of the benzene vapor condensed. dFinal mean cluster radius.
TABLE 11: Sensitivity of Model 3 Results to Input Parameters (Benzene Vapor in frac cond = 0.05 param Jmax at x X TD T 0.160 237 128 ab = 0.300' 1.12 x 1025 239 144 0.448' I .05 x 1025 0.157 0.150 237 153 0.600' 0.93 x 1025 6d = 0.oe I .06 x 1025 0.133 240 I50 0. I / 1.05 x 1025 0.157 239 144 0.28 I .06 x 1025 0.044 252 180 C,(Clust)h
Argon CarrieP) X
0.25 0.20 0.17 0.18 0.20 0.1 1
TD 218 223 223 226 223 240
frac cond = 0.50 T 131 158 172 165 I58 I94
X
0.50 0.34 0.28 0.3 1 0.34 0.21
( R ) , nm
7.7 9.3 11.0
9.5 9.3 12.6
For flow conditions and meaning of symbols,see Table I. bThermal accommodation coefficient of carrier. B a d e value = 0.448. c r = 2 X = 0.448, r = 1.1 X IO3. fa = 0.448, r = 2 6 = 0.1 nm. dTolman parameter for surface tension, eq 6 text, r adjusted to keep J,,, constant. f a = 0.448, r = 8.3 X IO-'. *Effect of variation upon results insignificant.
has an option for applying the Kelvin equation to clusters at all stages of growth, gradually shifting their vapor pressures toward those of the bulk as their radii increase. Models for Comparison. To test the principal assumptions outlined in the previous paragraph, three models of increasing complexity were considered. In the simplest, model 1, clusters are handled in separate groups identified by the position in the nozzle where they were nucleated, but certain properties are averaged. The kinetics of cluster formation is broken into two stages, A and B. In stage A where the nucleation rate rises rapidly and then falls precipitously when heat is given off to the gaseous medium, chster temperatures are taken as Tbal,calculated by equating Ei, and E,,,. This approximation assigns the same temperature to all clusters in a given volume element without regard to where they were nucleated. The resultant temperature is appreciably warmer than that of the surrounding gas. All of the clusters considered are generated in region A. In stage B the clusters continue to grow and continue to have the same temperature at a given position in the nozzle. The temperature, however, is now allowed to adjust itself in accord with governing relations 8 (based on a finite difference between the explicitly calculated quantities Ei,and Eout). Model 1 corresponds to the one adopted in ref 9. Model 2 differs from model 1 by allowing clusters born at different intervals in the nozzle to establish their temperatures T , from the start in accord with eq 8. The different groups of clusters are treated separately at all times, and their temperatures are not constrained to be the same. In both models 1 and 2, clusters suffer the abrupt shift of vapor pressure from that given by the Kelvin equation (for critical clusters) to that of the bulk (as soon as clusters begin to grow). This defect is remedied in model 3, where the Kelvin equation is applied at all times. On the other hand, in model 3 the vulnerability of the more volatile clusters to evaporation is so great at first that numerical instabilities associated with the finite time steps taken may cause temperatures and sizes to oscillate violently, leading to destruction of the clusters. To counteract this numerical artifact, a damping is introduced whereby neither are cluster temperatures allowed to fall below the temperature of the surrounding gas nor are cluster radii allowed to fall below the critical radius. To conserve mass and energy, any excess heat or matter withheld in one step is added to the next. It the case of model 3, the accumulation of new clusters is terminated at the position x, in the flow at which the calculated critical radius begins to increase. It is assumed that clusters with smaller radii than the current critical radius are unstable and evaporate.
X
As a compromise between numerical instability and computation time in models 2 and 3, time steps of 1 ns were taken. Such steps were sufficiently short to allow integration of the governing kinetic relations to be represented by simple summations.
Results Illustrative calculations were carried out on a system of benzene vapor seeded into argon carrier and expanded through a nozzle similar to "nozzle 6", a nozzle used extensively in experiments in this laboratory.28 Expansion conditions paralleled those of typical experiments. The assumed benzene mole fraction was 0.04, and the stagnation temperature and pressure were taken as 298 K and 2.7 bar. With the exception of the exact nozzle Tolman's parameter 6 (which was usually taken to be 0.1 nm), and the constant r, input parameters were the same as those in the analyses reported in ref 9. For comparisons in Table I, constant r was adjusted to a common value, making the maximum nucleation rate approximately IOzs clusters/(m3 s) for all three of the models. A somewhat higher value of r might be more realistic. Calculated temperatures of the gas and growing clusters at characteristic points in the flow are listed together with the final average radii of the clusters. Table I1 tabulates for model 3 the effects of varying the thermal accommodation coefficient, CY,and the Tolman parameter, 6, while adjusting r to preserve the maximum nucleation rate. Values for the cluster heat capacity constant C,, from 4 to 6 (the range of uncertainty) were also tested. Because results were insensitive to this parameter, they are not tabulated. When C,, was varied, cluster temperatures differed, at most, by a fraction of a degree, and cluster sizes by hundredths of a nanometer. Figure 1 portrays the distribution in temperature of clusters along their trajectories as calculated by model 3 and compares the results with those of model 2. A comparison between the distributions in cluster size yielded by models 2 and 3 is shown in Figure 2 . Results for model 2 differ too little from those of model 1 to justify plotting both. Discussion Qualitatively, the results of the computations obtained by all three of the models tested are much like those reported in a number of prior studies of condensation in nozzle f10w.11~15~30-40 At a (28) Valente, E. J.; Bartell, L.S.J . Chem. Phys. 1983, 79, 2683.
(29) Calculations of the nozzle shape in the present paper used the coefficients of ref 28. Calculations in ref 9 replaced coefficient do by 0.01 1 5 cm. Note a misprint in ref 9: the polynomial should read D ( x ) = Ed&". (30) Sherman, P. M.; McBride, D. D.; Chmielewski, T.; Pierce, T. H.; Oktay, E. Aerospace Research Laboratories Report ARL 69-0089 1969.
The Journal of Physical Chemistry, Vol. 94, No. 16, 1990 6471
Formation of Clusters in Supersonic Flow
in many previous studies, it is crucial for the computation 7 ofminor the most important information sought in the current research,
loot
‘0
0.2
0.1
0.3
0.4
x,cm Figure 1. Calculated temperatures as a function of distance along the nozzle from its throat, expansion conditions of Table I. Dashed line, cluster temperatures according to model 2. Solid lines, model 3: top line, maximum calculated cluster temperature; remaining solid lines, successively, correspond to the minimum temperatures of the warmest 2076, 40%, 60% and 80% of the clusters in the flow. The thickness of the lines at smaller distances illustrates the swings in temperatures computed for the embryonic clusters as they slowly grow large enough to overcome their enhanced volatility. The small discontinuity at x = 0.22 cm is at x, where r* began to rise, after which further nucleation was disregarded (see text). Dotted curve, gas temperature for model 3, increasing beyond x = 0.2 because of heat of condensation. I
P(R’
,
8
I
1 I I\\
the thermal history of the growing clusters. As will be seen, the principal difference between the models tested lies in the duration of the time lag. Models 1 and 2 differ in the method of calculation of cluster temperature during the period between onset and termination of nucleation and in the assumption of model 1 that all clusters at a given position in the flow have the same temperature. In the zone of nucleation, model 1 adjusts cluster temperatures to Tbl and then, beyond the zone, allows the clusters to update their temperatures TD in accord with eq 8. Model 2 lets the clusters seek their own temperatures TD a t all times via eq 8. Because the left-hand side of eq 8 is so much smaller than either of the two right-hand terms, the difference between models 1 and 2 is very small as confirmed by Table I. Near the end of cluster growth, a spread of perhaps 5 OC develops in cluster temperatures in model 2, but the mean value remains close to that of model I. Results for model 3 differ significantly. This model is distinguished from model 2 only in its consistent application of the Laplace pressure, which enhances vapor pressure-not only at the moment of nucleation but also at all times afterward. Therefore, while clusters are small, they are appreciably more volatile and, hence, slower to grow and evolve heat. As a consequence, the medium warms more slowly and nucleation persists longer, climbing to a greater maximum rate before peaking. Table I documents this characteristic aspect of model 3, including the resultant more prolific production of clusters. Cluster growth is delayed somewhat, causing the growth to occur in a cooler region of the flow, and the fully grown clusters are smaller because more compete for the same supply of vapor. Figure 2 compares the size distributions implied by the two models. The unrealistically abrupt decline for model 3 inside 5 nm is caused by the neglect of clusters nucleated beyond the position x, in the flow where r* begins to rise. As might be inferred from the preceding paragraph, clusters of model 3 are more sensitive to the condition of their surroundings than are those of the other two models. This is illustrated by the much larger initial spread in temperature of clusters shown in Figure 1. Recently nucleated clusters remain close to the gas temperature until they manage to grow large enough for effects of the Laplace pressure to abate. Later, mean cluster temperatures approach those of model 2. The model 3 results testing variations in the imperfectly established parameters characterizing cluster heat capacity (Can), thermal accommodation coefficient ( a ) ,and surface tension (6) are readily understood, qualitatively. Variations in the parameter C, had too little effect upon the outcome of computations to warrant tabulation, as mentioned in the previous section. This null result is a consequence of the kinetics of heat transfer explained above that makes the cluster temperature TD very nearly the same as Tbal,the temperature at which energy in-flow and evolution would be exactly balanced. The hypothetical temperature Thl is, of course, independent of heat capacity. Increasing a,however, has a significant effect: it makes things happen sooner in a warmer part of the flow. This, in itself, would tend to produce warmer clusters, but of course a higher thermal accommodation coefficient cools clusters more efficiently, countering the effect of the surrounding temperature. The balance struck is shown in Table 11. One clearcut response to an increasing a is in the speed with which nucleation is shut down. The faster is heat carried off by the carrier gas into the nucleation medium, the sooner is nucleation inhibited. Therefore, increasing a decreases the number of clusters and increases their final size. The possible magnitude of Tolman’s parameter 6 is much more uncertain, and the results of its variation are more complex. When 6 is varied from 0 to 0.1 nm, a range that seems most plausible from the fragmentary evidence available, the influence upon the run is modest (Table 11) if parameter r is adjusted to preserve J,,,, the maximum nucleation rate. On the other hand, if 6 is increased to 0.2 nm, there is a large decrease in the surface tension
I 5
IO
15
20
25
R,nm Figure 2. Probability distributions in cluster radius for models 2 and 3 with maxima at 1 1 and 7 nm, respectively.
critical supersaturation the onset of nucleation occurs, and the nucleation rate continues to rise steeply until heat evolved from the condensation terminates nucleation. The rise and fall are so abrupt that clusters are produced only over a limited zone. Since they all have comparable growing times, the spread in their sizes is modest, as shown in Figure 2. What quantitatively distinguishes the present results from almost all of the prior results is the difference between the temperatures of the clusters and the surrounding gas attributable to the collision-controlled pace of heat evolution from the clusters into the surroundings. Although the corresponding lag in heat transfer has been considered to be ~
(31) Chmielewski, T.; Sherman, P. M. AIAA J. 1970, 8, 789. (32) wu, B. J. C. NITS Report AD-775 257, 1974. (33) Wegener, P. P. In Nonequilibrium Flows; Wegener, P. P., Ed.; Marcel Decker: New York, 1969. (34) Wu, B. J. C.; Wegener, P. P.; Stein, G.D. J . Chem. Phys. 1978, 68, 308.
(35) Yang, S.N.; Lu, T. M. Phys. Reo. 1987,835,6944. (36) Yang, S.N.; Lu, T. M. J . Vac. Sci. Technol. 1987, 85, 355. (37) Be, S.H.; Yano, K.; Kawai, S.Jpn. J . Appl. Phys. 1982, 21, 1755. (38) Jaeger, H. J.; Wilson, E. J.; Hill, P. G.; Russell, K. C. J . Chem. Phys.
1969, 51. 5380. (39) Dawson, D. B.; Wilson, E. J.; Hill, P. G.;Russell, K. C. J . Chem. Phys. 1969, 51, 5389. (40) Yamada, 1.; Usui, H.; Takagi, T. J . Phys. Chem. 1987, 91, 2463.
6472
J . Phys. Chem. 1990, 94, 6472-6478
of critical nuclei and an associated catastrophic increase in nucleation rate. Nucleation peaks quickly, considerably earlier in the flow, leading to fewer, larger, and warmer clusters. It can be concluded that one of the approximations adopted in earlier papers in this series, namely, that of assuming all cluster temperatures are the same at a given distance along the nozzle, is not seriously in error. More serious is the abrupt assignment of bulk vapor pressure to newborn clusters. This crude approximation, common to almost all previous model calculations published, has the effect of making the clusters somewhat warmer and larger than they would have been with the same input pa-
rameters in computations with a consistently applied Laplace pressure. Fortunately, however, the effects are modest in magnitude, as are the effects of reasonable variations in the poorly known quantities that enter the calculations. This is particularly true for the derived quantities of greatest concern in the present research program, namely, the temperatures experienced by clusters as they pass through formative stages of their growth. Therefore, despite imperfections in the modeling, the present approach can provide valuable supplementary information in research on the phases and phase changes of molecular clusters generated in supersonic flow.
The Microcanonical Weight Function: Application to Molecular Dynamics Simulations Marek Litniewski Institute of Physical Chemistry, Polish Academy of Sciences, Kasprzaka 44/52, Warszawa, Poland (Received: March 10, 1989; In Final Form: January 22, 1990)
Integration of the classical microcanonical sum of states over the velocities makes it possible to derive the microcanonical weight function in a useful form. Expressions for partial derivatives of thermodynamic quantities (among others for d p / 8 q T ) have been drived by means of this function. The derived formulas are exact for particles as well as for molecules treated clasically. It was shown that differences inherent to the canonical or microcanonical approach as well as the influence of the linear momentum conservation rule could be neglected in the computer simulations. A computer test on 5 12 Lennard-Jones particles has been performed with the use of an algorithm characterized by small fluctuations of total energy. In a few cases evolutions have been carried out for time over t = 2004m/t)'/*.This allowed us to draw interesting conclusions about the measurements performed on a system containing a small number of particles.
1. Introduction
The possibility of comparing computer simulation data and measurements performed on real systems is of great importance. Such comparison is the best method to answer the question of how reasonable is the model of interactions between molecules used to simulate reality. The comparison of quantities easily measurable by molecular dynamics, such as bonding energy and pressure, may not lead to the appropriate answer (e.g., the most widely used potential for water, i.e., the CI potential, leads to drastically erroneous pressure'). A much better way is to compare thermodynamic quantities along with their partial derivatives. It is very useful to apply formulas that make it possible to evaluate the derivatives in a single computer simulation. The relations between the fluctuations and the derivatives of thermodynamic quantities, adequate for the canonical ensemble, can be easily obtained from the known form of the weight function. The classical (NEV; constant energy and volume) molecular dynamics (MD) method2' simulates the microcanonical ensemble, and therefore, these relations cannot be applied. A few MD methods have been presentede7 that make it possible to fulfill the constant temperature (NTV), temperature and pressure (NTp), or pressure and enthalpy (NpH) conditions. The analysis of the partition functions realized by the NTV methods has been ( 1 ) Reimers. J. R.; Watts, R . 0. Chem. Phys. 1982, 64, 95. (2) Evans, M. W., Ed. Dynamical Processes in Condensed Matter; Ad-
vances in Chemical Physics 68; Wiley: New York, 1985. ( 3 ) Allen, M. P.; Tildesley, D. J. Computer Simulations of Liquids; Clarendon Press: Oxford, 1987. (4) Andersen, H . C. J. Chem. Phys. 1980, 72, 2384. ( 5 ) Haile, J. M.; Graben, H. W. J. Chem. Phys. 1980, 73, 2412. (6) Haile, J. M.; Gupta, S. J . Chem. Phys. 1983, 79, 3067. (7) Nose. S Mol. Phys. 1984, 52, 255.
0022-3h54/90/2094-6472$02.50/0
presented by Nose.8 It has been shown that some of these functions were equivalent to the canonical one. In such a case the canonical relations can be applied. The formulas for some of the derivatives (dH/d7')lp and dp/d711s)for the NpH system have been presented by Haile and Grabems The motion of the particles in the NTV, NTp, or NpH methods is disturbed either through a change of the equations of motion or in a direct way (e.g., scaling the velocities). Therefore, these methods differ fundamentally from the NEV one which simulates reality in a direct way. There are two ways to derive the formulas for derivatives valid for the NEV system. One of them is to correct the relations for the canonical ensemble. The general theory of the differences between the two ensembles was discussed by Lebowitz, Percus, and Verlete9 They also derived the formulas for CYanddp/d71Y adequate for an NEV system consisting of molecules having three degrees of freedom. Cheung generalized the formulas to more than three degrees of freedom and presented themlo together with the formula for d p / d v s . Wallace and Straub" derived general relations between the fluctuations for the different ensembles. The formulas for the canonical ensemble can be transformed by use of these relations into a form adequate for the microcanonical ensemble. All the expressions derived in the above three papers (refs 9-1 I ) are exact in the thermodynamical limit (Le., when we neglect the quantities of the relative order of l/N. Another way, presented in this paper, consists of deriving an expression for the microcanonical weight function. Exact formulas for derivatives at constant energy and volume can be easily obtained, by a differentiation, using the weight function in partly ( 8 ) Nose. S . J . Chem. Phvs. 1984. 81. 5 1 1 (9j Lebowitz, J. L.; Percis, J. K.; Verlet, L. Phys. Reu 1967, 153, 250 (10) Cheung, P S . Y . Mol. Phys. 1977, 33, 519. ( 1 1 ) Wallace, D C ; Straub, G K Phys Reo 1983, 27, 2201
0 1990 American Chemical Society
