J. Phys. Chem. 1993,97, 323-331
323
Numerical Simulation of Particle Formation and Growth in Rapidly Expanding Axisymmetric Flows B. J. Jurcikt and J. R. Brock’ Chemical Engineering Department, University of Texas, Austin, Texas 78712 Received: July 17, 1992; In Final Form: October 14, 1992
Results are presented on a systematic investigation of particle formation and growth by homogeneous nucleation and condensation in axisymmetric free jets. The model is validated insofar as possible by comparison with available experimental data. Most of the previous models for particle formation in supersonic jets have been one dimensional. It is shown that one-dimensional models do not completely represent the physics of the flow and therefore cannot display a variety of phenomena found in condensing axisymmetric free jets. It is demonstrated that, in axisymmetric free jets, particle formation occurs initially off-axis due to the cooling in the radial expansion a t the nozzle exit. The region in the jet with the largest average particle diameters can be off-axis or on the centerline depending on the stagnation conditions or, equivalently, the dominant particle growth mechanism. At higher stagnation temperatures, largest diameters occur on the centerline, while for lower stagnation temperatures average particle diameters are highest off-axis. These observations on particle growth characteristics in the jet are from this model study; experimental data for confirmation are not available, An increase in stagnation pressure is equivalent to an increase in nozzle diameter with respect to the temperature field in the expanding jet. In scaled coordinates, the location where nucleation is initiated is invariant with nozzle diameter and the subsequent flow is not self-similar.
Introduction
only with the mass fraction of the vapor condensed and has not examined the effect of the expansion parameters on the particle size distribution. Yang and LuI6 examined the effect of source pressure and expansion angle of a free jet on the size distribution of metallic and semiconductor clusters. However, they used a 1-Danalysis and assumed a constant and small expansion angle, which is not physically reasonable. This approach assumes implicitly (incorrectly) that the cluster formation process and jet expansion are decoupled. The purpose of this paper is to present a general two-dimensional analysis of highly nonequilibrium nucleation and condensation that includes a description of the aerosol size distribution and to present some results using this model. We begin with a brief description of the important physics in the problem then describe the model and numerical method used tosimulate the flow. Thesimulationmethodis thenvalidated by comparison with published experimentaldata. Further results and conclusions on particle formation and growth in expanding flows are then presented. Unfortunatelythere are no experimental data on details of particle evolution in these flows for comparison with model predictions.
In expanding free jets or supersonic nozzles thevapor undergoes an isentropicexpansion,which results in a decrease in temperature, pressure, and density. For most vapors, the saturation vapor pressure decreases with temperature more rapidly than the pressure decreases in the expansion. The vapor, therefore, approaches and crosses the saturation line and becomes supersaturated during such an expansion,provided the isentropecrosses the saturation line and the pressure drop in the expansion is large enough. If the cooling is rapid enough, the gas or vapor becomes quite supersaturatedbefore any vapor condenses and smallclusters or particles are formed by homogeneous nucleation. The formation of aerosol particles in supersonic expansions has been the subject of a great many investigations.1,2These studies have used particle formation to study homogeneous nucleation?+ the structure of small clusters,6and the deposition of ionized cl~sters.~J Particle formation takes place in a variety of applications such as in wind tunncls,9JOrocket nozzle exhausts,II and the sampling of gas from high-pressure cylinders,l2 in all of which particle formation occurs by the same mechanism but is an undesired effect. Furthermore, particle formation in expanding flows may play an important role in the production of high-quality nanometer particles, which have applications in deposition of high-quality films and fabrication of quantum devices, electrooptic devices, and advanced materials including ceramics. Despite the enormous literature on systems in which highly nonequilibrium nucleation and condensation occur there have not been many attempts to model the flow field simultaneously with the particle formation processes. With only a single exceptionl3 the analysis of highly nonequilibrium condensation has been one dimensional.14J5 One-dimensional analysis fails to account for the two-dimensional expansion and shock structure that are inherent in supersonic jets. As a result one-dimensional analysiscannot show the effects of the large gradients that occur off-axis. While the analysis by Davydov” was two dimensional, only results along the axis were presented and no off-axis results were given. In addition, previous analysis has been concerned
meow If an expansion is sufficiently rapid, the flow will pass through the point at which the saturation temperature and pressure are reached without a condensed phase being formed. The vapor is then in a nonequilibrium(supersaturated) state as the expansion progresses because the existing vapor pressure is larger than the equilibriumvaporpressure. The rate of expansion and the kinetics of phase change therefore determine the extent of supersaturation that is reached. In the rapid expansions examined here, particles are formed by homogeneous nucleation and subsequently grow by condensation. The condensation process releases the heat of vaporization to the gas, raising the temperature, and therefore the saturation pressure, thereby decreasing the supersaturation, bringing the expansion path closer to the equilibrium line. Obviously, the ability to simulate the expansion flows in which homogeneous nucleation and condensation occur is limited by the accuracy of the expressions used for the nucleation and condensation processes. Stated succintly, the expressions used here are from the classical theory of homogeneous nucleation
’Current address: American Air Liquide, 5230s. East Ave., Countryside,
IL 60525.
(8
1993 American Chemical Society
324 The Journal of Physical Chemistry, Vol. 97, No. 2, 1993
Jurcik and Brock
With the stated assumptions,the governingequationsdescribing and transfer processes to a particle in the free molecular regime." the gas and vapor phases are written: The well-known inaccuracies of classical nucleation theory l 8 are accepted, primarily because there is no other reliableoption. The expressions used for the classical homogeneous nucleation rate, Jnuo and free molecularmass condensation rate are ~ e l l - k n o w n l ~ ~ ' ~ where bold quantities are column vectors. In eq 2, for axisymand will not be repeated here. metric coordinates the parameter Y = 1 and for Cartesian It is important to note that during an expansion, since vapor coordinates Y = 0. After the mapping defined above is performed, density and temperature, T,decrease with pressure, p, the rate the mapped equations are in the same form as the unmapped of condensation (being proportional to the product pT1I2) equations if the vectors F and G are redefined. The elements of decreases. However, there is not necessarily a commensurate the vectors E,F,G, and S (prior to performing the mapping) are decrease in the nucleation rate because the supersaturation ratio enters exponentiallyin the expression for Jnuc.This suggests that condensation may be more important at higher temperatures and pressures. Also, expansionsin which nucleation is delayed,thereby reaching relatively low temperatures and pressures, may not exhibit strong effects of condensation but may be dominated by pG= E= F= nucleation. We now briefly describe the model and numerical method for simulating particle formation in a supersonic expansion. This PlV* model is not specificto a particular gas/vapor system but is readily adaptable to different chemistries, as well as different nozzle (3) geometries. The restrictions used throughout are the following: In these equations p = the mass density of the carrier gas and (1) continuum mechanics describe the fluid dynamics of the jet; vapor (g/cm3);pl = themassdensityof thecondensing component ( 2 ) molecular transport effects are all negligible compared to in the vapor phase (g/cm3) (if the system consists of a carrier gas convective transport; (3) flow is two dimensional;(4) flow is steady; and vapor then po is the mass density of the carrier gas and p = ( 5 ) the axial Mach number 1 1 everywhere in the flow field; (6) po + pi); u, = the velocity in the axial direction (cm/s); uz = the the particles all follow the gas streamlines with no velocity lag. velocity in the axial direction (cm/s); P = the pressure (dyn/ The approach used is to march in the axial direction (which cm2); 7f = the specific enthalpy of the carrier gas and vapor acts like a time variable), finding the radial profiles of velocity, (erg/g); H = the specific enthalpy + kinetic energy (H= 7f + temperature, density, pressure, and particle size distribution from (012 + uZ2)/2)of the carrier gas and vapor (erg/g); N,= the ith the axis to the jet boundary (or nozzle wall for condensing flows volumetric property of the aerosol distribution (Ni); k = the rate in nozzles) at each axial location by solving the PDE's that describe of mass condensation per volume [g/(cd/s)]; P,= the voluthe flow in the radial direction with an explicit finite difference metric source of particle distribution property Nc [Ni/(cm3/s)]; technique. This approach is alwaysvalid when the Mach number and h = the volumetric heat released due to condensation [erg/ ( u / c ) 1 1 because no feedback from a downstream location can (cm3/s)]. It is assumed that the pressure exerted by the particles occur. The solution domain changes with axial location because due to their random motion is negligiblecompared to the pressure the jet boundary or nozzle (for in-nozzle expansions) changes in exerted by the gas molecules. Obviously, if there is only one the axial direction. The equations are solved in a rectangular component (Le., the carrier gas is also the condensing species) domain by mapping the equations from the r,z axisymmetric then the last equation is not included in the simulation and p = coordinatesystemtoant1,tcoordinatesystem. Themapisdefined PO. by the equations The equations presented above are not the usual conservation equations of mass, momentum, and energy because the source vector S # 0. The fact that S # 0 docs not imply that mass, momentum, and energy are not conserved. The reason for the nonzerosourcevector is that thesystemdescribed by eq2 describes Where rjb is the radius of the jet (or nozzle) boundary at an axial the vapor phase only. Within a differential volume, however, location and ZMd is the axial location of the Mach disk. The there may also exist a condensed phase in the form of particles, intersections of the r-z lines in the physical domain represent the and only the totals of the mass, momentum, and energy of the actual locations of the nodes used in a calculation. The two phases are conserved. The elements of the source vector intersections of the vt lines in the computational domain are the represent respectively the rate of loss (gain) of the conserved locationsof the nodes in the computational domain. The location quantities to the carrier gas and vapor. It is assumed in all these and properties on the jet boundary at an axial location are found equations that the volume of the condensed phase is negligible by using the solution from the previous z (or I ) ,and a method compared to the volume of the carrier gas and vapor. Therefore, of characteristics procedure is employed to calculate the location the calculated density of the carrier gas and vapor is q u a l to the and properties on the jet boundary. This procedure generates density of the carrier gas and vapor if there is no condensed phase the boundary conditions for the node on the jet boundary for the present. subsequent finite difference step. The details of the method of In the formulationpresented aboveno description of the particle characteristics solution at the jet boundary aregiven elsewhere.lg size distribution is specified and a wide variety of simulation Marching in space and solvingfor the radial profiles is referred methods could be used for particle formation and growth. In the to as solving the parabolized NavierStokes (PNS) equations. context of the classical nucleation theory, the simplest model is Dash and Thorpe20 simulated inviscid supersonic rocket exhaust the moment method of McGraw and Saunders24 in which the flows using a PNS method that was extended21g22to incorporate variation of theaerosol sizedistribution is followed exactly through turbulence and mixing. The basic numerical method is patterned the first three moments. The most complete description is the after these models and the fundamental difference is that here solution of the discrete equation in which a conservation equation nonquilibrium effects due to nucleation and condensation are is written for each particle ~ i z e . 2 ~In ~ ~between 6 these simulations accounted for. PNS solutions of nonquilibrium flows have been in rigor are the discrete-sectional method27and finite element obtained by Prabhu et al.23 where chemical reactions that take methods.28 However, for the simulation of the growth of aerosol place in a shock layer are simulated. particles in the two-dimensional expansion a moment method is
I}
The Journal of Physical Chemistry, Vol. 97, No. 2, 1993 325
Particle Formation and Growth the most computationally efficient method and the increase in information that could be obtained from the more complicated descriptionsis probably minimal. The reason for this is that the only important growth mechanisms in these simulations are nucleation and condensation. Thereforeall simulationspresented here will be based on the moment method of McGraw and Saunders,Z4which is an exact, closed-form moment formulation if evaporation and coagulation are negligible. The twocritical assumptions in the moment method of McGraw and Saunders24 are that evaporation and coagulation are negligible. Evaporationis not negligible when particles pass through a shock because of the large temperature jumps that occur in shocks. In our situation we are concerned with the growth and formation of particles and their propertiesupstream of the shock. Furthermore, in many studies a skimmer is placed upstream of the Mach disk location so that the particles do not pass through a shock. Coagulation is negligible in the expansion flows considered here because the time scale of coagulation is much greater than the time scale of flow. This point can be illustrated by the following argument. A characteristic time scale for the flow properties is ~ n , = D/u 1 - 100ps, while thecharacteristic time of coagulation is29T~~~ = 2/KN. With a particle diameter of 0.01 pm and T of 100 K, the coagulation constant K is 1W9 cm3/s, and, for a typical particle concentration, N 1O1O~ m - ~ , the characteristic coagulation time is -0.1 s. The time scales of coagulationare therefore much longer than the flow time scales and so coagulation effects were neglected. For extremely large nozzle diameters, the flow time scales are longer and coagulation could become important; these situations were not examined here, however. In the above equations it is implied that the values of the base variables p, u,, uz, P, H , Ni, and PI can be derived from the axial flux vector E. Since only the flux variables are calculated, two solutionsfor thevector of basevariablesexist, an axially supersonic and an axially subsonic solution. The two solutions correspond to the conditionson each side of a normal shock. Consequently, a d d e procedure was used to calculate the base variables from the computed axial flux values using the supersonic solution. The method of characteristics is used to calculate the location and propertieson the jet boundary. Since highly nonequilibrium nucleation and condensation processes are occurring on the jet boundary, the well-known characteristic equations of gas dynamics’o must be modified to take these effects into account. The basic properties of the jet boundary-namely, that the pressure equals the ambient pressure and the jet boundary is a streamline-remain the same, however, and are used here. The results from this derivationlg are along the lines defined by the equation
-
dr/dz = tan (e and the following equation holds:
- -
+x)
(4)
where 6 = the angle of inclination of the velocity vector with the axis, x = the Mach angle defined by the equation as sin-l(l / M ) = x, and V = the magnitude of the velocity vector (Y= (u? + U ~ ~ ) O The . ~ ) . variable I$in this equation is a term that includes the effects of the noncquilibrium nu,clcation and condensation. A seven-step procedurewas used for finding theconditions on the jet boundary. Details are a~ailable.’~ A similar procedure is used tocalculate the flow within nozzles since in inviscid flow the nozzle contouris a streamline. Therefore, for in-nozzle flows the angle of flow on the boundary is specified by the n o d e geometry and all that needs to be found is the pressure and particle size distribution at an axial location. A six-step procedure for nozzie flow was used.19
TABLE I: Summary of Simulation Conditions Used for Model Validation. conditions P,JP, = 8.03; pcJpm = 9.64 D = 0.0305 cm Po = 15.2 atm To = 295 K M = 1.0 PcJPm= 2000; p c I p - = 2400
ref 49
properties compared radial density profiles
45
correlation function for
~~
M = 1.0 centerline Mach number P c J P , = 91.68; peJpm= 117.22 46 centerline temperature and centerline velocity D = 0.635 cm Po = 10.2 atm To = 250 K M = 1.0 P o J P m 1OOO; PoJpm = 1000 50 centerline Mach number and two-dimensional structure M = 4.4 D = nozzle diameter; M = Mach number at exit; T = temperature; P = pressure; subscript o refers to upstream stagnation; e refers to exit condition, and
refers to the receiving conditions.
The calculation of the expanding flows starts either from the nozzle exit or the nozzle throat (for particle formation in nozzles). When simulating an underexpanded sonic jet the pressure at the nozzle exit (PJ > ambient pressure (P-), and, in order for the jet boundary to expand to the ambient pressure, an initial angle needs to be specified. This initial angle is found from a PrandtlMeyer expansion,3’J2 assuming a constant heat capacity ratio (y = Cp/Cu)and no particle formation within the first axial step. Due to the largeunderexpansionssimulated here (P,/P, = lOOO), a grid expansion technique is used in the near-field region. When extremely large underexpansions occur the Prandtl-Meyer angle is greater than 90’. For these simulations a large initial angle (typically 45-70’) is specified. At large underexpansions the flow during the Prandtl-Meyer expansion will no longer be governed by continuum m e ~ h a n i c sand ~ ~ so , ~the ~ location of the jet boundary is somewhat vague and the above assumption does not affect the core of the jet. Validation of the Model and the Numerid Mctbod. In numerical simulations it is important tovalidate the computational procedures and to ensure that the model used represents the physical system well. In this section the procedures used to validate the model and numerical method are described. We begin by verifying the model for noncondensing supersonic jets by comparing with published experimental data, empirical formulas, and other numerical simulations. The model is further validated by comparing with the experimental data of Dankert and Koppenwallner’s for the condensation of nitrogen in fret jets. The simulation of noncondensingjets is a good way of validating the numerical procedures because the physical system is less complicated and the complexities of the numerical procedures are reduced. Furthermore, since supersonicjets have been used extensively for many years there is a large literature for comparison. Table I lists the simulations performed for noncoildensing jets, the reference from which the conditions were taken, and which properties were compared. The results comparing the work of Ashkenas and Sherman36 and Moosmuller et aL3’ are centerline properties, therefore it is easy to see the quality of the simulationsfrom a direct comparison. These comparisons are presented in Figures 1-3; the calculated values are obtained from the solution of q 2 with no adjustable parameters. As can be seen from these figures the model represents the centerline data well for a wide range of operating conditions. Accurate centerline solutions imply that the rate of expansion is calculated correctlyand that the isentropic condition is maintained correctly. The rate of expansion on the centerline is determined by the jet boundary and barrel shock locations, so the total of these effects is simulated correctly.
Jurcik and Brock
326 The Journal of Physical Chemistry, Vol. 97, No. 2, 1993
OM)
0
2
1
6
8
1
0
1
2
00
02
06
04
08
I D
I2
r/D
Dialance from Nozzle lmml
Figure 1. Comparison of computed (line, no adjustable parameters) and experimentaP7 (data points) centerline temperatures.
Figure 4. Comparison of computed (line, no adjustable parameters) and experimental4 (data points) radial density profiles at z / D = 2.0.
ann -
7sn-
~-~
r - .
~
0
2
4
6
8
1
0
1
2
I
2
4
6
8
1
0
2
r 4
.
I
6
.
1
8
.
I
10
dD
Distance rrom nozzle Imml
Figure 2. Comparison of computed (line, no adjustable parameters) and e ~ p e r i m c n t a l(data ~ ~ points) centerline velocities.
0
0
1
2
zlD
F g r e 3. Comparison of computed (dotted line, no adjustable parameters) centerline Mach number with empirical ~ o r r e l a t i o n ~(solid ~ J ~ line).
The comparison with the empirical formula for the centerline Mach number, Figure 3, is also in good agreement with the calculations. The correlation function has been well established for its representation of experimental data for highly underexpanded jets.38J9 The correlation function was found by using a method of characteristics solution for an ideal gas expanding into a perfect vacuum and then the curve fit to be valid for z / D 1 5 on the centerline. The initial Mach number in the method of characteristics solution was taken as 1.01, which may explain the slight but consistent difference between the calculation and the empirical formula, since these simulations used a Mach number of 1.0. The other possible explanation for the slight differences is that the method of characteristics more accurately maintains the isentropic condition, while MacCormack’s method is known to have an inherent numerical viscosity, Le., numerical diffusion. The increase in entropy due to the numerical viscosity would decrease the calculated Mach number. As can be seen from Figure 3, however, the difference is slight and the calculation agrees well with the empirical formula and shows that the model can handle highly underexpanded jets and that the boundary procedure used for large Prandtl-Meyer angles is valid. The results presented above are on the centerline only and do not necessarily verify that the model is representing the physical system off-axis as well. In order to verify that the simulations accurately represent the behaviom off-axis some two-dimensional comparisons of noncondensing jets will be presented. Figure 4 shows the calculated (no adjustable parameters) and experimentally measured40radial density profile at a constant axial
Figure 5. Comparison of this calculation (line, no adjustable parameters) and Teshima and Sommerfeld’s calculations4~(data points) for the centerline Mach number.
distance. The density profile in this figure is normalized to the upstream density, po, and it can be seen that the calculated density profile agrees with the calculated profile extremely well. The isentropic core, barrel shock, and jet boundary are obvious in this figure and demonstrate the large gradients that occur within the underexpanded jets. Notice that the density of the gas closer to the centerline than the barrel shock is lower than the ambient density but increases rapidly at the barrel shock. A further test of the two-dimensional accuracy of the model was performed by simulating a underexpanded Mach 4.4 jet exiting from a nozzle. The Mach 4.4 jet simulated here was also simulated by Teshima and Sommerfeld41 using a time split piecewise linear method (PLM) and experimentally investigated by using laser-induced fluorescence. A comparison of the centerline Mach number calculated by the PLM and the PNS methods is presented in Figure 5 . The agreement between the two computational methods further proves the validity of the PNS solution. The drastic decrease in Mach numbers shown in Figure 5 is due to the intersection of the barrel shock with the jet axis. A Mach disk is not present in this jet so the flow just downstream of the Mach number decrease remains supersonic. The location of the intersection of the barrel shock with the jet axis was determinedexperimentally by Teshima and Sommerfeld41 and this location agrees well with the two calculations. The results presented above suggest that the simulations of noncondensing jets accurately represent the physical system. The accuracy of simulations in which condensation occurs has been checked by a comparison with published experimental data. Dankert and K~ppenwallner~~ present a table of experimental upstream stagnation conditions and the corresponding (condensation) onset points for an underexpandedsonic nitrogenjet. Eight of these conditions were simulated to test the accuracy of the condensingjet model. We do not show all these simulations but present two selections, Figures 6 and 7. The calculated expansion paths along the centerline are shown in Figures 6 and 7. All conditionssimulatedarelistedinTableI1.Whiletheexperimental upstream stagnation conditions are well defined in these simulations, it is less clear what calculation results should be explicitly compared with the experimental onset points. The reason for this is that an experimental onset point is not a well-defined
The Journal of Physical Chemistry, Vol. 97, No. 2, 1993 327
Particle Formation and Growth
0
IO
20
30
40
10
$0
70
80
PO
100
o
Temperature [KI
Figure 6. Calculated expansion path (no adjustable parameters) for condition 7, Table 11, giving comparison with experimental onset points.35 in 3
I
t
K
101
o
IO
20
30
4n
50
60
ao
70
90
IM)
Temperature [K]
Figure 7. Calculated expansion path (no adjustable parameters) for condition 2, Table 11, giving comparison with experimental onset points.35
TABLE II: Conditions Used To Verify the Nonequilibrium i ~ the ~ Nucleation and Condensation C d ~ u l r tand Experimentally Determined Onset Points42
1 2 3 4 5 6 7 8
0.50 0.50 0.50 0.50 0.50 0.50 1.o 1.o
3 3 3 3 3 3 0.75 0.75
120 125 160 175 190 200 125 150
34.5 29.9 4.85 2.20 1.30 1.07 2.54 0.76
36.2 36.2 27.6 24.1 22.5 22.4 26.6 22.7
quantity. From an experimental point of view an onset point can be viewed as the point where condensation processes begin to measurably affect the flow. However, for condensation to affect the flow there must be a substantial density of particles present, which implies that the nucleation processes must have occurred just upstream of the onset point. The result of this obscure experimental definition is that no explicit numbers can be compared and only qualitative comparisons are possible for onset points. The calculated onset points in Figures 6 and 7 can be found by the break in the expansion path toward the vapor pressure curve. Once the break has begun, if the condensation rates are large enough, the approach of the vapor pressure curve is quite sharp. This break leads to the nose profile shown in Figures 6 and 7. For the other expansion paths, the onset point cannot be easily determined, and examining other properties such as nucleation rate or condensation rate is more useful in determining the calculated onset points. On the basis of all eight simulations (conditions in Table 11), we infer that the simulations do a better job of representing the experimental data for the conditions with lower upstream stagnation temperatures, To.There are several possible explanations for this. At the lower To'sthe effect of condensation is stronger so the onset point is more sharply defined and therefore easier to measure with accuracy. Also, at the higher To's nucleation and condensation is delayed until lower temperatures and pressures are reached, therefore free molecular effects could be important in the clustering or nucleation process. The results
IO
20
in
40
so
60
70
80
90
100
Temperature [ K l
Figure 8. Calculated expansion paths for nine To's (120, 125, 130, 140, 150, 160, 175, 190, 200 K) and other conditions listed in Table 11.
presented in Figures 1-7 demonstrate that the model qualitatively and in some respectsquantitatively represents the physical system. The disagreement between some of the calculated results and the experimentally determined onset points could be due to the inaccuracies of the nucleation theory or to the free molecular effects that will be present at low temperatures and densities. Free molecular effects become important at low temperatures because at large distances from the nozzle (and therefore low temperatures) the flow becomes frozen.42 Neither free molecular effects nor a reliable nucleation theory (since one does not exist) could be incorporated into the framework of the model so it was decided that for simulations in which low onset temperatures occurred only qualitative results could be achieved.
Results The comparisonswith experiment presented above suggestthat the model represents the physical system reasonably well. In this section the effect of systematic variation of operating parameters on particle formationand growthwill be explored. The parameters that are varied for condensing underexpanded sonic nitrogen jets are upstream stagnation temperature, To, nozzle diameter, D, and upstream stagnation pressure, Po. We conclude with the simulation of a carrier gas (N2) and condensing vapor (H20) expanding through a supersonicnozzle, where experimentalresults are available. From our comparisons with experimental onset points, it is clear that as To increases the onset point is delayed to lower temperatures and the expansion path does not approach the equilibrium line. To examine this effect further three more simulations were run at the same conditions but with To's of 130, 140, and 150 K. The expansion paths for these nine temperatures for conditions of Table I1 are shown in Figure 8. Figure 8 shows that after the onset point the gas approaches the equilibrium line. The approach to the equilibrium line is sharper for the lower To's. This result is not surprising since at the lower To's the onset point is at a higher pressure; therefore the rate of condensation is larger and the approach to equilibrium is more rapid. Another interestingeffect seen in Figure 8 is that, after following theequilibrium line in the pressure range 5-1 5 Torr, the expansion paths then begin to deviate from the equilibrium line again. A second deviation from the equilibrium line has been observed experimentally previously by Duker and Koppenwallner43 and Koppenwallner and Dankert.I4 This leads to further large supersaturations and can lead to a second region of substantial nucleation rates. This isdemonstrated in Figures 9 and 10, which show the nucleation rates and supersaturation ratios for nitrogen expansion for the cases of Toequal to 130,140, and 150 K.Two regions of large supersaturation ratios and high nucleation rates can be seen clearly in these figures. The reason the expansion path deviatesfrom the equilibriumline as the expansionprogresses is that the condensation rates decrease with distance along the centerline. The kinetics of condensationare not sufficientlyrapid to maintain the equilibrium values at the larger values of z / D (equivalently, lower pressures, temperatures, and gas densities).
Jurcik and Brock
328 The Journal of Physical Chemistry, Vol. 97, No. 2, 1993
100.0
1017 1016
I
3
Y
?5
5
1015 1014
76.0
T
1013
60.0
10'2
loll
.o
-a
liS.0
ln'n 109
:1 0 8
2
10'
10 6 10 5 0
1
2
3
4
5
6
7
zlD
Figure 9. Calculated nucleation rates along the jet centerline. in
-
SI~gnmonT
i
l4OK
S l r g n ~ i m1
i
IUIK
Figure 12. Temperature distribution in a noncondensing nitrogen jet. Stagnation conditions are summarized as condition 1 in Table 11. Source terms are zero everywhere. ioo.0
0
1
2
3
4
5
6
7
ZlD
\
Figure 10. Calculated supersaturation ratios along the jet centerline.
-
26.0
J
inxin
P O
r/D
Upstream stagnalion Iempcraturr IKl
Figure 11. Final average particle diameters and concentrations as a function of To at jet centerline.
TABLE III: Stagnation COII~MOM Simulated To Examine the Effect of D rad P . ~~~
~
~~
condition D,cm Po,atm To,K condition D,cm Po,atm To,K 1 2 3 4
5 6 7 8
0.50 0.50 0.50 0.50 0.50 0.50 0.50 0.50
1.5 2.0 2.5 3.5 4.0 4.5 5.0 3.0
120 120 120 120 120 120 120 120
9 10 11 12 13 14 15
1.0 0.75 0.25 0.10 0.075 0.05 0.01
3.0 3.0 3.0 3.0 3.0 3.0 3.0
120 120 120 120 120 120 120
The effect of Toon the final concentration of particles and the average particle diameter is shown in Figure 1 1. From this figure it is apparent that the concentration of particles increases with To, while the average particle diameter decreases. The decrease in average diameter and increase in concentration is consistent with thescaling relations for condensingjets.Figure 11 shows that the average particle diameter decreases approximately linearly with To, while the concentration of particles increases with To approximately exponentially. From Figure 11, it is apparent that in an expansion large concentrations of ultrafine particles can be produced. From the P-Tdiagrams it is apparent that significant changes occur in the expansion process by varying Toonly. It is worthwhile at this point to examine this effect in a little more detail. The three-dimensional profiles of temperature, average diameter, and particle concentration will be discussed for a high Tocondition (condition 6 in Table 111) and a low Tocondition (condition 1 in Table 111).
Figure 13. Temperature distribution in a condensing nitrogen jet. Stagnation conditions aresummarized as condition 1 inTable 11. Source terms are nonzero everywhere.
Figure 12 shows the temperature surface for a noncondensing nitrogen jet with To equal to 120 K (with condition 1 of Table 111, calculated by setting all source terms equal to zero). Figure 13 shows the temperature surface also for condition 1, but with condensation occurring. A comparison of Figures 12 and 13 indicates that the strong condensation that occurs in the low-To jet occurs initially off-axis. The condensation region in this situation appears to be nearly annular, beginning almost immediately downstream of the nozzle due to the isentropic cooting in the Prandtl-Meyer expansion. The temperature rise due to condensation is larger nearer the axis, so that on the axis the effects of condensation are delayed longest and the effects are greatest. This effect can only be found from a two-dimensional analysis and has not been predicted previously. Furthermore, condensationc a u m the jet boundary to expand less than in the caseof the noncondensingjet. The temperaturesin thecondensing jet are greater than in the noncondensingjet and the barrel shock is stronger in the condensingjet. We have observed that at higher To the heat effects due to condensation are not noticeable and that qualitatively the expansion appears identical with that of a noncondensing free jet. The variation of average particle diameters for the low- and high-Toprofies exhibited somequalitativedifferences. The threedimensional profiles of average particle diameter are shown in Figure 14 for the low-Tooaseand Figure 15 for the high-T, case. For the low-To jet the average diameter increased with radial distance from the axis, while the higher Tojet has the largest average diameter on the centerline. The reason for this qualitatively different behavior is that as To is increased, nucleation
The Journal of Physical Chemistry, Vol. 97, No. 2, 1993 329
Particle Formation and Growth
0
2
I
4
z/D 0
Firmre 16. Particle concentrations on iet centerline and iet boundary for coidition 1 of Table 11.
Figure 14. Average particle diameters in a condensing jet for condition 1 of Table 11. The vertical axis is linear and ranges from 0.0 to cm.
dD
Figure 17. Average particle diameters on jet centerline and jet boundary for condition 1 of Table 11.
rP Figure IS. Average particle diameters in a condensing jet for condition 6 of Table 11. The vertical axis is linear and rangts from 0.0 to 10" cm.
-E 5 E
100-
.
m-
'
60-
and condensation are delayed until lower temperatures and pressures. The condensation rate for higher Tojets is therefore lower and is not sufficiently rapid to adjust to the changes in the fluid properties. For the higher Tojet the largest condensation rate is found on the centerline where the largest values of temperature, density, and pressure are found (within the barrel stock). The largest particles are therefore found on the axis for the higher Tojet. For the lower Tojet the condensation rate is sufficientlyrapid throughoutthe flow field to adjust to the property changes during the expansion. The main factor in determining the condensation rate in the low-To jet is the supersaturation ratio. The location where the lowest temperatures occur (and therefore highest supersaturations) therefore exhibits the largest condensation rates. As in the high- Tojet, the largest temperatures are found on the axis (again, within the barrel shock), so for the low-Tojet the highest condensation rates are found off-axis. The largest particles in the lower Tojet are therefore found off-axis. At high To the particle growth is 'kinetically" limited because the driving force of condensation (i.e., the supersaturation ratio) has no effect on the growth rate, while for the low-To case the flow is "mass transfer" limited because the supersaturation determines the region of fastest growth. Figure 16 shows particle concentrations at the jet centerline and at the boundary for condition 1 of Table 11. Condition 6 of Table I1 is used in the calculation of results shown in Figure 17 for the average particle diameter. These both imply that the flow on the jet boundary is essentially frozen because the condensation rates are not sufficiently large to grow the particles beyond their critical cluster size. Along the jet boundary the density and pressure are so low that the collision rates are extremely low; the nucleation rates and condensation rates are therefore almost negligible. Along the jet boundary, free molecular effects arc probably dominant andquantitative behavior should not be believed in this region. It is probably safe to say, however, that along the jet boundary the effects of condensation
E 40
-.
207 0
.
.. .
I
I
..
. .
I
. . . .
2
I
.
.
.
3
. 1
z/D
Figure 18. Comparison ofjet centerline temperature profiles in condensing and noncondensing jets for identical stagnation conditions corresponding to condition 1 of Table 11.
are kinetically limited becauseof the extremelylow collision rates and that the flow closely approximates frozen flow. It should be noted that P, = 0.1 Torr in all the simulations in Table I1 and that varying this parameter in' the range 0.01 Torr I P, I 0.5 Torr did not affect the results. It is surprising that the effect of condensation in an underexpanded jet is to decrease the expansion of the jet boundary. In fact, this variation of the jet boundary does not affect the solution on the centerline. This effect can be seen in Figure 18 where the centerline temperature profiles of a condensingand noncondensing jet of condition 1, Table 11, are presented. Figure 18 clearly shows that the difference in jet boundary location for the two situations does not affect the conditions on the centerline. The drastic divergence of the two solutions at z / D 1.2 is caused by the heat released by condensation and is not due to the differences in the jet boundary location. The effects of the upstream stagnation pressure,Po,and nozzle diameter, D, on condensing jets have been examined. Table 111 shows the conditions under which simulations were performed for this investigation. The expansion paths for the different Po and D simulations were studied. The expansion paths show that the effect of increasing Po is similar to the effect of an increase in D, in that the departure from the equilibriumline decreases with an increase in D or Po. This effect can be explained as follows: as the stagnation pressure increases, the condensation rate increases
-
330 The Journal of Physical Chemistry, Vol. 97, No. 2, 1993
Jurcik and Brock
TABLE Iv: Conditions for Water-Nitrogen Simulations condition To,K P,,atm 1 2 3 4
351 340 350 310
3.74 3.14 3.14 3.14
w,,
0.02 0.02 0.02 0.02
condition To,K Po,atm 5 6
7 8
370 350 340 351
3.74 3.74 3.74 3.74
w,
0.015 0,015 0.015 0.015
r/D
Figure 19. Effect of increase of nozzle diameter, D, and stagnation pressure, Po, on jet centerline particle concentrations. Solid line corresponds to condition 4 in Table 111, dotted line to condition 10 in Table 111, an,d dashed line to condition 8 in Table 111. 0
Temperature [K]
Figure 21. Comparison of calculated and experimentally" determined onset points for condensation of water in nitrogen carrier gas.
0
I
2
3
1
r/D
Figure 20. Effect of increase of nozzle diameter, D, and stagnation pressure, Po, on jet centerline particle diameters. Solid line corresponds to condition 4 in Table 111, dotted line to condition 10 in Table 111, and dashed liny to condition 8 in Table 111.
approximately linearly because the temperature profiles are approximately identical up to the onset point. The effect of an increase in D is a linear increase in the flow time scale, so the mass condensation rate increases approximately linearly at any set of conditions. An increase in Po is therefore similar in effect to an increase in D, because the break in the expansion path occurs at a higher temperature and pressure. The effect of these changes on particle properties is indicated by Figures 19 and 20, giving a comparison of the particle concentration and average particle diameter for an increase in Po and d over the jet of condition 8 in Table 111. From the temperature profile, the effect of an increase in either Po or D is to move the onset point closer to the nozzle exit. The temperature profiles are almost identical for the Po = 3.5 atm and the D = 0.75 cm cases. The particle properties are not the same in these situations, however, which implies that identical onset temperatures do not necessarily imply that the particle propertics are identical. As Figures 19 and 20 show, an increase in nozzle diameter increases the average particle diameter more than does an increase in stagnation pressure, while the particle concentrationdecreases less for an increase in stagnation pressure than for an increase in nozzle diameter (for the situation when the onset temperatures are the same). These results imply that, while it is possible to have a self-similar onset temperature by changingltheupstream stagnation pressure and nozzle diameter, the particle properties will not be identical. Carrier Cas and Vapor in a Supersonic Nozzle. All the simulations resented previously have been concerned with the nucleationa d condensationof nitrogen and have not considered the effect of a carrier gas. In addition the simulations have all been for free jets and have not considered particle formation in the flow of a gas in a nozzle. The model describedpreviously can describe a carrier gas and flow in a nozzle; so for the sake of completenis a few results from two component nozzle flow will be given. The physical system considered is the expansion of water and njtrogen through a supersonic nozzle. Table IV describes the stagnation conditions simulated. The parameter w, in Table IV is the initial mass fraction of water. The initial
np
conditions are found by assuming no condensate at the sonic throat, and the sonic conditions are found by assuming isentropic flow with a constant heat capacity ratio of 1.4. The nucleation of water in a carrier gas has been explored by Jaeger et aL4' and a comparison of the calculatedand experimental onset points is shown in Figure 21. The nozzle used in the simulations is defined by the effective area ratio as a function of distance given by Jaeger et al.4' The stagnation conditions for each onset point were not listed by Jaeger et ~ 1 . : ~so there is no way of knowing if the simulations performed here are a direct comparison of the experimental conditions. However, tne range of stagnation variables is correct, so a qualitative comparison is possible. The calculated onset points are easy to determine from the temperature profiles because the large heat of vaporization of water changes the temperature significantly. As shown in Figure 21, the simulations compare well with the experimental onset points and show that the model can describeflow in a nozzle as well as including a carrier gas.
concllrsions The results presented are the first systematic investigations of two-dimensional simulations of nucleation and condensation in expanding flows. From the results the following were obtained: (1) Particle formation occurs initially off-axis in all jets due to the cooling in the radial expansion at the nozzle exit. Similar behavior has been found in supersonicnozzle nitrogen flow where the pressure at the nozzle exit is greater than the chamber pressure.48 (2) The region with the largest particle diameter can be offaxis or on the centerline, dependingon the stagnation conditions. At higher stagnation temperatures the averagediameter is largest on the centerline, while for lower stagnation temperatures the average diameter is highest off-axis. (3) An increase in Po is similar to an increase in D for the temperature profile but not for the particle size distribution. (4) In dimensionless coordinatesthe location where nucleation begins is the same irrespective of nozzle diameter. Subsequent to particle nucleation the flow differs because of the variation in flow time scales. Acknowledgment. This material is based upon work supported by the National Science Foundation under Grant ECS-91 19043.
References and Notes (1) Wegener, P.P.;Wu, B.J. C. Adu. ColloidInterfaceSci. 19ll.7,325. (2) Wegener, P.P.In Nonrquilibrium Flow with Condensation;Wegener, P. P.,Ed.; Marcel Dekker: New York, 1969; Chapter 4, p 163. ( 3 ) Stein, G. D.; Wegener, P. P.J. Chcm. Phys. 1967, 46, 3685.
Particle Formation and Growth (4) Clumpner, J. A. J . Chem. Phys. 1971. 55 (lo), 5042. ( 5 ) Abraham, 0.; Kim, S.S.;Stein, G. D. J . Chem. Phys. 1981,75( I ) , 402. (6) Mark, T. D.; Castleman, A. W., Jr. Adu. At. Mol. Phys. 1985,20, 65. (7) Takagi, T. Z . Phys. D At. Mol. Clusters 1986,3,271. (8) Fukushima, K.; Yamada, I. J. Appl. Phys. 1989.65 (2), 619. (9) Willmarth. W. W.; Nagamatsu, H. T.J . Appl. Phys. 1975,23,1089. (LO) Hill, P. G. J . Fluid Mech. 1966, 25 (3), 593. (11) Crowe, C. T.; Willoughby, W. AIAA J . 1%7,5, 1300. (12) Wen, H. Y.; Kasper, G.; Montgomery. D. J . Aerosol Sci. 1987,19 (1).
(13) Davydov, L. M. Fluid Dyn.-Sou. Res. 1971,I (I), 90. Translated from: Mekhanika Zhidkosti i Gaza 1971,3,66. (14) Koppenwallner, G.; Dankert, C. J. Phys. Chem. 1987,91, 2482. (15) Clarke, J. H.; Delale, C. F. Phys. Fluids 1986,29 (3,1398. (16) Yang, S.N.;Lu, E. M. J. Vac. Sci. Technol. 1987,B5, 355. (17) Hidy, G.M.;Brock, J. R. The Dynamics ofderocolloidal Systems; Pergamon: Oxford, 1970. (18) Zettlcmoyer, A. C., Ed. Nucleation; Marcel Dekker: New York, 1969. ~. (19) Jurcik, B. J. Ph.D. Thesis, The University of Texas at Austin, 1990. (20) Dash, S.M.; Thorp, R. G. AIAA J . 1981, 19 (7),842. (21) Dash, S.M.; Wolf, D. E. AIAA J . 1984, 22 (7),905. (22) Dash. S.M.; Wolf, D. E.;Seiner, J. M. AIAA J. 1985,23 (4).505. (23) Prabhu, D. K.; Tannehill, J. C.; Marvin, J. G. AIAA J. 1988,26(7), 808. (24) McGraw, R.;Saunders, J. H. Aerosol Sei. Technol. 1984, 3, 367. (25) Cqrls, J. C. Ph.D. Thesis, The University of Texas at Austin, 1988. (26) Brock, J. R. In Theoryof Dispersed Multiphase Flow; Meyer, R. E., Ed.; Academic Press: New York, 1983;p 135. (27) Gelbard. F.; Tambour, Y.; Seinfeld, J. H. J. Colloid Interface Sci. 1980,76 (2),541. (28) Tsang, T. H.; Brock, J. R. Aerosol Sci. Technol. 1983,2, 31 1. (29) Friedlander, S.K. Smoke, Dust and Haze; John Wiley and Sons: New York, 1977.
The Journal of Physical Chemistry, Vol. 97, No. 2, 1993 331 (30) Ferri, A. In High Speed Aerodynamics and Jet Propulsion; Scars, W. R., Ed.; Oxford University Prcss: Oxford, 1955;Vol. VI. General Theory of High Speed Aerodynamics; Section G.19. (31) Shapiro, A. S.The Dynamics and Thermodynamics of Compressible Fluid Flow; Ronald: New York, 1957. (32) Landau, L. D.; Lifshitz, E.M.Fluid Mechanics, 2nd ed.; Pergamon: Oxford, 1987. (33) Bird, G. A. Rarefied Gas Dynamics Symposium; Academic Press: New York, 1980; Vol. 2, p 681. (34) Campbell, D. H. J . Spacecraft 1989,26 (4),285. (35) Dankert, C.; Koppenwallner,G. RarefiedGas DynamicsSymposium; Academic Press: New York, 1978;Vol. 2, p 1107. (36) Ashkensas, H.; Sherman, F. S.Rarefied Gas Dynamics Symposium; Academic Press: New York, 1966; Vol. 2, p 84. (37) Moosmuller, H.; Herring, G. C.; She, C. Y. Opt. Lett. 1984,9(12), 536. (38) Abraham, 0.; Binn, J. H.; DeBocr, B. G.; Stein, G. D. Phys. Fluids 1981, 24 (6),1017. (39) Matsuda, T.; Umeda, Y.; Ishii, R.; Yasuda, A,; Sawada, K. Mem. Fac. Eng., Kyoto Univ. 1987,49 (l), 84. (40) Kutatdladze. S.S.;Novopashin, S.A,; Perepelkin, A. L.; Yarygin, V. N. Sou. Phys. Dokl. 1988,32 (7),548. (41) Teshima, K.; Sommerfeld, M. Exp. Fluids 1987,5, 197. (42) Cattolica, R.; Robben, F.; Talbot, L.; Willis, D. R. Phys. Fluids 1974,I7 (lo), 1793. (43) Duker, M.; Koppenwallner, G. RareJed Gas DynamicsSymposium; Academic Press: New York, 1980;Vol. 2, p 1190. (44) Hagena, 0. F.; Obcrt, W. J . Chem. Phys. 1972,56 (9,19 793. (45) Hagena, 0.F. Z . Phys. D.: At., Mol. Clusters 1986,4, 291. (46) Palopezhentsev, S. A.; Yarygin, V. N. J . Appl. Mec. Tech. Phys. 1984,24 (2), 145. (47)Jaeger, H.L.; Willson, E. J.; Hill, P. G.; Russell, K. C. J. Chem. Phys. 1969,51 (12), 5380. (48) Bailey, A. B.; Price, L. L. Phys. Fluids 1985,28, 1583.
