THE UNIVERSITY OF MICHIGAN COLLEGE OF ENGINEERING Department of Aeronautical and Astronautica!l Engineering Aircraft Propulsion Laboratory DIGITAL COMPUTER ANALYSIS OF CONDENSATION IN HIGHLY EXPANDED FLOWS by James. L. Griffinr ORA Project 05169 under contract with; Aeronautical Research Laboratory Contract AF 33(657)-8867 Wright Patterson Air Force Base, Ohio administered through: OFFICE OF RESEARCH ADMINISTRATION ANN ARBOR August 1963

This report was also a thesis in partial fulfillment of the requirements for the degree of Aeronautical and Astronautical Engineer in The University of Michigan, 1963. Reproduction in whole or in part is permitted for any purpose of the United States Government.

ACKNOWLEDGMENTS The author is deeply grateful to Prof. Pauline M. Sherman for suggesting this topic, and for her guidance and assistance throughout the course of the investigation. The author wishes also to express his appreciation to Dr. Thomas C. Adamson and Mr. Kenneth Sivier for their helpful suggestions and to Mr. Leonard Harding for his assistance in programming the computer. Additional thanks go to Mrs. Eunice M. Sherbrook for typing the manuscript and to my wife for typing the rough copy. Professional Committee: Associate Professor Pauline M. Sherman Assistant Professor Martin Sichel

DIGITAL COMPUTER ANALYSIS OF CONDENSATION IN HIGHLY EXPANDED FLOWS SUMMARY The application of the IBM 7090 Digital Computer to the theoretical prediction of condensation in highly expanded flows is presented. The equations of the spontaneous nucleation theory of Frenkel (13) are combined with the steady one-dimensional diabatic flow equations for a solution of the expansion of a pure vapor. A digital computer program for the solution of these equations is compiled and presented. The theoretical prediction of the condensation of nitrogen is compared to experimental results and variations in specific heat, latent heat of vaporization, surface tension, and rate of expansion are investigated. The theoretical calculation is applied to metal vapors and the results for copper and zinc vapors are presented. The results indicate that the theoretical solution gives a reasonable prediction of the condensation in highly expanded flows. The degree of supersaturation increases with an increased rate of expansion and for a proper set of initial conditions "condensation free" flow is obtainable. The rate of expansion and the surface tension are the most critical parameters in the equations for condensing flow. Variations of specific heat and latent heat of vaporization show only minor effects on the end result.

TABLE OF CONTENTS Page TABLE OF SYMBOLS v I. INTRODUCTION AND SURVEY OF THE LITERATURE 1 II. THERMODYNAMIC FUNDAMENTALS 5 A. Isentropic Expansion of a Vapor 5 B. One-Dimensional Steady Flow Equations with 7 Condensation III. PREDICTION OF THE ONSET OF CONDENSATION 11 IV. NUCLEATION THEORY 16 A. Condensation Nuclei 16 B. Kinetics of Nucleation and Condensation 19 C, Limitations 21 D. Growth Estimates 22 V. ANALYSIS 27 A. Basic Assumptions 27 B. Supersaturated Expansion 29 C. Saturated Equilibrium Expansion 34 VI. DIGITAL COMPUTER PROGRAM 37 VII. CALCULATIONS AND RESULTS 40 A. Validation of the Program 40 B. Justification of Assumptions 43 Co Variation of Parameters 45 D. Application to Metal Vapors 47 VIII. CONCLUSIONS 50 REFERENCES 52 FIGURES 1-13 55 APPENDICES e A. Derivation of the Nucleation Rate Equation 68 B. Digital Computer Program 71 C. Input Data 89 iv

TABLE OF SYMBOLS A nozzle cross-sectional area A* throat area B, C saturation curve constants for In P B - C specific heat at constant pressure E energy per molecule f(n) distribution of n-sized drops g fraction of the mixture, by mass, which is in the condensed phase h enthalpy J nucleation rate per unit volume k Boltzmannvs gas constant L latent heat of vaporization M Mach number m mass of an atom or molecule mfh rate. of mass flow n number of molecules in a liquid drop N nucleation rate NA Avogadros number N number of n-sized drops n p free stream pressure p saturation vapor pressure for a plane surface of liquid Q heat transferred from the droplet to the vapor R universal gas constant r drop radius v

LIST OF SYMBOLS (continued) S drop surface area T temperature Ts surface temperature of the liquid drop t time U velocity Ufg change of internal energy due to condensation v specific volume V volume of one molecule W work of the formation of the droplet x station in the nozzle measured from the throat a accommodation coefficient a probability that one molecule leaves a unit surface area of a droplet of n molecules number of molecules striking a surface per unit area per unit time Bin 1 probability that one molecule condenses on a unit surface area of a droplet of (n - 1) molecules ratio of specific heats Ep fractional deviation from isentropic of the pressure change over increment Ax ET fractional deviation from isentropic of the temperature change over increment Ax molecular weight p density of the mixture a surface tension ~1/2 ~ nozzle diffuser half angle thermodynamic potential vi

LIST OF SYMBOLS ((continued) Subscripts and Superscripts' K )s refers to droplets of critical size K( )o refers to stagnation conditions () refers to vapor phase K) )liq or L refers to liquid phase ( )sak; refers to saturation point K)i refers to incremental step in which droplets were originally formed K)j. refers to incremental step presently being evaluated g AT refers to incremental step at temperature T which is being evaluated in the search for onset of condensation (() refers to an infinite surface of the liquid phase Go

Io INTRODUCTION AND SURVEY OF THE LITERATURE The phenomenon of the spontaneous condensation of a gas or vapor into the liquid phase has been the subject of numerous theoretical and experimental investigations for the past threeaquarters of a century. The condensation of an expanding vapor was first noted in 1887 by R. von Helmholtz- however, this phenomenon did not enter the field of aeronautics until about 1940 when its occurrence was noted in the supersonic wind tunnel. In a supersonic nozzle, the gas or vapor undergoes an isentropic expansion which results accordingly in a decrease in the temperature and pressure. Generally the saturation vapor pressure decreases with temperature more rapidly than does the static pressure of the isentropic expansion, hence saturation conditions are approached. If the rate of expansion is very rapid the flow will pass through the point at which saturation temperature and pressure are reached without any condensate being formed. The vapor is then in a ncn-equilibrium state, supersaturated, in which the existing vapor pressure is higher than the corresponding equilibrium saturation pressure. The first person to rnote that a supersaturated vapor can exist was R. von Helmholtz, who showed that saturated steam expanding through an orifice into the atmosphere would remain clear for some distance before becoming cloudy~ Stodola I1) furthered this by devising experiments on superheated steam which woud detec the onset of condensation and provide evidence of supersaturation0 1 c

-2In the late 1940's and the early 1950s an intensive series of investigations were conducted to determine the condensation of air, water vapor in air, and nitrogen. Many of these were conducted at the California Institute of Technology, and worthy of note are the investigations of Head (2) on water vapor in air, Willmarth and Nagamatsu (3) on nitrogen, and a summary report of GALCIT work by Nagamatsu (4). A series of investigations of air in high speed wind tunnels conducted during this era include the works of Stever and Rathbun (5), Wegener and Reed (6), McLellan (7) and Bogdonoff and Lees (8). There was a lag of interest in this area in the middle to late 1950's. However, with the advent of the nuclear space propulsion system and the use of metallic vapors as the working fluid for hypersonic wind tunnels, the interest has been renewed. The Office of Naval Research sponsored the recent review of existing theories by Courtney (9) in an effort to stimulate this interest. Many of these theories are built around the spontaneous self -nucleation concept, and they base the collapse of the supersaturated state upon the determination of the rate of droplet growth to the critical drop size. The complexity of the nucleation phenomenon requires that many simplifying assumptions be made in order to evaluate the problem by hand calculations. Courtney (9) attempts to obtain a more complete evaluation of the nucleation theories by carrying out the calculations on a digital computer. A recent study of air in hypersonic wind tunnels by F. L. Daum (10) indicates that, as yet, there is no proven available theory for predicting the amount of supersaturation expected in a hypersonic nozzle from which may be determined the point of onset of condensation.

-3Until recently none of the investigations had included the problem of metal vapors. The investigation by Hill, et al, (11) in 1962 adapted the existing theory to the case of metal vapors and provided results which indicate that this is a valid extension of the theory. The development of the nucleation equations and flow equations did not prohibit their application to metals, but since the need had not existed the application was not attempted. The authors noted have not mentioned any adaptation of the digital computer to the hypersonic flow problem. Hill suggests that this is the next logical step and that a computer program should be used to investigate the effect of variations in surface tension, nozzle length, area distribution, inlet conditions, and other parameters as well. Purpose of the Present Investigation At the inception of this study there was no known application of condensation theory to the case of metal vapors. The Hypersonic Wind Tunnel Group of The University of Michigan was undertaking the experimental investigation of the condensation of metal vapors in highly expanded flows, and thus arose the need of a method for the theoretical prediction of the onset and magnitude of condensing flow. Of greater importance was the need for an understanding of the effect a variation in any of the various parameters might have on the end result. Thereby arose the need for a digital computer program from which numerous results can be obtained for multiple combinations of the various parameters. The extent of the present study is too (a) review the literature to determine which of the more promising theories will most readily adapt to the IBM 7090 digital computer,

(b) write the program and test against some known experimental results to validate the method, Kc) obtain theoretical results for several metal vapors under reasonable hypersonic test conditions,'d) vary surface tension, specific heat, heat of vaporization, and nozzle parameters to show the magnitude of influence on the problem. The first approach was one toward a general problem with all quantities variable and with multicomponent working fluids. This rapidly exceeded the capabilities of the computer, and of the author, and dictated the simplified approach which is employed. The system will be developed for onehditnensional expansion of pure vapors which obey the perfect gas law and which have a constant ratio of specific heats. Initially the specific heat, latent heat, and surface tension will be taken as constant and later will be varied to observe the effect0

Ito THERMODYNAMIC FUNDAMENTALS A, Isentropic Expansion of a Vapor The assumption that the vapor follows the perfect gas law is, in general, well justified and the pressure and temperature can be related by. T:7 1 p11~~~~~ T 1~ ~~(2 1i) p T 0 0 Any vapor whose latent heat is large, i. e., L > C T, will approach p saturation as the temperature is decreased by an isentropic expansion, and if the expansion is rapid enough there will be a thernmal lag so that a supersaturated region will be reached, The Clausius-Clapeyron equation can be used to relAae thie variation of the saturation vapor pressure with temperature to tom latent lieat of vaporization. dpoo L(T) (2. 2) where po is the saturation vapor pressure corresponding to a particular temperature and for an infinite droplet radius, L(T) is the latent heat of vaporization, and v and vL are the specific volumes of the vapor phase and liquid phase, respectively. Since tile liouid icci(eu)ies a very sj nlal- -)ortion of tile total volull-ie, we will assume v - v v and v Assuming L - constant, v L v v p equation (12. 2) becomes V

dPm L p Li.o (2 3) dT vT 2 R v T where, p. the molecular weight of the vapor R - universal gas constant. In order to compare the slope of the isentrope with (2. 3) we will differentiate (2. 1) logarithmically2 dp: d/ dT (2~ 4) p -1 T from which, dp p p(25 dT R T by virtue of v C /C and R/ tA C - C o The ratio of slopes is p v p v therefore dPoo dp L (2. 6) dT dT CT Thus for the case of L/T > Cp the saturation vapor pressure decreases with decreasing temperature more rapidly than does the vapor pressure in an isentropic expansion, and the supersaturation ratio, p/po, continually increases until saturation and, eventually, supersaturation are obtained. Head (2) gives slope ratios for several vapors and shows a typical value of 2 for nitrogen nozzles operating at normal reservoir conditions. Similarly, he shows that if a vapor had a small latent heat of vaporization, L/T <Cp, an isentropic compression is necessary for the approach to saturation conditions.

-7Upon reaching the saturation point the fluid has a choice of two routes by which to continue the expansion. One, the saturated equilibrium expansion, involves a gradual process of condensation which commences at the saturation point and maintains thermodynamic equilibrium throughout the expansion. The other, a supersaturated expansion, is characterized by a delay in the onset of condensation to some point downstream of the saturation point at which the vapor temperature is well below the corresponding saturation temperature. Before continuing further with these alternatives it will be necessary to develop the diabatic flow equations which account for the effect of the condensate. B. One-Dimensional Steady Flow Equations With Condensation The equations have been derived by many authors; however, we are interested here in a specific case which includes the following assumptions: 1. The vapor phase is a pure vapor and may be treated as a perfect gas. 2. The condensed mass is a liquid, is uniformly distributed throughout the gaseous components, and has the same speed and temperature as the streamno 3. The volume of the condensed phase is negligible compared to the total volume. 4. The nozzle flow is frictionless and is without heat transfer. 5. The saturation curve for the vapor may be adequately approximated by the Clausius-Clapeyron equation. 6. The latent heat of vaporization and the specific heat at constant pressure remain constant.

Assumption 2 requires some disecussion, in that the crystallization of the liquid drops or the sublimation directly to the solid state is probable at the lower temperatures. In the area concerned, the available experimental data on the properties of the vapor and the liquid phases is extremely limited, and for the crystalline phase it is practically non-existento Highly extrapolated curves are the only sources, in many instances, so the liquid assumption will be as good as the available data. The zero drag approximation is supported by Stever (12) in citing lnumerous investigations which verify high fractional condensation rates. From the equations cf Frenkel (13) it is noted that a small critical drop radius is necessary for a high condensation rate. In high speed flow the drops are in the nozzle for only a short period of time, so there is only time for a small increase in size due to growth. Thus, for high speed expansions with supersaturation, the drop size will be small, and the logical assumption is that the condensate and the vapor travel at the same velocity. Only the major steps i, the development of the flow equations for a pure vapor are presented in this section. A more complete development is presented in Wegener (14). The continuity equation can be written in its usual form, since the volume of the liquid is negligible and the vapor and drops travel at the same velocity. With p taken as the densiity of the mixture and rh as the total mass which passes a given location in unit time~ pUA = m = constsant (continuity) dp +.dA 0 (2~ 7) p U A As a result of the assumptions 2 and 3 the equations of momentum and state ares

-9(momentum) dp = UdU (2. 8) P (state) dp dp dT dg (2. 9) p p T (1 - g) From the first law of thermodynamics, with no external heat sources, the energy equation is d (h + U2/2) 0. (2. 10) The density of the mixture, p, used in the above equations is defined by the equation, p t pv + P' The density p L is the mass of liquid condensate per unit volume of the vapor. This requires an assumption that the liquid droplets be uniformly distributed throughout the vapor so that the densities of the vapor and of the liquid may be referenced to the same volume. g is the mass fraction of the working fluid which has condensed into the liquid phase. By use of the assumption that the vapor and the liquid are at the same temperature, and for constant C and L equation (2. 10) becomes, d (U /2+ C T- gL) = 0 (2. 11) p and (energy) UdU + C dT - L dg 0 (2. 12) p By introducing the sound speed by the perfect gas formula in 2 R terms of vapor density, A = T, an alternate form of (2. 12) is - 1) M2dU dT L dg0 (213) (UM T C T p

-10The four equations, (2. 7), (2. 8), (2. 9), and (2. 12) involve the six unknowns, p, p, T, U, g, and A. For this particular problem, it will be assumed that g and A are determined from other sources and p, p, T and U are obtained from the above equations. Area will be known for a given nozzle, or the area distribution must be determined from a given Mach number or pressure distribution. The question of determining the condensation, g, depends on the assumed flow conditions. If it is assumed that thermodynamic equilibrium exists at all times the Clausius-Clapeyron equation can be used. This method will be developed in Section V. However, in the event supersaturation takes place, the determination of g is not so simple and nucleation theory must provide the missing equation. After devoting one section to the determination of the onset of condensation, and one to nucleation theory, the flow equations will be married to the condensation equations for a complete solution.

IIIo PREDICTION OF THE ONSET OF CONDENSATION A review of the literature shows that indeed there is no proven available theory which is valid for predicting the onset of condensation in a supersaturated flow. A solution of the equations for predicting this point always requires an assumption as to the amount of condensate formed, the droplet formation rate, or the deviation from isentropic flow. A marriage of the flow equations (thermodynamic) and the condeir sation equation from nucleation theory provide a means of corMput.ng the flow properties for a given area distribution. A plot of a flow paramr, eter during the period of condensation shows a sharp deviation fro m a plot of the same parameter undergoing isentropic expansion. Experimental investigations have determined that static pressure, p, and pitopressure, pi, are strongly affected by condensation and a jump in these quantities usually is used to detect the onset of condensation. Here we are defining the point of onset of condensation as that poinit along the flow at which the initial condensate formed is sufficient to cause the computed value of p or T to differ a detectable amount from its corresponding value for the isentropic expansion. The complete solution could be solved for very small steps starting at the saturation point and continuing throughout the length of the rizzle, From a plot of actual pressure, p, compared to the pressure distr-ibution for an isentropic expansion the point of onset of condensation cou!d be determined. However, this would require an excessive number of calculations, and for certain combinations of nozzle parameters and stagnation conditions the vapor may traverse the nozzle without ever condensing and several minutes of computer time would be wasted. -11-~

-12Therefore, it would seem expedient to establish a criterion for estimating the point of onset of condensation and then to continue the stepwise computation of the full set of equations from this point onward. Since the expansion of the vapor, in which there is no condensed liquid, will be assumed isentropic, any deviation from isentropic conditions should serve as an indication of condensation. As previously stated, the static pressure is sensitive to condensation and this will be used as an indicator. From the calculations of Hill (11) the static temperature, T, also appears to give a strong indication. At certain times it is quite possible that the temperature might give a better indication than the pressure, and this factor now will be investigated. The flow equations from thermodynamics, (2o 7), (2. 8), (2. 9), and (2. 13), relate isentropic conditions when there is no condensate present g = Ag = 0, and relate actual conditions when the vapor is condensing. A proper combination of these equations should produce a parameter which indicates the deviation from the isentropeo Since we are interested in temperature and pressure we will recombine the flow equations to obtain expressions for pand dT dp dU dA -p + dU dA 0 (2, 7) p U A dp = - pUdU (2. 8) dp dp dT dg -- + - (29) p p T 1 -g - 1) M2 dU dT L dg 0 (2 13) Substitute equations (2. 7) and (2. 13) into (2, 9) to obtain

-13dp _ dU dA 2 dU L dg dg p — U A 1)M C T (1 -g) dp 1)2 dA I L 1 [1+ 1) M2 + - Id (3 1) [pM -iU A1 -T g) _ (1 - g) ( - 1) - adg (3, 1) p from (2. 8) and the equation of state, p = p (1 - g) T, we obtain dU (1 - g)dp U- y M2 p Insert this into equation (4. 1) and rearrange dtp mndp (1ta 9) 1) 2r L 1 d d-g -~ M2'~ dA! L d: ~P [+ 2 A +C T 1-g pM p. P dp WM dA L 1 [ M2 (I - g- (1 - g)(- 1) M2] (. p, / For this particular problem we are interested only in the point of onset of condensation. Therefore, we assume isentropic conditions exist up to the point and the only condensate present will be the dg formed over the incremental step at which we are considering the onset to take place. The calculations will be accomplished in a step by step manner on a digital computer, so we write all increments as A and set g -0 to get 2 Ap yM AA L( P 2 A +C T) From (3. 3) it can be seen that over steps in which there is no condensate (Ag = 0) the stepwise change in pressure will be isentropic. However if Ag / 0 then the step change in p will differ from isentropic. We therefore rewrite (3. 3) as

-14P 2 A1 A Ap yM AA 1I AgA (3 4) p 2 IA AA C T M-p We now can define p = AA C T g 3 5) as the fractional deviation from isentropic of the pressure change over an increment Ax. The magnitude of E which can be considered a significant deviation is yet to be determined. Now we examine a similar analysis of the temperature as a possible indicator. From equations (2. 13) and (2. 8) we obtain dT 1 1 dp L ( -g)-+ dg (3. 6) P By substituting the value of Ap/p from (3. 4) the equation becomes AT Y- YM2 LA A {L A 1+ + 9g T 1M2 1 A AA C T ) C T M2- p P AT (M - 1) M2 AA A2 -~1+ -_ + 2AgC T 1M2 A1 A C[ Y- T1 Again we define T [ J1 ~~~~~$Ag (3.8) T AA y - 1 i" CiT as the fractional deviation from isentropic of the temperature change over Ax.

Note: A comparison of (3. 5) and (3. 8) shows a difference only of the coefficient of L/Cp T Thus as - 1/M21 )-, i.e., as M-. 1 ETapproaches the value of E p However, an examination of eT and E for given conditions of L/CpT, AA, Ag produces. (a) for condensation in the supersonic portion (M > 1) ET > epp (b) for condensation in the subsonic portion (M < 1)-: e > CT" As an example we consider reasonable values of L/C T z- 3 and M 10'' A A the result gives E (2) Ag and e (9. 5) l Therefore. a p AA T AA dual criteria will be selected so that eT will determine the onset of condensation when M > 1 and e will be used when M <' 1.

IV. NUCLEATION THEORY A. Condensation Nuclei In the process of condensation, the saturated state is referenced to a plane liquid phase and the vapor is considered to condense on a flat liquid surface or on a "cold" container wall such as the moisture formed on the outside surface of a glass of ice water. However, in highly expanded flows the liquid surface is not present and the "cold" wall is not available, due either to the "hot" boundary layer surrounding the flow in a nozzle or to the total absence of a wall in the case of the free jet. The vapor then must look for other surfaces on which to condense. These surfaces appear either in the form of impurities, such as particles of dust, etc., or as small drops formed of clusters of molecules of the vapor itself united by statistical collision. Obviously, for a truly pure vapor there would be no foreign nuclei; however, dust particles in clean air have been estimated upward of 103 particles per cubic centimeter, so this possibility must be considered. Head (2) has shown that these foreign nuclei do play an important role in the onset of condensation. However, for highly expanded flows where the stay time is very short, Stever (12) has shown that the high rate of condensation could not be caused by condensation onto the foreign nuclei alone. Stever (12) quotes an example from Oswatitsch's study of water vapor which considers 10 foreign nuclei per cubic centimeter traveling in a gas for a distance of 10 centimeters where there is super cooling of 300C. A slightly supersonic velocity of 3. 3 x 104 centimeters per second would traverse this 10 centimeters in 3 x 10 seconds. Using the faster growth formula for droplet growth gives a -16

-17radius of growth of about 3 x 10 centimeters in this distance. The amount of water condensed over this span then is 10 grams per cubic centimeter, which releases an amount of heat that is very slight and which has an insignificant effect on the flow quantities. The actual condensation process has a marked effect upon the flow properties; therefore it must be concluded that the number of drops far exceeds the number of foreign nuclei. We must then resort to the spontaneous formation of the droplets onto which condensate growth occurso There is an opening in the above argument, in that the growth rate onto the foreign nuclei could be grossly in error. However, Head (2) found droplet formation rates on the order of 10 drops per cubic centimeter per second in his investigation of water vapor; thus, the high particle formation rate tends to support Stever's argument For this analysis it will be assumed that the condensation phenomenon is one of spontaneously formed nuclei and is triggered by the attainment of critical size droplets onto which the supersaturated vapor condenses. The concept of the critical drop size as developed in Stever (12) stems from the energy required to evaporate a drop of liquid which contains n moleculeso E(n) =nE - 47r a (4o 1) In (4. 1) Ec is the energy required to remove by evaporation a single molecule from an infinite surface of liquid and o is the surface tension. To evaporate one molecule would require 2cV AE(n) = E - q (4. 2) Go r

where Vliq is the volume of a molecule of the liquid. An application of Boltzmann's law gives vapor pressure AE kT p =K e and the equilibrium vapor pressure for a plane surface E co ekT Combining these gives the Thompson-Gibbs equation 2Vl or rkT P e l4q 4. 3) P0 Expressing the molecular volume Vliq as pNAPL and rewriting it as a log function, (4. 3) becomes in 2p/p c2 _~ 2(4, 4) 00 rpLkNAT rpLRT 4) Note that with the occurrence of r in the denominator of equation (4. 4) the equilibrium between the vapor and the liquid drop is unrstable. For a vapor which is saturated with respect to a given drop size) if one molecule sticks to the drop the radius increases and the existing vapor pressure is greater than the new equilibrium vapor pressurethus, condensation will continue. The reverse is true for evaporation, in that if one molecule evaporates the existing pressure is lower than equilibrium and evaporation continues. Thus arises the concept of the critical drop radius, r'* for which r < r* the droplet vaporizes and for r > r* the droplet grows. Solving (4. 4) gives

-19r*= 2((4~ 5) PLRT In p/po0 Hill (11) points out that in a constant environment, if the critical radius was used to specify the initial drop size then the growth rate would be zero and there would be no condensation because at critical size the decay and growth are equally probableo Hill bases his use of 1. 3 times the critical radius on Becker and Doering's conclusion that at that size the probability of decay is close to zero. In this analysis the vapor will be flowing and any critical drops formed wall not remain in a constant environment very long since the quasi-steady steps for the calculations will be very small~ Therefore, the commonly used r* of equation (4. 5) will be used as the initial drcp size. B. Kinetics of Nucleation and Condensation As implied previously, the spontaneous nucleation approach will be the only one considered in this analysis. There have been volumes written in the development of this area and no attempt at reproduction will be made here. The general approach and final result of the most frequently referenced authors will be stated, however. In highly expanded flows the state of the vapor can change very rapidly from the unsaturated state to the highly supersaturated state. Frenkel (13) treats the equilibrium distribution of drops in an unsaturated state as a dilute solution of different substances in the vapor as a solvent. By classifying the different solutes by the number of molecules, n, they have in a drop, he obtained the expression for tihe number of drops containing n molecules, N - N e [4 6W na v

-20where Nv number of molecules in the vapor state ~liq' v = thermodynamic potential of single molecule 4 r 2a = surface energy term. This droplet distribution is based on the assumption of equilibrium (Ov < bliq); however, it is generally assumed that the same distribution can be used for supersaturated conditions {(liq < 0v), Note that for the supersaturated condition the number of drops increases with increasing radius and the distribution is unstable. Thus, it becomes essential to determine the rate at which the distribution changes from the stable distribution of the unsaturated state to the unstable distribution of the supersaturated state. This requires a determination of the rate of formation of critical sized drops, J, in drops per cubic centimeter per second. Head (2) states the equilibrium equation of Volmer AW 47 r *2 kT 3kT (4 7) J=Koe — Koe (4~7) where AW is the total work of the formation of the droplet. The determination of the constant, KS has been developed by several authors, including Volmer, Becker and Doring, and Frenkel, All assume a quasi-stationary phase transition process in which the number of molecules which condense are instantly replaced in the vapor such that the number of vapor molecules is maintained constant. Volmer considered only the number of drops reaching critical size, the number of vapor molecules striking the surface of the critical size drop per unit timed and the number of molecules in a drop of critical sizeo

-21Becker and Doring derive the nuclei formation frequency purely from kinetic considerations and take into account the evaporation from the surface of the drop, as well as the condensation onto it. This approach is considered by most authors to be better than others; however their solution is rather difficult to apply. Frenkel uses the method of Zeldovich to integrate the differential equations of Becker and Doring and gets a much simpler derivation of an equation that is almost identical to that of Becker and Doringo The ratio of the result of Frenkel to that of Becker-Doring is (n*)3 and this is considered minor. n* is the number of molecules in a critical sized dropo A partial derivation of the nucleation rate equation is presented in Appendix A. This derivation is necessitated by an apparent error of rT in the original equation presented by Frenkel (13). For this analysis, equation (A. 7) will be applied as the nucleation rate equation. p2T 1L e2r 3kT4 p e 4. 8) L A where J -- number of critical drops formed per unit volume per unit time. C. Limitations 1.' Generally a as a function of temperature is not known and this is a large factor in the result. Surface tension appears to the third power in the exponent; hence a small error will have a large effect upon JO

20 The rapidity of expansion will show marked deviations from the assumed quasi-stationary condition, 3. The condensing flow may have sublimation directly from vapor to solid which is not considered. 4. At high degrees of supersaturation the drop radius will be small, and Head (2) has shown that the variation of Ac with r is large for drops of few molecules. Stever and Rathbun (5) have included this effect into the Frenkel equation to obtain r* 8 IF, rudr - 1 r*21 (kT)2 r' m r r = Bogdonoff and Lees (8) have an approach which is entirely different from that of Head and Stever and which produces quite different results. Until these theories have been developed further, it is not practical to include either theorem in this analysis. From the previous discussion of the assumptions made in the formulation of the nucleation theory and of the limitations above, it is apparent that the best one can expect in an analysis such as this will be an order of magnitude estimate of the condensation parameters and a qualitative analysis of the flow conditions. Do Growth Estimates Heretofore, in order to calculate the formation of the critical drops, we have rconsidered the temperature of the liquid drops to be tpie sinme as the surrounding vapor0 In actuality, the impinging and vapcr'zing (f4

-23molecules onto and away from the droplet leave a surface temperature that approaches the saturated equilibrium value for the surrounding vapor. This temperature difference allows for growth of the particle. The growth of small drops, i. e,, r < < mean free path, occurs by the molecular transport of molecules onto the surface, some of which condense, and by the conduction away of latent heat of vaporization by the evaporation of a fewer number of molecules which leave at higher kinetic energies. Oswatitsch laid the basic groundwork for growth in a flowing vapor. Stever (12) follows the method of Oswatitsch to develop an equation for the growth rate of the small drop. The development presented by Stever makes use of the kinetic theory of free molecule flow (See Patterson (15)). Consider a drop of radius r, temperature Tliq, and saturation pressure (ps,)liq with a surrounding vapor of pressure p and temperature T. The heat transfer from drop to vapor per unit time is Q ~ pi (Tli- T)a where a is an accommodation coefficient which varies from 0 to 1, depending on the vapor and liquid surface. Also from kinetic theory the heat transferred to the drop surface due to condensation is Lfp v/kT/2r m where p fikT/27r m is the amount of vapor (in the form of molecules) impinging on a unit surface and f is the fraction which condenses. The energy being transferred to the surface is equated to the energy being transported away to obtain an expression for f;

-242a (Tliq - T) LT the mass MD of the drop is then increased by MD 1 2 c / (Tli2 _ L,7 r a-p (T.- T) dt L V mT l T 4 3 2 Since MD 4 r3 p and dMD 47r PL dr, d 2rr a p (Tliq T) dt f, LPL m Nl If the calculation is made in quasi-steady steps, time can be expressed as a function of step length and Ar 2- L L (T T) (4. 10) rL pL mT liq U where the accommodation coefficient, ac is a constant, approximately equal to 1 and Tliq is the temperature of the surface of the liquid drop. The accommodation coefficient a is expressed by Kennard (16) as a (E E ) ( gEi - E ) 1 w i r where E. - Energy brought up to a unit area per unit time by the impinging molecules Er - Energy carried away from the molecules leaving a unit surface per unit time E - Energy that the departing molecules would carry away if they carried away the same mean energy per molecule as does a stream issuing from a gas in equilibrium at the surface temperature.

-25Patterson (15) lists coefficients of the range from. 8 to 1. 0, and Stalder (17) found that for free molecule flow an accommodation coefficient of 0. 9 checks well with test results. Thus we arrive at the statement, "approximately equal to 1.' This only applies for cases where free molecule flow can be assumed to prevail. On the other hand, where r approaches the magnitude of the mean free path the heat removal I.s controlled by normal gaseous heat conduction. For this large drop Stever (12) gives the growth equation for period of time, t, /2k (Tiq T) or = X L 4 11)p Other authors have made similar analyses, and Wegener (14) records the quasi-steady solution of Buhler (18) which includes an express sion for the condensation coefficient. 1 dr p 2(T sT) dr p viL i PI (2zr/ 15 _(4. 12) Pdx L pLM r2LL (-L where Ts, the surface temperature of the drop, may be chosen as the equilibrium saturation temperature of the vapor with respect to the drop. As in the previous case, Wegener states (4. 12) is valid only for small drops and the equation recommended for large drops is merely a revised form of (4. 11). Hill (11l) notes that typical calculations of condensation nuclei indicate very small sizes which are generally much smaller than the mean free paths; thus only the free molecule equations need be considered.

-26His analysis produces a growth rate expression of dr _ 3p [RTD - RTJ (4. 13) d 21 Tg-L 2RT in which Ufg.. change of internal energy due to condensation. The term TD is droplet temperature, and Hill presents calculations to show that the rate of drop temperature change is enormous compared to the rate of growth. For water vapor the temperature approaches wet bulb temperature in approximately the same time that it takes the number of molecules in the drop to increase by a factor of 10. Hill uses this in expaiatiorn of the assumption used by Oswatitsch that the drops always have wet bulb temperature during the growth period. After a displacement, the drop temperature will rapidly approach the temperature at which d T/dt - O at a much greater rate of change than dr/dt. Hill uses this to form an expression for TD which depends only upon the instantaneous values of r, p and T. e3xRT 1t1 -4p PLRTr, TD TtE -9, X4 14) 3 RT o p RTDrTi/ T where t is the fraction of impinging molecules which condense. From the investigations on water vapor, nitrogen, air and a few metal vapors it appears that for highly expanded flows the small drop size is the general rule. The calculations will then be initiated by applying (4. 10) as the growth rate equat;ion and choosing drop temperature as the equilibrium saturation temperatiure of the vapor. The accommodation coefficient in (4. 10) will be varied to study the effect and a wll be placed into the computer program as a data input so that accurate values, when knrown, can be used for the calculationr.

V. ANALYSIS The basic equations are established in Sections II through IV. The proper combination of these equations produces a theoretical solution for the hypersonic expansion of a pure vapor in which a supersaturated state occurs. Another solution is obtained for which there is saturated equilibrium expansion in lieu of the supersaturated portion. Both solutions, however, stem from the same basic model. The basic model considered here consists of assuming. 1.o "One-dimensional nozzle" with geometrically specified inlet and diffuser. 2. A pure vapor with specified stagnation conditions and known: constant specific heats, latent heat of vaporization; surface tension, for infinite radius as a function -of temperature; and saturation curve. 3. An isentropic expansion to t!he point of onset of condensation and then diabatic flow through the region in which condensation is taking place. A. Basic Assumptions In addition to the assumptions of the basic model, this analysis also embodies the multitude of simplifying assumptions which have been included in the derivation of the basic equations. The type of equations and the specific assumptions made are, 1. Flow equations a. Mass flow is constant b. Flow is one-dimensional and steady state co Volume of condensed phase is neglgihble when compared to the total volume 27=

-28d. Nozzle is frictionless with no heat transfer across nozzle wall. eo Condensed mass is liquid, is uniformly distributed throughout the gaseous components, and has same speed and temperature as the stream 2. Prediction of onset of condensation equations a. This embodies same assumptions as flew equations bo No condensate is formed prior to the predicted point of onset of condensation 3. Nucleation equations a. Drops are formed by spontaneous nucleation only bo Only drops reaching critical drop size continue to grow c. Number of molecules of vapor is maintained constant d. Saturation pressure corresporeds stc aturat~orn pressure of a droplet of infinite radius e. The ordered velocity of the partU(~.es is neglected f. Equilibrium values of the probabilities that one molecule will leave from or condense on a unit surface and equilibrium particle distributions are assumed valid for the non-equilibrium case 4. Growth rate equations -a. Drop radius is much less than mean free path b. Accommodation coefficient, a, is approximately equal to 1. 0 c. Drop temperature is assumed equal to the saturation temperature that corresponds to the pressure of the vapor

-29B. Supersaturated Expansion The supersaturated expansion passes the saturation point without any condensate being formed; therefore, a four step solution is required: 1. Isentropic equations from Shapiro (19) and NACA Report 1135, are used to calculate the expansion of the vapor to the point of condensation. P -PoT (5. 1) 0 M 1I T(- (5. 2) p=pO1 +Y 1 M21 (5. 3) y+ 1 A* 2 21M 2K2 1) A = -A 1 + i M (5. 4) U=[2 Cp (T- T)]2 (5. 5) At any particular value of T the above equations are readily evaluated for given p0, To, po, y, Cp and A*. The station at which these values apply is then determined from the given nozzle geometry, 2. The saturation point is established by an iterative balance of the saturation curve equation with the isentropic flow equation.

-30The saturation curve can generally be expressed in the form In poo - B - X in which B and C are constants whose value depend upon the choice of units for p and T. The equations which require joint solution are T (5.6) sat B - In psat T =1 Psat = ( Tot) (5, 7) The isentropic relations above are then utilized to calculate M P T A p and U sat' sat' sat' sat Psat' sat" 3. The point of onset of condensation is determined by taking steps, AT, of supersaturation and computing the flow parameters from the isentropic relations for each successive step. The isentropic values of each incremental step in temperature are introduced into the nucleation equations for critical drop size and formation rate. r* L Tlnpp (4. 5) T PLRTInp/Poo 4 - 4cr*2 Jr ke T (4. 8) TkTJPL NA Note that choice of a very small AT will cause the first step to remain very near the saturation point and the resultant P/Poo - 1 will cause r* to approach an infinite radius. r* falls off very rapidly, however, as the degree of supersaturation increases, and for most cases a choice of AT > 50k is sufficient.

-31The growth rate is zero for the critical droplet at the instant of formation; therefore the nucleation rate for this increment is NT T JT AT (5 8) The mass fraction of condensate formed in this increment is 4pL r 3 AgT L-3i N- T N T (5. 9) The epsilon equations from the onset of condensation development, (3. 5) and (3~ 8), are now used to determine if this amount of condensate is sufficient to cause an appreciable deviation from the isentropic expansion. e = 1 Ag (3. 5) p AA C T P E j - (3.8) T AA l- 1 C T If the e is less than the assumed critical value then the same calculations are carried out at step T + AT until the criticall value of E is equalled or exceeded. If exceeded, a bracketing procedure is employed to improve the est mat~e All quantities computed for Ag at values of e eCritical are discarded, and it is assumed that there has been no condensate formed prior to e E ec tical' The choice of Ax in equation (5. 8) will effect directly the magnitude of N and of Ag in equation (5. 9). AA in equations (3, 5) and (3. 8) is also a function of Ax so that in the evaluation of E the effect of Ax is nullified~ Calculations have been

-32performed for values of Ax ranging from. 01 to 1. 0 cm to verify that the quantity has no effect on the prediction of the point of onset of condensation. Since the Ax that will be used for incrementing the calculations through the condensing portion of the nozzle will be a required input to the computer program, it is advantageous to use that same Ax in evaluating (5. 8). 4. The condensing portion of the flow requires a joint solution of the nucleation and growth equations and the diabatic flow equations. The calculations are performed for increments of Ax starting from the onset of condensation. Values are estimated of p, T and U at the point i and Agi is calculated from the equations r RT lnpp (4, 5) i PL PRT n p2/ p J. = ( P ) 1 ( ~A~u)1/2 3kT (4. 8) N =J A. A x, (5.8) 1 11 1 A. j (T T)- (4. 10) J Lp L s LT r..=r. r5.10) rij ri (j - 1) j 10) j-1 Ag-N EAZ N r. A r + N. r* (5. 11) L~gj =g 1 1 Ig g.=gi 1 +Agj (5. 12)

-33The subscripts i and j are numbered increments of x starting from the condensation point. i acts as a label for each group of droplets of a particular size and denotes the particular increment in which these drops originated as critical sized drops. j denotes the increments presently under consideration. Thus a designation r13 denotes droplets which were formed as critical drops in increment 1, have undergone the growth Ar2, and are now undergoing the growth Ar3. Using the above approximate values of Ag and g and the unknown values of AA and A the equations from Section II are solved for p, T, p and U. p U A @+ AU~+ =AA 0 (2~,7) Ap _ UAU (2.8a) P R (1 - g) R T Ap Ap AT Ag2 p p -- ( (-g9) UAU + C AT- LAg =.: 0 (2. 12) The solution of the flow equations is accomplished by the following iterative procedure: a. assume AU and compute U. - U. +H AU b. solve (2. 12) for AT and compute T. - T. + AT solve 2 7) f c. solve (2. 9) for Ap/p d. solve (2. 9) for Ap/p e. solve (2. 8a) for AU and check with value from step (a). Improve the assumed AU and continue until a solution is obtained.

-34These corrected values of p, T and U are now used to resolve the equations for Ag, and the process is continued until a set of flow parameters which satisfy both sets of equations is obtained. This iterative procedure is continued for each increment Ax up to a specified ax, usually the end of the nozzle. This set of equations could be evaluated from the saturation point onward in order to obtain a more accurate prediction of the onset of condensation. It is obvious, however, that the magnitude of the calculations involved is formidable; and, the slight error induced by attempting to predict the point of onset of condensation, does not justify the refinement. This is especially true since the approximations in the basic equations are so gross. C. Saturated Equilibrium Expansion The saturated equilibrium expansion is a gradual process with condensation starting at the saturation point and equilibrium flow being maintained throughout. This occasionally occurs for slow expansions of certain vapors and for cases in which the impurity content is high. In this analysis the saturated equilibrium expansion is presented only as a comparison with the supersaturated case. The calculations are greatly simplified in that the mass fraction of condensate can now be specified by the Clausius-Clapeyron equation dpTLM p (2. 3) dT R 2

-35This equation and the flow equations (2. 7), (2. 8), (2. 9) and (2. 12) usually are integrated and expressions are written for each parameter (see Wegener (14) or Daum (10)). By proper manipulation and combination the integrated equations can be reduced to one involving g. C T T g LT - n T t+ ( T (5.13) CL Tsat T s However, in this analysis the flow equations are programmed in the incremental form for the supersaturated case, so it is convenient to write (5. 13) in the form C T. T. Ag = L Ti T + -l 1 (5 14) L i Tsat sat i and to use this for the Ag expression in a solution similar to that used for supersaturation. 1. The saturation point is established in the same manner as for the supersaturated case. Small increments of length, A x, are started at the saturated point. 2. An incremental increase in the mass fraction of condensate, Ag, is assumed for the step Ax. With assumed value of Ag and known value of A the flow equations are solved for p, p, and T. Ap AU AA0 (2.7) S + / + A=0 (2 7) p U A AAp _ _ UAU, (2. 8a) P (1 - g) T

-36Ap A ~AT Ag + - (2.9) p p T (1 - g) UAU + C AT - LAg = O (2.12) Equation (5. 14) is solved to check and correct the assumed value of Ag and an iterative procedure is established to obtain a satisfactory solution of the five equations. Again this is continued in a step by step process until the solution is complete.

VI. DIGITAL COMPUTER PROGRAM A major task in this study has been the adaptation of the mathematical analysis of Section V to the digital computer. The description contained herein is presented for the future user of the program and in no way covers the myriad of detours encountered between the input and output statements. The Computer Center of The University of Michigan has developed a Michigan Algorithmic Decoder (MAD) computer program which compiles much faster than the commonly used FORTRAN program. The computer program for this investigation has been written in the MAD language and can be translated into FORTRAN if desired. The IBM 7090 computer program for the condensation of a pure vapor is presented in Appendix B. The program has been subdivided into a main program and six subroutines. The subroutines enable easier modification of the program and provide a greater flexibility in its use. Statements as to the function of each of the calculations, the input parameters, and the obtainable results are provided in the program at the beginning of each subroutine. Appendix B also includes a list of the program symbolic names assigned to the parameters used in the equations of Section V. A brief summary of the subroutines follows: 1. Isentropic flows (IFLOWS)-Given T, this program solves the isentropic flow equations for p, p, U, M, and A. 2. Condensation Flow (CFLOW)- For given g and A, this program solves the diabatic flow equations (2. 7 through 2. 12) for p, T, p, M, and U. -37

-383. Nucleation Theory I (NUCLEl)-From nucleation theory equations (4. 5, 4. 8, 5. 8 and 5. 9), NUCLEI computes the condensate formed due to the nucleation of critical sized drops only, i. e,, growth is not included in this portion. This subroutine is used to determine condensate formed for that portion of the main program which determines the point of onset of condensation. 4. Nucleation Theory II (NUCLE2)-This subroutine evaluates the growth equation (4. 10). NUCLE2 is used in conjunction with NUCLEI to determine the amount of condensate present for all calculations downstream from the point of onset of condensation. 5. Nozzle (NOZZLE)- The geometry of the nozzle is placed in this subroutine, from which A, AA, or x can be determined. 6. Vapor Parameters (FRHOL, FL, FSIGMA)-The parameters PL' L and a are placed in this subroutine as functions of T or as constants and are evaluated as called upon by the Main Program. The Main Program consists primarily of a proper utilization of the subroutines. The operations by the Main Program appear in the following order: 1. Read in the input data, 2. Use IFLOWS to expand isentropically from stagnation conditions, 3. Match isentropic equations to Saturation curve to determine saturation point. 4. Use IFLOWS for expansion past the saturation point and calculate Ag from NUCLEi. 5. Evaluate epsilon equations (3. 5 or 3. 8) to determine onset of condensation~

-396. Downstream of the onset of condensation, a joint solution of NUCLEI, NUCLE2, and CFLOW is obtained at each increment, Ax, and the parameters x, N, Ax, AA, A, P, p, T, U, M, r*, r Ar, Ag and g are stored for later print out. The saturated equilibrium expansion is compiled in a very simple program and need not be presented at this time. The Nozzle and Vapor Parameter Subroutines are included in a Program Common such that any of the included parameters, when used in the Main Program, are automatically evaluated and returned. Conditional statements are inserted in the program to account for the possibility of condensation occurring in either the converging or diverging portions of the nozzle, the possibility of reaching the end of the nozzle without any condensation having occurred, and the failure of the iteration of NUCLE1 and NUCLE2 with CFLOW to converge. The development of the equations is such that any consistent set of units can be used for the proplem. However, it was found that minimum scaling was necessary if the problem was worked in a cgs system of units. The problem, as written, is compatible with the cgs system of units and a list of the proper input data is presented in Appendix C.

VII. CALCULATIONS AND RESULTS The calculations performed in this study are divided into four parts: 1. A test of the program by hand calculations and a test of the method by comparing the results for a nitrogen expansion with the experimental results of Willmarth (3). 2. Justification of the assumptions of small radii, no drag, and accommodation coefficient approximately equal to one in the basic assumptions. Verification of the prediction procedure for the onset of condensation. 3. A variation of latent heat, specific heat, surface tension, and nozzle geometry for the nitrogen expansion to observe the induced effect of each parameter. 4. Application of the theory to metal vapors by computing the expansion from an initial pressure of approximately one atmosphere of copper, tin, lead, and zinc vapors. A. Validation of the Program The input data for the calculations performed in this study are tabulated in Table C-1 of Appendix C. A reference to Nitrogen (1), Zinc (1), etc., indicates the set of input data from which the indicated results were obtained. Any additional modifications of the input data, such as a change in a, L, or nozzle geometry, are listed directly on the plot of results. An explanation of the input data required, the proper units for each parameter, and the references (20-24) from which the vapor properties were obtained is also presented in Appendix C, -40

-41Initially the computer program was divided into small segments and hand calculations were performed for Nitrogen (I) to verify the computer results. With the assurance that each small segment was working properly, the entire program was assembled and hand calculations were performed to spot check the calculations of the saturation point, the point of onset of condensation, and the first few incremental steps of the condensing flow regime. With this degree of confidence in the program, the next step was to check the theoretical results with known experimental results. The Nitrogen (l) input data is designed to approximate an experimental test of the condensation of nitrogen performed by Willmarth (3). Nitrogen properties from NBS Circular 564 (24) and the International Critical Tables (23) required extrapolation below 70~K for the surface tension and the latent heat of vaporization. These values are generally valid only to within ~10% and, as will be shown later, this becomes extremely critical in the case of the surface tension. The vapor pressure curve for nitrogen is reasonably well known, and it was approximated by the Clausius-Clapeyron equation (see Fig. C-1). For comparison with experiment the calculations for nitrogen were carried out in full from the saturation point onward and no attempt was made to predict the point of onset of condensation. A comparison of the experimental results of Willmarth (3) to the theoretical prediction of the condensation of Nitrogen (I) is presented in Figure 1, The theoretical curve for Nitrogen (I) shows a degree of supersaturation of approximately 3. 5 degrees (20%o) lower than is indicated by experimental results. This is expected in that the amount of impurities in the experimeintal nitrogen was unknown and this advanced nucleation would cause an early break, However, of greater importance is the degree of

-42uncertainty in the extrapolation of the available surface tension data. In anticipation of the effect of a on the theoretical results, a second extrapolation was chosen for a. The a equation for Nitrogen (II) remains within the ~ 10% accuracy of the original surface tension data and, in actuality, gives a value for a which is only 7% less than the a from Nitrogen (I) at the onset of condensation. The theoretical prediction for the condensation of Nitrogen (II) compares most favorably with the experimental results in Figure 1. The strong influence shown by only a 7% change in a will be discussed more fully at a later point in the discussion. Figure 2 presents a more complete comparison of the theoretical results for Nitrogen (II) as compared to the experimental results of Willmarth (3). The isentropic and saturated equilibrium expansions are also presented in this figure. The theoretically calculated results presented in Figure 2 employ the full calculation of the nucleation and flow equations for the entire length of the nozzle. This is significant in that the flow properties prior to the "condensation shock" lie on the isentrope and the flow properties after the "shock" tend to approach the saturated equilibrium expansion. This is the expected result, and it gives further confidence in the program. The experimental points of Willmarth (3) fall on the isentrope as well as along the "condensation shock", which indicates a proper estimate of the physical expansion and validates the program. The close correlation between theoretical and experimental results within the "shock" itself lends further support to the basic nucleation theory of the condensation process.

-43Willmarth (3) obtained the pressure results, p, po and po' by direct measurement. The effective area ratio, (A/A*)eff, temperature, T, and mass fraction of condensate, g, were then calculated from the measured quantities by a joint solution of the flow equations. B. Justification of Assumptions 1. Drop Size A typical distribution of maximum and minimum droplet radius is presented in Figure 3 for the Nitrogen (I) expansion. Note the very short span over which the nucleation of critical sized drops, J curve, takes place. The average critical sized drop has a radius of 4. 5 x 10 cm, and the largest particles -7 which leave the nozzle have a radius of 8 x 10 cm. The radius of the nitrogen molecule from Loeb (25) is approximately 1i. 5 x 10-8 cm and calculated mean free paths vary from the order of 10 cm at the onset of condensation to 10-3 cm near the nozzle exit. According to Stever (5) the drop radius is larger than the minimum required to contain the 10-12 molecules necessary for the application of macroscopic theory. Therefore it is not unreasonable to assume that small drops are the general rule, that the no drag approximation is valid, and that since r << mean free path the samll drop growth equation can be used. 2. Accommodation Coefficient The effect of a variation in a is presented in Figure 4. Small variations near the value of a = 1. have very little effect. However, since the growth equation is a direct function of a, the curve for a- 0. 1 differs markedly from a = 1. O. Fortunately, the known values for a are from 0. 8 to 1, 0 so that the choice of a = 1. 0 is a

-44reasonable one. The computer program is so constructed that accurate values of a, when known, can be inserted as input data. An additional feature of Figure 4 is the indication that growth has little effect up to the "break" point but thereafter controls the remainder of the condensation process. This is further brought out by the plot of J in Figure 3, which indicates the rapid nucleation of a large number of particles to start the "shock" and then a continuing growth of these particles until saturated equilibrium conditions are met. 3. Point of Onset of Condensation From the previous discussion and from Figure 3, it appears that little or nothing happens in regard to condensation until a degree of supersaturation is reached for which the critical drop radius is very small. As r* from equation (4. 5) becomes small 2 the exponential term, which contains r,9 ceases to be the predominant term in equation (4. 8) and J increases rapidly. + 2411a (4. 5) pL RT in p/p. 47T r*2 J = T PL NA e 3kT (4. 8) The epsilon of Section III should then predict the approximate point of onset of condensation, Several calculations for various values of E are presented in Figure 5. e = 0 indicates that the full calculation of the nucleation and flow equations was carried out from the saturation point onward. The choice of epsilon results

-45in only a minor error when compared to the possible 10% error in surface tension. Epsilon values of 10 3 and 10 give almost identical results, thus E =. 001 appears to be a suitable choice for the expansion of nitrogen. Also in agreement are HeadVs (2) calculations for water vapor in which he indicates that the attainment of J = 1016 droplets per centimeter cubed per second signifies the point at which the condensation shock will break away from the isentrope. Several 16 calculations on the computer indicate that- J 10 at the breakaway point for nitrogen also. At present, however, this author sees no justification for predicting that this will be true for other vapors as well. C. Variation of Parameters 1. Surface Tension As indicated previously, the variation of a is the most critical aspect of the theoretical prediction of the condensation phenomenon. From equations (4. 5) and (4. 8) it can be seen that a enters the exponential term of J as a third power. Thus for a 10% change in a the nucleation rate changes by several orders of magnitude. Figure 6 and Figure 7 present several variations in a to show the large effect that this parameter has on the degree of supersaturation, Later results will show that the effect of other parameters is minor when compred to surface tension. There appears here a possible method of improving the surface tension data. The results of a closely controlled experimental vapor expansion could be used as a reference plot and surface tension varied on the computer until a matching plot is obtained. Providing

-46the other parameters are reasonably well known, this should provide a better estimate of a than merely extrapolating a curve far away from the few experimentally determined values that are available. 2. Specific Heat Cp is reasonably well known and only varies a few percent over the range considered in this analysis. For this analysis Cp is chosen as a constant, and from the plot in Figure 8 the constant assumption for Cp appears reasonable. The effect of a 10% variation in Cp is negligible when compared to an equivalent uncertainty in a. 3. Latent Heat of Vaporization Latent heat is generally a weak function of temperature. For Nitrogen (T) the variation in L is only 9% over the range of temperature from saturation to the nozzle exit. From the plot of variations in L presented in Figure 9 it is evident that the constant latent heat assumption in the basic equations is a reasonable assumption and any variation in L is negligible when compared to the effect of an equivalent variation in,. 4. Nozzle Geometry The condensation process is a function of time, so it is expected that the rate of expansion should have some effect on the "condensation shock. " Theoretical calculations were performed for the Nitrogen (I) expansion in a two dimensional nozzle for which the exit half-angle was varied from 9. 75 to 45 degrees. Nozzle lengths were selected so that a Mach 8 nozzle was used for each run. The results of varying the rate of expansion for Nitrogen (I) are presented in Figure 10 and Figure 11, In Figure 10 it is noted

-47that the faster expansion causes the vapor to supersaturate to a lower temperature before the onset of condensation. The "condensation shock" is thinner for the fast expansion, it is (in physical dimensions) 3 cm for the 0 1/2 450 expansion versus 5 cm for 7501/2 the 0 1/2 = 9. 75 expansion. 0 1/2 is the nozzle exit half-angle. This is somewhat misleading because, for Figure 11, in which all nozzles are referenced to a non-dimensional length for a Mach 8 nozzle, the faster expansion tends to smear out over a broader relative length. The plot of the percent of condensed vapor in Figure 11 indicates that the higher the rate of expansion the less moisture is formed in the nozzle. It appears that for sufficiently high rates of expansion a condensation free flow can be obtained. This has been borne out by "condensation free" experiments in the past. D. Application to Metal Vapors As stated in Section V, the theory as developed herein should apply to metal vapors. The primary interest here is in the finding of a metal vapor that will condense in a highly expanded flow and that is suitable for laboratory experimental use. Hill (11) has predicted that sodium and potassium will condense; however, these vapors introduce critical handling problems. Copper, zinc, tin and lead appear to have possible experimental use, with copper being the preferred vapor of the four. Theoretical calculations were performed for a two-dimensional nozzle with a constant slope diffuser for the metal vapor input data listed in Appendix C.

-481. Copper Fortunately, copper vapor is theoretically predicted to condense for a fairly rapid expansion in a two-dimensional nozzle with a 20 degree diffuser half-angle. The results are presented in Figure 12 for an expansion in a Mach 10 nozzle from an initial pressure of approximately one atmosphere and an initial temperature of 3000~K. The "condensation shock" occurs very near the throat and the change is very rapid with the'"shock" spanning less than a centimeter of nozzle length. The average critical drop radius is approximately 3 x 10 cm and the maximum droplet formation rate is of the order of 1019 drops per centimeter cubed per second. Again the onset of the condensation appears to be marked by a value of J 10 16 or 1017 drops per centimeter cubed per second. Copper could be a suitable choice for experimental work in rapid expansions because the large fraction of condensate formed in this theoretical solution indicates that higher rates of expansion and higher initial temperatures and pressures probably could be used without resulting in a "'condensation free" expansion. 20 Zinc Two computer runs were made for zinc, and neither run predicted any detectable amount of condensate. Zinc (I) was computed for a 45 degree diffuser half-angle and initial pressure and temperature of approximately one atmosphere and 1500 K respectively. The results of this calculation are presented in Figure 13 where a sufficient number of the calculated points are plotted along the isentrope to show that the expansion is essentially isentropic. This indicates a negligible amount of condensation. The tabulated calculations from the computer give a maximum g of 1. 48 x 10 at the nozzle exit.

-49The maximum nucleation rate obtained is J = 2 x 10 drops per centimeter cubed per second and the critical drop radius is -8 2 x 10 cm. It was anticipated that a slower expansion with a lower starting temperature would predict a detectable amount of condensate. Zinc (II) was run on the computer for a 20 degree diffuser and a pressure and temperature of one atmosphere and 13000~K respectively. The results in Figure 13 again plot along the isentrope and indicate no condensate. The maximum g of 4. 8 x 10 5 for the Zinc (II) run indicates some movement toward condensation. The maximum nucleation rate of J = 3 x 1014 drops per cubic centimeter per second is obtained at a vapor temperature of 375~K and the critical drop radius is again 2 x 10-8 cm. The J - 10 value with no condensate formed lends support to the 16 J 101 criterion for the onset of condensation. A slower expansion of zinc probably would produce a condensation shock- however, this would be out of the desired rapid expansion regime. From the results of the Zinc (I) and Zinc (II) runs it appears that for rapidly expanding flows zinc vapor may not condense. 3, Tin and Lead Unfortunately, the computer runs for tin and lead failed to converge in the CFLOWS Subroutine. The reason for this is not immediately apparent, and a further investigation will have to be made before any definite conclusions can be drawn. However, the few points that were calculated show promise that both tin and leaOd will condense in a rapid expansion.

VIII. CONCLUSIONS An adaptation of the digital computer to the analysis of the condensation of a flowing vapor is possible, and the program for a pure vapor is presented herein. Theoretical predictions of the condensation of nitrogen, copper and zinc indicate the following: 1. The nucleation theory of Frenkel (13) provides a reasonable approximation to the rate of formation of condensate in a flowing vapor. 2. The theoretical prediction for the condensation of nitrogen compares surprisingly well with the experimental results of Willmarth (3). 3. A small drop size is the general rule. The critical drop radius -8 is of the order of 10 cm and droplets of maximum growth rarely exceed 10 6 cm. 4. The nucleation rate controls the onset of condensation, and a prediction procedure based on this parameter might prove more useful than the one proposed herein. 5. The growth equation plays a predominant role once condensation has started and a more accurate expression would be beneficial. 6. The surface tension is the most influential parameter. Generally the value of a is uncertain and this leads to the largest source of error in the program. 7, Values of specific heat and latent heat of vaporization are not critical, and treating these parameters as constants is a reasonable assumption. 8. The condensation process is very sensitive to rates of expansion, and for rapidly expanding vapors a "condensation free" expansion is possible. -50

-519. Copper vapor will condense in rapidly expanding flows and appears suitable for experimental work in the condensation of metal vapors. 10. Zinc vapor does not condense for high rates of expansion from one atmosphere pressure, and runs for lower pressures and temperatures are recommended.

RE FERENCES 1. Stodola, A. Steam and Gas Turbines, McGraw-Hill Book Company, 1927, New York 2. Head R. M, Investiations of Spontaneous Condensation Phenomena, PhD Thesis, California Institute of Technology, Pasadena, California, 1949. 3. Willmarth, W. W., and Nagamatsu, H. T., "Condensation of Nitrogen in a Hypersonic Nozzle," GALCIT Hypersonic Wind Tunnel Memorandum No _6, January, 1952. 4. Nagamatsu, H. T., "Summary of Recent GALCIT Hypersonic Experimental Investigations, " GALCIT Report, January, 1954. 5. Stever, H. G., and Rathbun, K. C., "Theoretical and Experimental Investigations of Condensation of Air in Hypersonic Wind Tunnels, " NACA Tech Note 2559, 1951. 6. Wegener, P., Reed, S., Jr., Stollenwek, E., and Lundquist, G., "Air Condensation in Hypersonic Flow, " Journal of Applied Physics, Vol. 22, p. 1077, 1951. 7. McLellan, C. H., and Williams, T. W., "Liquefaction of Air in the Langley 11-inch Hypersonic Tunnel, " NACA Tech Note 3302, 1954. 8. Bogdonoff, S. M., and Lees, L., "Study of the Condensation of the Components of Air in Supersonic Wind Tunnels, " Princeton Univ. Aeronaut. Eng. Lab. Rept. 146, 1959. 9. Courtney, Q. G., "Kinetics of Condensation from the Vapor Phase," Texaco Experiment Incorporated TM 1250, 15 July 1961, Texaco Experiment Incorporated, Richmond, Virginia. 10. Daum, F. L., "The Condensation of Air in a Hypersonic Wind Tunnel,"IAS Paper No. 63-53, Institute of the Aerospace Sciences, New York, 1963. 11. Hill, P. G., Witting, H,, and Demetri, E. P., "Condensation of Metal Vapors During Rapid Expansion,'" ASME Paper No. 62 - WA - 123, 1962. -52

-53REFERENCES (continued) 12. Stever, H. G, "Condensation in High Speed Flows," High Speed Aerodynamics and Jet Propulsion, Vol. III, Princeton Univ. Press, Princeton, N.J., 1958, pp 526-573. 13. Frenkel, J., Kinetic Theory of Liquids, Clarendon Press, Oxford, 19569 pp 366-400. 14. Wegener, P. P., and Mack, L. M., "Condensation in Supersonic and Hypersonic Wind Tunnels," Advances in Applied Mechanics, Vol. V, Academic Press, 1958, pp 307-447. 15. Patterson, G. N., Molecular Flow of Gases, John Wiley and Sons, New York, 1956, pp 166-170. 16. Kennard, E. H., Kinetic Theory of Gases, McGraw-Hill Book Company, 1938, pp 311-315. 17. Stalder, J. R., Goodwin, G., and Greager, M. 0., "Heat Transfer to Bodies in a High-Speed Rarefied-Gas Stream," NACA Report 1093, 1952. 18. Buhler, R. D., Condensation of Air Components in Hypersonic Wind Tunnels. Theoretical Calculations and Comparison with Experiment, PhD Thesis, California Institute of Technology, Pasadena, California, 1952. 19. Shapiro, A. M., The Dynamics and Thermodynamics of Compressible Fluid Flow, The Ronald Press Company, New York, 1953. 20. Elliott, J. F., and Gleiser, M., Thermochemistry for Steelmaking, Addison-Wesley, Reading, Mass., 1960. 21. Weatherford, W. D., Tyler, J. C., and Ku, P. M., "Properties of Inorganic Working Fluids and Coolants for Space Applications, WADC TR 59-598, December, 1959. 22. Lyon, R. N., Liquid Metals Handbook, U. S. Government Printing Office, 19 54. 23. International Critical Tables, National Research Council of the U. S. A., McGraw-Hill Book Company, New York, 1933. 24. Tables of Thermal Properties of Gases, National Bureau of Standards Circular 564, U. S. Government Printing Office, 1955.

- 54REFERENCES (continued) 25. Loeb, Leonard B., Kinetic Theory of Gases, Dover Publications, 1961, Appendix 1.

-5560K 7 / 40 K Saturation Point M = 4. 61 (4. 60)* 20 K N~ 10 K!.! I I_ -. - 8K M = 5. 89 i(5. 84)- / f i J i-N.trogen (II) 6Experlment - M-KL 6. t t/ M 6/e ~~~~~~~~~~~.2_........................... 2 7sentrpe / NLtrogenl (I) - -Satur:tion Curve 2 K (a=2.8 2T - *Mach number in parenthesis indicates exper imental result 1 K- I LL______ ______ of Willmarth,(3S) 20 25 30 35 40 45 50 55 60 65 70 T, OK Figure 1. Theoretical Results for Nitrogen Compared to Experimental Results of Willmarth (3).

-5613, ~~~11 __. i~ Expermentai Results of 1 Willmarth (3) 58o - - At4.i4' — Saturat re Expansion 0. T eorE~tical Cah'.ulations Isentrope 50 4. 6 (4. Theretial C lculatio s 54 50 426 20 24 28 Distance From Throat, x, cm.e Figure 2. Theoreretic i-Toheoretihal Caluondensations 4 i ~iO..Cal ulations 9f Willmal.th (,3 0 4 8 12 16 20 Distance From Throat, x, cm. Figure 2. Theoretical Prediction of the Condensation of Nitrogen (II) Compared to Experimental Results.

-5719 -4 _10 -f.. __. __ __.... _L __., I - — t- I —i -....... 108 - _ j i,1 10o5 ____ I__ _ 1 I __ a 1 7 6 _, i 1 _ L A C _ _ _ _ _ / - t - _ _ _ _ _ | _ _~ O t a) t Xl 0 0 10'~C 1017 i, 2,: 6....20 24 I I I I I x cm - _ _ _ _ iso1 I i {___ I -1 4 8 12 16 20 241 X, cm Figure 3. Typical Expansion of Nitrogen (I) Presenting Variation of Nucleation Rate, Maximum Droplet Radius and Critical Droplet Radius Along Nozzle Axis.

-5850 48 - 46 44 at= i.0 42 /Y = -a=. 42 / -— = o. 9 o 40'"~-a =(. 8 0 40 38 /7 36 /34 32 -= i sentrope - 30110 9 8 7 6 5 _ 6-" 51 1 -j I // 0 -= -b -. -,a= I.9 -a= i. I 3 2 1.. 10 12 14 16 18 20 22 24 x, cm Figure 4. Theoretical Prediction of the Effect of a Variation in the Accommodation Coefficient on the Condensation of Nitrogen (I).

-5997 6 x,1~~~~~~~~~~~90 00 45............ ~0,or-K 50~.cO a E =.100 48 /E:.001 46 44 42 0 401 E;~ 38 ~~0001 36 0(10,. 34 32 Is ntrope 30 28 6 8 10 12 14 16 18 20 x, cm Figure 5. Effect of Choice of E as a Predictor of the Onset of Condensation for Nitrogen (T).

-6013 —— r_ 11 sentl ope9 07 04 ~.90 5c 5 1 I I 1 7- 1` 6~t''- Saturation Cui've 3 58 Sauration Point 54 S SLtura ion Curv 50 46 "e _ _ _'1 1~A ll11 42 / 90 o or a O 42 E? / / 38 34 30 Ise tro e 26 10 Nitr ogen (I) 8 E = 001 Bo 6'~90u/ 1 __ /90 aI /I 1 /' 2 / 2 0 0 4 8 12 16 20 24 28 x, cm Figure 6. Effect of + 10% Variation in Surface Tension on the Theroetical Prediction of the Condensation of Nitrogen.

-61)- l -18 60 l Em 121. t0o 23. 20 26. 50 Isentrope Saturation lurve I 10..90 25 35 45 55 65 75 T, OK Figure 7. Effect of Variations in Surface Tension on:the Predicted Supersaturation of Nitrogen.

-6252 S52stura tion Curv - 50 90ca C 48 Ol ~ 1.0 C 46 SL~~~L/P 44- _ __/ I/ /1. C 42 I ~ ~./ 0~/ 38 ~~~ / 36 IsentropI 34 1 e 32 30 -- Nitrogen (I) 7,,,, C 1. 107 N 1 07 / 6 __ 10 ~ 4 ~ 82 47 __ E=.i 5 4 3I A /. 5 c 2 /P: 1/ ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~00 06 8 10 12 14 16 18 20 x, cm Figure 8. Comparison of 10% mnd 50% Changes in Specific Heat to a 10% Change in Surface Tension for Nitrogen (I).

-6352 Satrati n 50 Cirve''~ T =56. 6 __ at 48/ 9.6 1.I i I I I /I I I 1 Y 46__ 44 M 42 0 1.~~~~~~~~~~~~~~~~1 5 40 38 36 34 Isentrope 32 30 8~~~~~~~~~~ Ni rogen (I ___ _'_ 9 L =2.3 x 0.90 = 01 5 - 4P 3 L~~~~~~~~~ ~~2~~~//_1. L 1- I I'1-P~~~~~~~L1 6 8 10 12 14 16 18 20 22 X, CmH1 Figure 9. Comparison of 10% anra 50% Changes in Latent Heat of Vaporization to a 10% Change in Surface Tension for Nitrogen.

-648 ~~~~~~~E1.O0l 2 04 3 _shr~ _,b 2 - 1 8~~~~~~~1125 NitrIsen ro 5~~~~~~~~~~~~~2 A* =..005 50 48 / 46 5~~~~~~~~~~~. = 9.-75 ___ __ / _0121/2 44 0 42 34 ~ Ise 32- - 50 ~~~~~~~1/2 0 2 4 6 8 0 12 14 Is6 x, cmI Figure 10. Theoretical Expansion of Nitrogen (I) for a Two-Dimensional Nozzle with Varying Exit Half-Angle, 8 1/2 46 61./2:4 40 38 1 0,~ [sent ropel 36/\% 34!' Isent~ope 32 0 2 4 6 8 10 12 14 16 x,q c m

-6560 56_ f Satur ttion 54 50 46 42 / // Nitrolken (I 30 0 0.2 0.4 0 ac6 0.h8 1.0ozz x/L 60 I~= 6 cn 0 0.2 0.4 0.6 0.8 1.0 x/L Figure 11. Theoretical Expansion of Nitrogen (1) in Mach = 8 Nozzle of Lengths 35, 16. 5, and 6 cm. 10 — 10~~~~~ / /~~~Oo 6 / /~~~~~~~~~~~~~~~~~~~~~~~~~~~~o 8 Nozzle of Lengths 35, 16. 5, and 6 cm.

-663000! I t3000 1 I Saturation 2800__ _ _ 2600.. 2400 2200 -. 1 __ 2000 -j 1800........t 1600_ _ 1400; Isen ropel 44.2 T1= 00D ~2 0 2op 4e01xitHaAge A =.645c 2m 2 0IseItrope T0 0 12 -2 0 2 4 6 8 10 12 x, cm Figure 12. Theoretical Prediction of Condensation for the Expansion of Copper Vapor in a Two-Dimensional Nozzle.

-67106 Zinc (I) 6 K _-_ P -= x10 d ynescm A r T = 15000K E /'IO = =450 / 10 de/e I 0O T =18000K 200 Fise 13 hertiaPeicfIse trope (I) E nn 10 Ise Zint (II) Vapor Indicating Copden Satdrione trope I) / Figure 13. Theoretical Prediction of the Expansion of Zinc (I) and Zinc (II) Vapor Indicating No Condensation.

Appendix A DERIVATION OF THE NUCLEATION RATE EQUATION This is in no way intended to be a complete derivation of the nucleation rate equation. The full derivation can be found in complete detail on pages 366 through 400 of Frenkel (13). A scant derivation is presented in Stever (12), pages 532 to 540. There appears to be an error of 7T in the final equation presented by Frenkel, so enough of the derivation will be reproduced here to clear up this discrepancy. The solution is a statistical problem of kinetic theory, In an equilibrium distribution of droplets where kv < bliq the following relation must hold. N S n N Sn - 1 N (A. 1) n n n n-1 n-1 n-1 where 0v and qliq are the thermodynamic potentials of the vapor and liquid drops, respectively. N and N are the number of droplets n n-i containing n and (n - 1) molecules, respectively, S and S. are n n-1 the surface areas of droplets containing n and (n - 1) molecules, respectively. an is the probability that one molecule leaves a unit surface area of a droplet containing n molecules. n 1 is the probability that one molecule condenses on a unit surface area. Assume an and n 1 to remain valid for the non-equilibrium case ($liq < ~v) and modify the distribution N to a non-equilibrium distribution f. The equation will no longer balance so (A. 1) become s -68

-69J -f S n -f S a X0 (A. 2) n n-1 n-1 n 1 n n n' where Jn represents the increase per unit time of the class (n) due to condensation of the vapor on the surface of droplets in class (n - 1) over the number, which, owing to evaporation, pass from the class (n) to the class (n - 1). The rate of change of the number of drops of a given class is determined from (A. 1) and (A. 2) to be f n = i J (A. 3) at n n-1 Frenkel employs the method of Zeldovich which utilizes the feature of the sharp maximum of A4b at n* to solve the equation. By considering only values of n > 10 say, quantities appearing in the equation will only vary slightly when n is changed by 1 and can be treated as functions of a continuous variable. Equation (A. 3) can then be integrated and the result is presented as equation (27b) on page 396 of Frenkel (13). J = N D(n xp _(A. 4) n )(ex kT) 21 kT (A. 4) where N - = number of molecules in the supersaturated vapor n kT 1'3 2/3 2/3 t2/3 D(n*) = S(n*) / A diffusion coefficient = (4) / 3/ V iq n* liq /3 - - P 4' 2 2 3 (=liq -v)/n* substituting these values into (A. 4) produces the equation J= 2r2 (A. 5) ()exp- 4 r)7

-70or, writing in the form of Frenkel's final equation (28a) n*J 2 - 4 7 ar *2 v) liq - O n*J 2r*2 ( - n* (A. 6) n 3kTex 3kT In comparing (A. 6) with Frenkelvs equation (28a), it is noted that Frenkel does have an extra factor of 7 in the coefficient as noted by Courtney (9) and others. The desired equation for this analysis is obtained by substituting a V liq the value of (Oliq - v) from the critical radius equation r* -liq1 4 3 and the value n*.. Vi r* into equation (A. 5) Vliq 3 J 2r*2 ( P ex4 3U 2 V liq Vliq3 l2*kT 3kT / mr* 3 2r*2 2T 2r*~'7 1rro 2c(T 4:!:r r*2 J(j )V liq exp.- (A. 7) kgT liq v 71 m UT To eliminate the necessity of calculating molecular values of V and m N the equation for computirng J in this liq lN' analysis is written J ) 1kT 41N r exp 3k) (A. 8)

Appendix BIBM 709-0 DIGITAL COMPUTER PROGRAM FOR CONDENSATION IN HIGHLY EXPANDED FLOWS Written in Michigan Algorithmic Decoder (MAD) Language This program is valid for a pure vapor with constant specific heat and constant latent heat of vaporization. The cgs system of units is utilized throughout and cgs values for the universal gas constant, Boltzmann's gas constant, and Avogadro's number have been included within the program. The symbolic names of parameters used in the program are listed and defined at the beginning of each applicable section of the program. The units for and input data required for utilization of this program are presented in Appendix C. -71

-72$COMPILE MAD,PRINT OBJECT,PUNCH OBJECT MAIN O01MAIN 000 R AERONAJTICAL ENGINEERING R NOZZLE CONDENSATION RESEARCH R A R --- MAIN PROGRAMR RRESEARCH — J. GRIFFIN PROGRAMMING — L. HARDING JULY 1963 PROGRAM COMMON GAMMA, MU, B, C, CP, L RHOL, SIGMA R STANi READ DATA RTHIS DATA SHOULD INCLUDE R R PZERO = INITIAL PRESSURE R RHZERO = INITIAL DENSITY...R-TZEO = INITIAL TEMPERATURE - _ - R PE = ISENTROPIC EXIT PRESSURE. THIS IS USED AS AN — R......- -. -i-NI-T-IAL A-PW ROX IMATION TO THE SATURA'TION PRESSURE. R R GAMMA = RATIO OF -THE SPECIFIC HEATS R MU = MOLECULAR WEIGHT OF THE VAPOR -— R-CF-C — -— = SPECIF-IC HEAT AT -— CONSTANT PRESSURE R L = LATENT HEAT (NOTE VAPOR PARAMETERS SUBROUTINE) -R —S-IGMA -= -SURFACE T'-SION (NOTE VAPOR PARAMETERS SUBROUTINE) R R C = SATURATION CURVE CONSTANT R B = SATURATION CURVE CONSTANT — R -A-LEP'H-A -.ACCOMODATION COEFF-ICIENT R -R — DEL-X.=...IN-C-fR-EMENTAL STETP IN X. THE NOZZLE CONDITIONS WILL R BE COMPUTED AT XCON+DELX,XCON+2*DELX,... WHERE R XCON IS THE CONDENSATION POINT. R XRANGE = LENGTH OF INTERVAL, STARTING AT THE CONDENSATION..-.. -R POINT, IN WH.ICH THE FLOW CONDITIONS WILL BE R COMPUTED. R..XPMTNT =-.THE- CO-NE-SAT iON POINT IS ASSIGNED A SUBSCRIPT OF R ZERO, ONLY THOSE POINTS WITH SUBSCRIPTS XPOINT, * 2*XPOINT,... WILL BE PRINTED. R TRANGE = PR-OGRAM WILL LOOK FOR A CONDENSATION POINT IN THE R INTERVAL (TSAT,TSAT-TRANGE), WHERE TSAT IS THE -W COMPUTED SATURATION TEMPERATURE. R DELTAT = INITIAL STEP SIZE USED IN LOOKING FOR CONDENSATION R POINT IN THE INTERVAL (TSAT,TSAT-TRANGE). R PERCNT = IF EPSLON(T) EXCEEDS PERCNT IT IS ASSUMED THAT R CONDENSATION HAS STARTED. R PRINT COMMENT $1 NOZ 1ZLE CONDENSATION RESEARCH $ PRINT COMMENT 0 INITIAL CONDITIONS$ PRINT RESULTS PZERORHZERO,TZERO,PE.... P-R-NT-CO-MMENT...........CONSTANTS DEPENDING ON THE VAPOR$ PRINT RESULTS MU,CP,GAMMA,SIGMA,L PRINT COMMENT - $0 SATURATION CURvE-C-ONSTANT-S PRINT RESULTS C,B PRINT "'-IT-" PROGRAM P-AWAETERS$ PRINT RESULTS DELX,XRANGEXPOINT,ALPHA -EX ECOUTE-NO-ZEZ-L PW-E ~ -THRO-TAT$-,'- AS -fAR ) -

-73FIND SATURATI N P INI BY' MATCHING ISEN'FROPIC PRESSURE'-OLUTION Tu ThE SATURATION CUkVE PRESSURE. Z =GANrIA/ (A C, M A-. ) Si,1 T=C/ (.'JT -ZLOG.(P)) PSA T=PZE7RO-( (T/TZiRO).P. Z ) TEST=.A\BS. (i.-PSAT/P)-.001 WHEN-VF? TECT.L. C., IR<AN-,TF - R'I10 2 P=P+ 5* ( PSAT-P) TRANSFER TO S1 S2 TSAT=T EXECUTE IFLOWS. (ASTAR,CP,GAMMA,PZERO,RHZERO,TZERO) EXECUTE IFLOW. (ASAT9,SATPAT,RHTOS'ATTSATUSAT) Z=MSAT EXECUTE NOZZLE.($INVERSS,Z,ASAT) WHENEVER Z(I).E.l.,TRANSFER TO 53 PRINT COMMENT $4 THE SATURATION POINT FOUND IS NOT INSIDE 1 THE NOZZLE, THE PROGRAM WILL CONTINUE.$ S3 XSAT=Z MDOT =ASAT*RHOSAT\ USA T PRINT COMMENT $2 COMPUTED SATURATION CONDITIONS$ PRINT RESULTS PSAT,RHOSATTSAT _ PRINT RESULTS MSAT,USAT PRINT RESULTS XSAT,ASAT PRINT RESULTS lDCT R li FIND CONDENSATION POINT BY MATCHING ISENTROPIC FLOW R SOLUTION TO THE NUCLEATION THEORY SOLUTION GIVING THE R AMOUNT OF CONDENSATE AS A FUNCTION OF AN ASSUMED FLOW,R DELT=DELTAT ITERF=O. _ TMIN=TSAT-TRANGE WHENEVER T,,l I N.L O.,Ti'I N=O..P-RINT CO~MMENT $2 PARAMETERS DETERMINING CONDENSATIO 1N POINT SEARCH$ PRINT RESULTS TSAT,T-fMIN,I)ELT PRINT RESULTS PERCNT R 6i4 T=TSAAT ITER=O. -- T=T-DEL rT WHENEVER T o.L TMIN.... -Z_"...' 56 ITERF=ITERF+1. DELT=. 5*DELT WHENEVER ITEkF.L.10.,TRAN'aFER [Ft S4 WHENEVER Z.E.O.oTRANSFER TO S7 PRINT COMMENT $0 TEN PASSES OVER TEMPERATURE RANGE SPECIFIED, 1 NO CONDENSATION POINT INDICATED.$ TRANSFFR TO START END OF CONDITIONAL EXECUTF IFLOW. ( A,M,P,RHO, T,U) EXECUTE NOZZLE. ( $ I NVEti', Z, A ) WHENEVER Z(1).E. O.

TRANSFER TO S6 S7 PRINT COMMENT $0 TEN PASSES OVER THE FULL RANGE OF THE NOZZLE 1, NO CONDENSATION POINT INDICATED.$. —-- TRANSFER TO START END OF CONDITIONAL EXECUTE NOZZLE.($AREA$,X-DELX,Z) R L=FLT. tT) J.. TDROP=C/ (R-ELOG. (P) ) RHOL=FRHOL. ( TDROP) SIGMA=FSIGMA. (TDROP ) EXECUTE NUCLE1. (A,DELX,DOT,P T,RAD, NDOT DELG') R -'R.COMPUTATION OF EPSILON ( T ) R EPSLON= (GAMMA- (1. /(M*M)) )/(GAMMA-1.) WHENEVER M.L.1.,EPSLON=1. EPSLON=A*(1.-(L/(CP*T))*EPSLON) EPSLON=EPSLON* ( DELG/DELG/DELA ) EPSLON=. ABS. EPSLON.-WHENEVER EPSLON.L.PE'ECNTTRANSFER TO S5 T=T+DELT DELT=.5*DELT I TER= I TER+ 1. WHENEVER ITER.L.LITER-T RANSFER TO 55 R R LITEi'R HALF-INTERVAL STEPS HAVE. BEEN PERFORMED AND R A CONDENSATION TEMPERATURE HAS BEEN FOUND. R PRINT COMMENT $0 HAVE FOUND TEMPERATURE THAT SATISFIES THE CO 1NDIT IONS FOR THE CONDENSAT ION POINT-$ PRINT RESULTS T,DELT TCON=T EXECUTE IFLOW.(ACONMCON,M PCONRHOCON,TCONUCON) Z MCON EXECUTE NOZZLE.($INVERS$,Z,ACON) -WHE-NEIE-ER Z ( 1 ). E. TI,TANSFER STO8.. PRINT COMMENT $4 THE CONDENSATION POINT FOUND IS NOT INSI 1DE THE NOZZLE, THE PROGRAM WILL CONTINUE.$ 58 XCON=Z L=FL.(T) TDROP=C/(B-ELOG. (PCON)) -HR-O-L-=-FRHOL H.TO P ( T.DROP ) SIGMA=FSIGMA. (TDROP) — EX E-XCU-TE NUCL-.-' -( (A-, CON-, DELX,MDOT,PCON,T TCON, RAD,NDOT,DELG ) G=DELG PR'NT COMM'ENT $4 COMPUTED CONDENSAT ION CONDI T IONS$ PRINT RESULTS PCON,RHOCON, TCON -P RI N-TT R- E-S-L T S- M-CON-",-UC ON PRINT RESULTS XCON,ACON PR I NT-RESULTS SIGMA,-fHOL, L,CP PRINT COMMENT $ NUCLEATION THEORY QUANTITIES$ 1P}i -N kRESU-LT-S RAD.~NDO) T,DELG,GD E L G, G RI N I T IAL I ZE FOR- NOZZLE COMPUTATIONS START I NG AT THE CONDENS

-75RATION POINT. THE VECTOR X CONTAINS THE POINTS ALONG THE RAXIS OF THE NOZZLE AT WHICH THE QUANTITIES DESIRED ARE TO RbE COt',,IPUTED. NOTE IHAT I 1' IS, NOT NECE"SARY THAT THE SPACING RbE UNIFORM SINCE A VECTOR DELX IS PROVIDED. AT EACH POINT ROF THE MESH WE WANT TO CALCULATE R P = VAPOR PRESSURE R RHO = MlIXTURE DENSITY R T = TEMPERATURE R' = FLOW VELOCITY R rM1 = MACH NUM1BER R DELG = CONDENsAfE INCREASE 0INCE IHE LAT MEsH POINT R G = FRACTION, BY MASS, OF THE VAPOR THAT IS IN THE R C~CONDENED oTATE L) LLR(I) = THE INC REASE IN RADIUS, SINCE THE LAST POINT, OF R THE DROPS INITIALLY FORMED AT POINT I R RAD(I) = RADIUS OF DROPS INITIALLY FORMED AT X(I) RNDOT(I) = THE NUMBER OF DROPS PER SECOND THAT ARE NEWLY R FORMED FROM X(I-1) TO X(I) R SIGMA = SURFACE TENSION R RHOL = LIQUID DENSITY R CP = SPECIFIC HEAT AT CONSTANT PRESSURE R L = LATENT HEAT RALSO THE IDENTIFYING QUANTITIES R X = POINT OF EVALUATION R DELX = DISTANCE FROM LAST MES'H POINT R A = NOZZLE AREA R DELA = AREA INCREASE FROM LAST MEsH POINT RONLY THE POINTS, THEIR INCREMENTS AND THE INTEGER INDEX N RARE KEPT IN -CORE, ALL THE' REST OF THE INFORMATION- IS MAINRTAINED ON TAPE. ALL QUANTITIES FOR POINT N, HOWEVER, ARE RREQUTRED TO OBTAIN' THOSE FOR POINT N+ —N, HENCE THEY' ARE RTEVPORARILY IN CORE AS UNINDEXED QUANTITIES WITH THE ABOVE P N A t,1LE 50 R R...I'NITIALIZE TO CONDENSATION POINT R T=TCON P=PCON RHO=RPHOCON U=UCON M=MCON X=XCON A=ACON EXECUTE NOZZLE. ($AREA$,XCON-DELX,A(1)) DELA=A-A(i1). N=O PPOINT=O EXECUTE IOCTRL, ]~-E L 0 =."O000 ]_*UON.-..... XLI MIT=XCON+XRANGE R R S29 U(1 )=U+DELU WHENEVER U(1).LO0. ~00 i~R-i -lq T C4-T-$z THE INITIAL -APPROX-IMATION TO THE VELO 1CITY AT THE NEXT MESH POINT IS NEGATIVE,$.............TRANsER TO S34 END OF CONDITIONAL PR

-76R BEGIN ITERATION ON DELG UNTIL NUCLEATION THEORY AND R FLOW EQUATIONS AGREE. A MAXIlIUM'OF GIl' iR iT AT IuNS ARE R ALLOWED PER ",iESH, POINT.'3-] — ~0 ITER=GITER N = N + 1 w'HEN-VEFR N.G.1000 0 0 PRINT COM',.MENT $4 NUM'BER OF fM ESH POI NTS LXCFE,L T1fi PRO. 1 G RAMED L I,"I IT.$ TRANSFER TO S34 END OF CONDITIONAL.... L)ELX (N) =L)FLX X(N)=X(N- )+DELX(N) WHENEVER XG.XLIMI T... PRINT CCO'i.'ENT $4 THE NEXT POINT WOCLD E XCE-ED T H E. Rit,' 1 DESIRED$. TRANSFER TO S34 END OF CONDITICNAL EXECUTE NOZZLE. ( $AREAX ( N ),A(1) ) WHENFVER A(1)F.F. PRINT COC'4ENT $4 NOZZLE SUROUTINE ERROR INDICATION - 1ME-SH POINT DOES NOT LIE INSIlE)F NOZZLE.$ PRINT RESULTS X(N) TRANSFER TO S34 END OF CONDITIONAL *,,, DELA(1)=A( 1)-A DELG ( 1 ) =DEL G( 1 )=G+DELG( 1 ) 531 EXECLTTE CFLOW. (A(1),DELA(1),AG(1),DELG(1) P,RHu,T,U,M) L=FL.(T(1) ) TDROP=C (B-ELOG. (P1)) ) RHOL=FRHOL. (TDROP) SIGMA=FS IGMA. (TDROP) TEST=G( 1 ) EXKCU UTE NUCLE2. ( X ( N ),DELX ( N ),A ( 1 ),DELA ( 1 ) P( 1),RHO( 1) T ( ) 1U( 1 ),fM ( 1 ),MDOT,ALPHA,N,RAD,DRAD,NDOT,DLG ( ) ) G(1) =G+DELG ( 1 ) TEST=.ABLS. (TEST-G(1)) WHENEVER TEST.L..00001, TRANSFER TO S32 ITLR= I T-R-1. WHENEVER ITERGE,O.,TRANSFER TO S31 PRINT COMMENT $.......... $ S3 2 -- X X( N ) DELX=DFL X( N ) A=A (1) DELA=DFLA( 1 ) P=P(1) RHO=RHO (i) T=T (1) DELU=U(1)-U U=U (1) M=M (1) DELG=DELG( 1) G=G ( 1 ) THROUGH S 93,FOR I=0,1,I.,I.(N-1) S33 RAD ( I ) =RAID ( I ) +RAD ( I ) EXECUTE IOC T R L.

-77TRANSFER TO S29 R 534 TRANSFER TO START R R R INTFRNAL FUNCTION FOR OUTPUT OF THE FLOW QUANTITIES$ AT R- MFSH POINT X. LATER THIS SECTI'CON MAY ALSO SAVE THIS R INFORMA4TION FOR FURTHER PROCESSING. R INTERNAL FUNCTION ENTRY TO IOCTRL. WHENEVER N.E.PPOINT PRINT COMMENT $4 NOZZLE CONDITIONS -- $ PRINT RESULTS N,X,DELX,A,DELA,P,RhO,T,U,M,RAD(N),NDOT(N), 1DELG,G WHENEVER N.G.O PRINT RESULTS RAD(O)...RAD(N-1) END OF CONDITIONAL''1PPOI NT=N+XPOINT END OF CONDITIONAL FUNCTION RETURN.R END OF'- FUN'CT ION......... R INTEGER I,NPPOINT DIMENSION DELX(i000OO)DRAD(1000),NDOT(1000),RAD(1000),X(1000) "-DIMENSION A( 1),DELA(1),DELG(1),G( i ),M( l,P(1),RHO( 1),T(1), 1 U(1),Z(1) VECTOR VALUES XPOINt = 1. VECTOR VALUES LITER = 10. V.ECTOR VALUES GITER = 10. R END O —,-b'' F P R O G'R M.......

-78-COMPILE MAD,PRINT OBJECT,PUNC-H OJECT IFLO IFLOOOIFLOWOO* R AERONAUTICAL ENGINEERING R NOZZLE CONDENSATION RESEARCH R R ISENTROPIC FLOW R RRESEARCH — J. GRIFFIN PROGRAMMING — L. HARDING JULY 1963 R R ISEN'TROP IC FLOW EQUATIONS FOR COMiPUT ING PRESSURE, MACH NUMBER RAREA, VELOCITY AND DENSITY GIVEN THE TEMPERATURE AND THE RFOLLOWING PARAMETERS R ASTAR - THROAT AREA OF NOZZLE — R..CP.-SP"ECIFIC HEAT AT CONSTANT PRESSURE R GAMMA - VAPOR EXPONENT R PZERO - IN I TIAL PRESSURE R RHZERO - INITIAL DENSITY - "TZERO - I N I T I AL TEMPERATURE R. —-E-X-T-E-RNAL F'UN'CT ION (.A1 A2,A3,A4,AA6 ) R R THIS ENTRY PROVIDES FOR PICKING UP THE PARAMETERS LISTED R ABOVE, THEY ARE THEN SAVED INTERNALLY FOR LATER USE. ENTRY TO IFLOWS. — A-S T A R =-A 1 CP=A2 GA M M A = A 3 PZERO=A4 P-'Z E R O = A. 4 -.RH7ZEROA5 - - -- TZERO=A6 -FUNCF-CTI N R E TUR -- R THIS ENTRY COMPUTES THE UN..KNOWN QUANT I T I ES R A] = NOZZLE AREA R -- ---- A2 -= MAC-H-'-NUMBER R A3 = PRESSURE -R -. -.- A-4- -T —- EN -SY —. ( -R HO- ) R A6 = VELOCITY R...IN..TERMtS" OF THE PARAMETERS SAVED- FROM THE LAST CALL. R OF'IFLOWS' AND -5 = EM PE1.ATU-R.. -A 5 = ---— T -M P -E-R A T-U R E R NOTE THAT THE LIST OF ARGUMENTS IS ALPHABETICAL. ENTRY TO IFLOW. T-A5 - R PRESSURE E —EXP=6GAMMA-/ ( GAMMA- 1. ) Z=(T/TZERO).P. EXP - -A-3 = —P ZE RO*-Z R VELOCITY AND MACH NUMBER.'.'Z.'';'' (Tf Z E RO-T T''.-... A6=SQRT. (CP*Z) -- — Z=Z'T-GAM MA-T'.. A2=SQRT. (Z/T) R DENS ITY AND NNO-ZZLE — AREA - - EXP=1./ ( GAMMA-1. ) z = I;* T A2' 7EXI -.______ _.. A4=RHZERO/(Z.P. EXP) E X PR' —-.-5'-( GMMA +.-1 —- EXP....

-79Z=2-.'*Z/( GAMMA+ 1. ) Z=(Z *P. EXP)/A2 A 1=ASTAR * Z R FUN —CT ION RETURN R END -OF FUNCT 0 N- -"

-80$COMPILE MAD,PRINT OBJECTPUNCH ObJECT CFLOW1 CFLOWOO* R AERCNAUTICAL ENGINERkINO R NOZZLE CONDENSATION RESEARCH R R CONDENSATION FLOW R RRESEARCH — J, GRIFFIN PROGRAMMING — L, HARDING JULY 1963 R RCOMPUTATION OF VALUES OF PRESSUREDENSITYTEMPERATUREMACH RNUMBER AND VFLOCITY THAT SATISFY THE FLOW EQUATIONS THAT RAPPLY AFTER THE ONSET OF CONDENSATION. IT IS ASSUMFD THAT WF RHAVE COMPUTLD THE CONDITIONS AT X AND DESIRE THOSE AT X+DELX RTHE ARGUMENTS REQUIRED ARE R DELA = AREA INCREASE FROM X TO X+DELX R A = AREA AT X+DELX R DELG = CONDENSATE MASS INCREASE FROM X TO X+DELX R G = CONDENSATE MASS AT X+DELX R P = PRESSURE AT X R RHOC = DENSITY AT X R T = TFMPERATURE AT X R U = VELOCITY AT X R M = MACH NUMBER AT X RTHE PARAMETERS CPL AND MU ARE REQUIRED FRC'OM PROGRAM COMM',ON. RTHE FUNCTION FL.(T) MUST BE USED HOWEVER SINCE ITERATION ON RT IS NECESSARY. THE UNKNOWN QUANTITIES ARE RETURNED IN P(1), RRHO(1),T(1) U(1) AND M(1) I. IS ASSUM1itED THAT U(1) CONTAINS RAN INITIAL APPROXIMATION TO THE VELOCITY AT THt NEXT MESH RPOINT ON ENTRY. R EXTERNAL FUNCTION (AgDELA,G,DELGPAtRHOATA1,UA,MA) R PROGRAM COMlA(MNA GAiMMA, MU, B, C, CP, L, RHOL9 SIG I A R ENTRY TO CFLOW. U=UA (1) DELU=U-UA MDQG=DELG/( 1.-G) DQA=DELA /A L=FL ( TA ) R.T -tHR'COUGH S2 FOR I=1 1IG1 100 DQU=DELU/U R COMPUTATION OF DELT AND T FROM THE ENERGY EQUATION R DELT=(L*DFLG-U*DELU) /CP T=TA+DFL T L=FL ( T) R COMPUTE DELRHO/RHO FROM CONTINUITY EQUATION DQRHO=DQA+DQU D-QRHRO=-"D Q RHO. R R COMPUTE DHEL-PJ — -FRMM TEQ-U-TI OrF STATE R " DQP -DZR HO+-(-DELT/'t )-MQ t-"''F ~I ~' Q - R R WE CAN NOW COMP-UTE THE VELOC'ITY CORRESPONDING -TO

-81R THESE DIFFERENCE QUOTIENTS FROM THE MOMENTUM EQUATION DELU=- ( DP*R*T* (1-G) ) / (MU*U) NEWU=UA+DELU TEST=.ABS. (1.-NEWU/U) U=NEWU 52 WHENEVER TEST.L.EPSLON, TRANSFER TO S3 R R ITERATION LIMIT HAS BEEN REACHED WITHOUT CONVERGENCE R PRINT COMMENT $'CFLOW' - ITERATION FAILED TO CONVERGE$ EXECUTE ERROR. R R CONVERGENCE CRITERIA SATISFIED R S3 TA(1)=T UA(1)=U DQRHO= 1.-DQRHO DQP=1.-DQP PA(1 )=PA/DQP RHOA( 1) =RHOA/DQRHO NEWT=R*GAMMA*T/MU MA(1 )=U/SQRT. (NEWT) FUNCTION RETURN R VECTOR VALUES R=8.314E+07 VECTOR VALUES EPSLON=1'E-06 INTEGER I R END OF FUNCTION

-82$COMPILE.qAD,PRINT OBJECT,PUNCH O03JECT NUCL1001NUCL100* R AERONAUTICAL ENGINEERING R NOZZLE CONDENSATION RESE-ARCH R R NUCLEATION THEORY - i R RRLSEARCH — J. GRIFFIN PROJGRAFMING — L. HARDING JULY 19o3 R RCOMPUTATION OF THE MASS DUE TO INITIAL CONDENSATION AS GIVEN RBjY NULCLEATION THEORY. TH TE FOLLO(WING OUANTITIES MviUST BE AVAILRABLE. R(1) CONSTANTS R NA = AVOGADROS NUMBER (M OILECULES/GRAM MOLE) R K = BCLTZO ANN CONSTANT (DYNE*CM/DEGREFEK) R P = UNIVERSAL GAS CONSTANT. (DYNE*-CN/G1/G RAN MOLE*DEGREEK) R(2) PARAMTERS RF B = SATURATION CURVE CONSTANT R C = SATURATION CURVE CONSTANT R DELX = DELTA XFOR VOLUME R M~',?DOT = RHOSAT*.-ASAT*'USAT R MU = MOLECULAR WEIGHT OF VAPOR R(3) VARIAE:LES K ~ A = NOZZLE AREA R P = PRESSURE R T = TEMPERATURE R RHOL = LIQUID DENSITY DEPENDENT ON TFMPERATURE R SIGMA = SURFACE TENSION DLPENDLNT ON TEMPERATURE R EXTERNAL FUNCTION (A1,A2,A3,A4,A5,A6,A7,Ab) R PROGRAv COMMON GAMMA, MU, B9 Cq CP, L, RHOL, bIGMA R RCOMPUTE THE RADIUS OF THE DROPS CONDENSING UNDER THESE RCONDITIONS, THEIR NUMBER AND HENCE THE CONDENSATE MASS RFORMED IN A VOLUME ELEMENT OF LENGTH DELX. THE PROGRAM RCOMMON PARAMETERS RHOL AND SIGMA ARE ASSUMED TO CORRESPOND RTO THE DROP TEMPERATURE FOR THE GIVEN FLOW CONDITIONS. R ENTRY iTo NUCLE1. A=A1 DELX=A2 MD T = A3 P=A4 T=A5 R R COMPUTE RADIUS CM* OF THE DROPS THAT ARE CONDENSING RADIUS= ( 2.*SIGMA*MU) / (RHOL*T) RADIUS=RADIUS/(R*(ELOG.(P)-B+C/T)) R R COMPUTE THE'E NUM-BER OF DROPS'- OF THIS SIZE THAT ARE R CONDENSING PER CENTIMETER CUBED PER SECOND. R TEXP=RADIUS*RADIUS/K..T E X —-T-, 1 88 7 9 0 2_- * S MA TE'xP")/'........... N=SORT.(SIGMA*,vU/NA(1 ))......N ='~-''NR",: (""[-i ~t'i.............................................................................................. N=P*N*EXP. (TEXP) /(K*T*RHOL) N=.7978846* *N

-83NDOT=N*A*DELX R R COMPUTATION OF PERCENT OF LIQUID MASS, DELTAG DELG=4. 1887902*RHOL*NDOT*RAD IUS DELG=bELG*RADI US*RAD I US/MDOT R A6=RADIUS A7 =NDOT A8=DELG FUNCTICON RETURN VECTOR VALUES NA=6.396*023,6.027E+03 VECTOR VALUES K=1.379E-16 1.379E-06 VECTOR VALUES R=8.314E+07 R END OF FUNCTION

-84$COMPILE MAD,PRINT OBJECT,PUNCH OBJECT NUCL2001NUCL200* R AERONAUTICAL ENGINEERING R NOZZLE CONDENSATION RESEARCH R R NUCLEATION THEORY - 2 R RRESEARCH — J. GRIFFIN PROGRAMMING — L. HARDING JULY 1963 R RCOMPUTATION OF THE QUANTITIES ARISING DUE TO CONDENSATION RUNDER THE FOLLOWING CONDITIONS R X = NOZZLE COORDINATE R DELX = INCREMENT SINCE LAST MESH POINT R A = NOZZLE AREA AT X R DELA = NOZZLE AREA INCREASE SINCE X-DELX R P = PRESSURE AT X R RHO = MIXTURE DENSITY AT X R T = TEMPERATURE AT X R U = VELOCITY OF FLOW AT X R M = MACH NUMBER AT X R MDOT = TOTAL MASS FLOW RATE R ALPHA = ACCOMODATION COEFFICIENT R N = MESH POINT NUMBER SINCE CONDENSATION. THE POINT R OF CONDENSATION IS MESH POINT ZERO. RTHE QUANTITIES TO BE COMPUTED ARE R RAD(N) = RADIUS OF NEW DROPS FORMED AT X. RAD(I) IS RADIUS R AT X-DELX OF DROPS INITIALLY FORMED AT MESH POINT R I AND IS NECESSARY INPUT ) RNDOT(N) = NUMBER OF DROPS OF RADIUS RAD(N) FORMED BETWEEN R X-DELX AND X. ( NDOT(I) IS THE NUMBER OF DROPS R INITIALLY FORMED BETWEEN MESH POINTS I-1 AND I AND 7~R IS NECESSARY INPUT. RDELR(O) = DELR(I) IS THE RADIUS INCREASE BETWEEN X-DELX AND X R THRU OF THOSE DROPS INITIALLY FORMED BETWEEN MESH RDELR(N-1) POINTS I-1 AND I, THESE QUANTITIES MUST ALL BE R COMPUTED. R DELG = CONDENSATE MASS INCREASE FROM X-DELX TO X DUE TO R NEWLY FORMED DROPS AND GROWTH OF OLD DROPS. R RTHE ARGUMENTS ARE THE QUANTITIES GIVEN ABOVE IN THAT ORDER. RNOTE THAT THE ZEROTH WORD OF THE VECTOR ARGUMENTS ARE TO BE RGIVEN, FURTHER THESE VECTORS CONTAIN REQUIRED INPUT AS WELL RAS SPACE TO PUT THE DESIRED OUTPUT. FURTHER THE USUAL RPROGRAM COMMON PARAMETERS ARE ASSUMED TO CORRESPOND TO THE RARGUMENT TEMPERATURE. R EXTERNAL FUNCTION (XDELX,A,DELA,P,RHOT,U.M,MDOT,ALPHA,N, T RAD,DELR,NDOT,DELG ) R PROGRAM COMMON GAMMA, MU, B. C. CP, L, RHOL, SIGMA INTEGER N R ENTRY TO NUCLE2. R EXECUTE NUCLE1. (A,DELX,MDOTP,TRAD(N),NDOT(N),DELG) R R WE HAVE NOW ONLY TO COMPUTE THE RADIAL INCREMENTS AND R THEIR CONTRIBUTION TO THE CONDENSATE MASS INCREASE. THE R DROP TEMPERATURE IS TAKEN AS THE SATURATION CURVE R TEMPERATURE CORRESPONDING TO THE GIVEN PRESSURE.

R TDRCP=C/(B-ELOG.(P)) DELR= 7978.846 * SQRT.( 1.379*.6027/(T*MU)) * ( TDROP-T ) DELR=DELX*(ALPHA/L)* (P/RHOL)*(DELR/U) THROUGH S1,FOR I=11,9.IG(N-1) S1 DELR(I )=DELR R R COMPUTATION OF RESULTANT TOTAL CONDENSATE MASS INCREASE R DELG2=O'.THROUGH S2,FOR I=O,1,I.G.(N-1) S2 DELG2=DELG2+ RAD(I)*RAD(I)*NDOT(I)*DELR(I) DELG=DELG+12.5663706*RHOL*D ELGG2/ MDOT R FUNCTION RETURN R INTEGER I R END OF FUNCTION

-863Cl P I LE 9'AC,PRINT INC JECT JO Z L E,P 1 N JECT NO ZLOOlNZLtO R AERONAUTICAL ENGINEkERIN R NOZZLE CONDENSATION RESEARCH R R NOZZLE R RRESEARCH — J. GRIFFIN PROGRAMM1,ING — Le. HARDING JULY 1963 R RORIGIN OF THE COORDINATE SYSTEM FOR THE NOZZLE IS THE CEINTER ROF THEF THROAT, I.E., A(O) = ASTAR. THE VARIALEL X INCREAS.S RPOSITIVELY INl THIE DIRECTION OF FLOWV' AND IS NEGATIVE UPSTREAM RFROM THE THROAT. THE NOZZLE PARAPMETERS ARE R ASTAR = THROAT AREA R X, MIN = X-COORDINATE OF INTAKE R XY AX = X-COORDINATE OF EXIT R I = INTAK." -E HALF ANGLE IN DEGREES R OUTA,.!iG = EX IT S I r)DE HALF ANGLE IN DEGREFS RTHF EfITRY THROAT RFAD S AND PRINTS THESF PARAMETERS. R XT NAL FIJN! T ION (A 1,A2,A3) ENTRY TO NOCZZL E - WHFNEVER Al.E. $THROAT$ R R RLTURN THROAT AREA AND INITIALIZE IF NECCESSARY READL DATA PRINT COMMENT $0 NOZZLE PARAMETERS$ PRINT RESULTS ASTAR, INANGOUTANG PRINT RESULTS XMINXMAX A= I NANG/RAD IAN INTAN=SIN.(A)/COS. (A) A='JU TANG/RAD I N OUTTAN=SIN. (A) /COS. (A) A3=ASTA. R FUNCTION RETURN R OR WHENEVER A1.E. $AREA$ R R CCOMP'JTE THE AREA AT A2 AND RETURN IN A3 WHEN:VER (XMIN.L.A2).AND. (A2.L.O.) A=INTAN* (,AR S. A2)'R WtHENEVFR (A2.G.O.).AND. (A2.L.XMAX) A= U T TAN - A 2 OTHERW I SE A3=. WHENEVER A2.F.O.,A3=ASTAR FUNCTION RETURN END OF CONDITIONAL A3=ASTAR+2. *A FUNCTION RETURN OR WHENEVER Al.E. $INVERS$ R COMIPUTE THAT POI NT ON THE X-AXIS OF THE NOZZLF THAT R HAS AREA'A3'. THIS POINT nFPFNDS ON THE MACH NUMBER R WHICH MUST f._ IN A2 ON' ENTRY, SEE THE PROGRAM FOR R THE PARTICULARS. IF A2(1) IS 1. ON RETURN THEN THE R S( OLUTION IS IN A2, IF IT IS ZLERO TH1EN THE SOLUTION IS

-87R STORED IN A2 BUT DOES NOT LIE INSIDE THE NOZZLE. A=.5*( A3-ASTAR) WHENEVER..'A.L.O. A=O. X(1)=l. MACH=A2 WHENEVER MACH.L.1. X=-A/INTAN WHENEVER XMIN.G.X,X(1)=O. OTHERWISE X=A/OUTTAN WHENEVER X.G.XMAXX(1)=O. END OF CONDITIONAL A2=X - A2(1)=X(1) FUNCTION RETURN R END OF CONDITIONAL R INTEGER Al DIMENSION X(1) VECTOR VALUES RADI'AN=57.2957795 R — E-ND OF FUNCT ION

-88$C0MPILE NAD,PRINT OBJECT,PUNCH OBJECT VAPOROO1VAPOROO* R AERONAUTICAL ENGINEERING R NOZZLE CONDENSATION RESEARCH R _ _ R VAPOR PARAMETERS R RRESEARCH — J. GRIFFIN PROGRAMMING — L. HARDING JULY 1963 R RCOMPUTE AS A FUNCTION OF TEMPERATURE THE FUNCTIONS R RHOL = LIQUID DENSITY R L = LATENT HEAT R SIGMA = SURFACE TENSION RNOTE THAT EACH OF THE QUANT I T I ES MAY BE COMPUTED SINGLY AND RTHAT THEIR VALUES ARE STORED IN LOW CORE. THIS PERMITS THEM RTO BE HELD CONSTANT WITHOUT INITIALIZING THESE SUBROUTINES. RANY ONE, OR ALL, OF THESE PARAMETERS MAY BE HELD CONSTANT B'Y RSIMPLY RETURNING THE VALUE FROM PROGRAM COMMON. RIT IS TO BE NOTED THAT SIGMA AND RHOL ARE FUNCTIONS OF THE RDROP TEMPERATURE, NOT THE VAPOR TEMPERATURE. IT IS THE DROP RTEMPERATURE THEN THAT MUST BE USED IN CALLING ON THESE TWO RENTRI'ES,',WHEREAS THE VAPOR TEMPERATURE IS USED WHEN COMPUTING RTHE LATENT HEAT. EXTERNAL FUNCTION (T) PROGRAM COMMON GAMMA, MU, B, C9 CP, L, RHOL, SIGMA R ENTRY TO FRHOL. R - COMPUTE THE VALUE OF THE LIQUID DENSITY FUNCTION RETURN VALUE R ENTRY TO FL. R CO.....CMPUTE THE VALUE OF THE LATENT HEAT FUNCTION RETURN VALUE R ENTRY TO FSIGMA. R - COMPUTE THE VALUE OF THE SURFACE TENSION FUNCTION RETURN VALUE R END OF FUNCTION

APPENDIX C INPUT DATA 1. Summary of input data and units required. 2. Table of input data used to obtain results plotted in Figures 1 through 13. 3. Plots of vapor properties obtained from References 20 through 24. -89

Summary of Input Data Requirements Parameter Definition Units PZERO Initial pressure (Po) dyne/cm RHZERO Initial density (pO) gm/cm TZERO Initial temperature (T) OK PE Initial approximation to saturation pres- dyne/cm2 sure (Pe) GAMMA Ratio of specific heats (y) N. D. MU Molecular weight of the vapor (Mi) gm/gmol CP Specific heat (Cp) dyne-cm/gm-~K L Latent heat (L) dyne-cm/gm SIGMA Surface tension (a) dyne/cm RHQL Liquid density (L) gm/cm B } Saturation curve constants {and OK ALPHA Accommodation coefficient (a) N. D. DELX Incremental step in X (Ax) cm XRANGE Length of interval, starting with conden- cm sation point, over which flow conditions will be computed XPOINT Points which will be printed out, i. e., N. D. every fifth, tenth, etc., point. TRANGE Temperature range over which program 0K will search for condensation point DELTAT Initial step size in condensation point 0K search over TRANGE PERCNT Value of epsilon for determination of onset of~ecimal condensation traction ASTAR Throat area (A*) cm XMIN X-coordinate of intake cm XMAX X-coordinate of exit cm INANG Inlet half angle degrees -90

Summary of Input Data Requirements (continued) OUTANG Diffuser half angle (0 1/2) degrees NA AvogadroVs number (in program) molecules/gmol 6. 027 x 1023 R Universal gas constant (in program) dyne-cm/gmol~K 8. 314 x 107 K Boltzmann's constant (in program) dyne-cm/~K 1. 379 x 10-16

-92Table C-1. Input Data N2(I) N2(II) Cu Zn(I) Zn(n) P 08.45x10 8.45 x 106 1 x 106 1 x 106 p0.00952.00961 2.54 x 10 4 3.93 x 10 4 6.04 x 10 4 0 T 298 295 3000 1500 1300 P 2000 2000 166. 7 166. 7 166. 7 e y 1.4 1.4 1.667 1.667 1. 667 A 28. 02 28. 02 63. 54 65. 37 65. 37 C 1.07 x 107 1.07 x 107 3.26x 106 3.18 x 106 3.18 x 106 P 9 9 10 10 10 L 2.3 x 10 2.3 x 10 4.94x1010 1.82x1010 1.82 x 10 a 25. 8-. 22T 24. 25-. 22T 1300 1000-. 26T 1000-. 26T PL 1.177-. 00476T 1.177-. 00476T 8.9-. 0007T 7.4-. 0008T 7.4-. 0008T C 882 882 34, 200 13, 800 13, 800 B 25. 7 25. 7 25. 46 25. 3 25. 3 a 1.0 1. 0 1.0 1.0 1.0 A*.0645.0645.645.645.645 X -5 -5 -5 -5 -5 min X 25 25 58 21.1 58 max INANG 45 45 45 45 45 OUTANG 9.75 9.75 20 45 20

2 P (dynes/cm cD 0~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ 0 0~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~0C Ii cD 02 "''0 0 0 H ill IIIIIIIIIIIIII...... I~~~d) I CD~~~~~~~~~~~o 0 I I II I il I l 1 II I 11 11 1 T F T 1,1 i 1 i i i IN z~~~~~~~~~ H HT o 0I I-'.~ ~ ~ ~ ~ ~~~~~IHl ~il IIII! 11111 NIIII I I 1 A1 (2~~~~~~~~~~~~~~~~~~~~~~1111 iliIiIIIHl c',1~ ~ ~ ~ ~ ~ ~~~~~~I11U11 11 ~ k11 I1:111MI11 E 111FVN I l 11 llii 1I,!ff - 01~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ fIIIIHl 01- 1t I I III II1 M 1-11u. 1 II T l Il 111l l, ~;1;V 02~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~'!

-9420 15 or = 25 8 - 22 T'k or ~= 24. 25 -.22 T C.) 10 Surface Tension for Nitrogen References (23) and (24). 0 30 40 50 60 70 80 90 100 110 120 1. 0K.9 0 Ill 1Iirlrirf~li~lllllllilllillli[Zlllillllll~ll~ T~]l~iiiiiiilitiitililiiiili~lili[ 0.~~~~~~~~0 570.046( 8 * 80 ~ V * 70 Figue C2. urfce Tnsin ad Lqui Density versus firoe Tempeatur fo Itrogencs 23 nd(2 i l~_~111111111ltllllllltllllllllllll1!III1I I I11. ~~~~~~~~~~ 00K illll IllllllllllllJ~Llllllllllllllllllitlllll i IllI ll l! Fir I-2 Sufc Tnio and LIqi Dest IversuIs J!I L.IJJ.

108 7 Ref. (20) Thermochemistry for Steelmaking 1 H (21) Properties of Inorganic Worki Fluids 106 (22) Liquid Metals Handbook. 105 - 10 10 103 C F I I I I I lausius-Clapeyron Curve M atch In P B - C/T 10 Metal C B Temp of Match __4,11f I I 11 I II vCu 34,200125.4613000 and 2000 10 IV I I I 1 HIPt~~t gU~~~ I | A|. El. | I ~Sn 32, 400125. 2030000 and 20000 _ _ IIIIIIII I I. I IIIIlIPb 121,240124. 0712000 and n1200 E1~z II.... I I1 1 1 1 111 1 1 1 1HI1 1Ilr Zn 13,800125.3011000 and 600~ IE I W I I I I I I VI I Y 11 I11,640 23 8411000~ and 600~ IIIII)III ItII IIIII I. II IIl K 9 900 23 4611000 and 500~ 10 1-W-6.4- I LJ 200 400 600 800 1000 2000 4000 6000 Temp, OK Figure C-3. Vapor Pressure Curves

-9611 _::: 1.0 7O 5~firt e i+Ht X = (2)Liquid Density H 1 1 ~I T I TI T t I II IZn, Sn Ref. (20) Thermochemistry for 1. 0 TTr I I I l I I I L r ww Steelmaking g 0+W0W00@+ 4 4 | i (22) Liq id Metals Handbook rs 8 111 1111!11111 i 11 11 IL AIlllH+ I~-~LI~T F~;btttFST~llf6~FfF~~STliltQ. 0. 6 0 500 1000 1500 2000 800 Kt'Ktt~i~tttttttff~Rit~t~HtttttttS~t~ f190 180 700 170 a ~~i- ~6~~~~~~~~~~~60 ah tillirr[!77 rT1 11 illllliTFI IIII IIII I Ilflill ll 400 0500 1000 1500 2000 Figure C-4. Metal Vapor Properties X I II | II I I I I I I | I I I I 1 I \ I I i I I I I I I i I I | I I I I I I I II I I I I X I — - - - - - - - -T T T fl T T T ItT S1 1 i 500~~~~~~~~~~~~~~+. 1IIIIIIIII b''llllllllllllll~l~lillljll -11111! Ifflkl~il~l~lwllll!Wlllll~lll~llllllll~!!l!l~l!!ll!!llllllllllllllllilSuri!1!!! Illlllllllllllll~lllllllllllllllllilllllldlilllll facelllll Tensllll Illlllllllllllllllllllll>LL I I 771 1 111111!1!1111!!!!!111!1111!11111111 LLL llllllllM L~lIIlll~lllllllilllllllllllllllllullllllllllllll I filll llll Na~~~

.4 o 0 0.)'-4 I III X 1: 0 - I fllll~ CIm I I Ln I mEl ml o I2 o ~ilT Bo Io: o o d 0 ~~ Ifl CV X0 2/T I d TOW 2A~ I IU I fill~~~~~~~~~~~~~~~~~~~ f I I~~~~~~~~9 TTI. f I~ ~ ~ ~ ~~~O~ Q,~~~~~~~~~~~tIll l i l t l i l t I l l I l l l l I r, ror~~~~~~~~~~~~~~ snaJ d3~ ~ ~~~~~~~~~~~L