condensation of water. Unresolved issues

Evaporation/condensation of water. Unresolved issues

Kazys Almenas

Evaporation/condensation of water. Unresolved issues I. Phase change at low pressures, laminar conditions

VYTAUTAS MAGNUS

UNIVERSITY V E R S US AU R EUS

UDK 556 Al-146

This monograph has been reviewed by: Prof. Dr. Bal Raj Sehgal, Emeritus Professor of Nuclear Power Safety, Swedish Royal Institute of Technology, Sweden Prof. Dr. Jurgis Vilemas, Vytautas Magnus University, Kaunas, Lithuania Dr. Tony Pollman, Chair for Research, Naval Postgraduate School, Monterey, CA, USA

Approved and recommended for printing by the Council of the Faculty of Natural Sciences of Vytautas Magnus University on 20 October 2014 (Protocol No. 04-01).

ISBN 978-609-467-098-5 (Print) ISBN 978-9955-34-531-2 (Print) ISBN 978-609-467-097-8 (Online) ISBN 978-9955-34-530-5 (Online)

© Kazys Almenas, 2014 © Vytautas Magnus University, 2015 © “Versus aureus” Publishers, 2015

1. Introduction This study considers very familiar, but still not completely understood phenomena – the evaporation and condensation of water. Mankind has used and studied these processes since antiquity; from a practical point of view these studies have been successful. Empirical and semi-empirical relationships have been developed which make it possible to predict the rate of condensation/evaporation processes for various applications over specified ranges of conditions. The properties of water and its phase changes have been analyzed by meteorologists, biologists, physicists and practically all categories of engineers. This has led to a fragmentation of the accumulated information, the development of discipline specific nomenclatures and empirical relationships which makes the transfer of information between disciplines more difficult. To summarize this brief overview – in most of its numerous applications the phase changes of water can be predicted and controlled successfully, and yet, a fundamental understanding of the mechanism of these processes is lacking. To be specific: The kinetic energy spectrum of evaporating and condensing molecules is not known. Very high condensation rates can lead to the disruption of the liquid-vapor interface and sudden decreases of pressure. The mechanisms which generate these “condensation shock” or “condensation implosion” events are not understood. Water is the most abundant substance on the surface of our planet. It is a deceptively simple molecule, just an atom of oxygen and two of hydrogen, yet its

5

properties are so complex that no completely satisfactory theoretical model of the water molecule is available. It is a polar molecule, and has an exceptionally high dielectric constant. Because it is so small and has a high H to O atom ratio, it forms strong cohesive bonds, has a high surface tension, and an exceptional capacity for dissolving ionic compounds. As a result, all of the liquid water that we encounter is in fact a more or less concentrated ionic solution, sea water being the most abundant example. Other prominent examples include the various fluids found in the bodies of all living organisms including our own. Tap water, as well as distilled water in glass or plastic bottles, all contain sufficient contaminants to influence some of its properties. That is one of the reasons why in spite of an abundance of experimental studies, relatively few experimental studies of truly pure water are available. The first part of this study assesses empirical data and theoretical approaches relevant for condensation/evaporation processes occurring at low vapor pressures. Information available in several disciplines is considered. Classical theoretical and experimental studies based on kinetic theory are augmented by information from the fields of thermodynamics, physics, chemistry and potentially relevant computational approaches developed in neutron transport theory and molecular dynamics. It is an interdisciplinary study, therefore definitions and nomenclature which are self-evident for specialists of a given field have to be defined and to some degree explained. The presentation of an inherently complex phenomenon can be overwhelmed and eventually confused by detail, this is especially true if one attempts to deal with all of the complexities simultaneously. The phase-change of water is a complex issue, in order to keep the review manageable, first the simplest theoretical model is considered, additional phenomena and corresponding analytical approaches are introduced sequentialy. For this purpose a sequence of four „control volumes “, identified as CV-1 to CV-4, is employed. The schematics of these control volumes depict one dimensional representations of the flow geometry and define the assumptions and boundary conditions used by the associated theoretical models. An advantage of such an approach is that it approximately parallels the historical development of the analytical models. The literature concerned with these issues is very voluminous, some reviews of condensation/evaporation processes encompass well over 100 positions. In this study only selected papers, illustrating key developments, are referenced.

6

2. CV-1 Schematic Fig 1 represents the “classical” 1-D control volume employed for evaluating liquid-vapor phase transfer on the basis of kinetic theory. Boundary conditions are the bulk temperatures and pressures in the liquid and vapor regions. The interface is a zero thickness surface at which an abrupt change from liquid to vapor density occurs. For the conditions considered in this study the bulk vapor velocity resulting from condensation or evaporation is much smaller then average molecular velocities, bulk flow is laminar, therefore energy transport in both phases occurs by thermal diffusion. This implies that for non-equilibrium conditions, a step temperature difference must exist at the interface. As shown in the schematic, the evaluation of net mass and energy fluxes across the interface requires the specification of Tfi, Pfi on the liquid side of the interface and Tgi, Pgi on the vapor side.

Pg Tg Evap.

Cond.

Tgi Pfgi Tfi, Pfi

Tf FIG 1. DV-1 Idealized phase-change control volume

7

The specification of these parameters is the major cause of uncertainty. In the liquid phase the intermolecular distance is on the order of 0.1 nano-m. A direct measurement of the diffusive temperature gradient in the immediate vicinity of the interface is thus not possible, and Tfi has to be assumed or inferred from other measurements. The pioneering work in this area harks back to the beginning of the previous century; Knudsen (1) showed that the evaporation rate of Mercury into a low pressure environment could be evaluated successfully by using kinetic theory and the assumption that the kinetic energy spectrum of the evaporating molecules has a Maxwellian distribution characterized by Tgi (2,3,4,5). Numerous investigators applied this conclusion to other volatile materials including water. It has remained the dominant starting point for the development of analytical approaches for a very wide range of phase change problems. The use of this distribution is questioned in this study, therefore the reason why it has been employed is presented in the words of previous authors. In an extensive study of condensation published in 1970 Shankar (4) writes: “We shall assume that the liquid emits vapor molecules in a Maxwellian distribution corresponding to the liquid temperature Tl” In a thesis evaluating alternative analytical approaches for non-equilibrium condensation evaporation problems Bond (8) in 2000 provides the following justification: “This is justified by assuming that at equilibrium (no net phase change conditions), the mass flux is zero, the temperature T across the surface is constant, and the vapor pressure is equal to the saturation pressure Psat(T). We then see that in equilibrium Pevap = Psat(T) This should be true for small perturbations from equilibrium, we can then replace Pevap with Psat(T) It is standard practice to assume the liquid phase is never far from equilibrium; this allows the molecules leaving the interface to be described by the Maxwellian. “ Applying these assumptions the absolute rate at which water molecules evaporate or condense is determined by integration over a half-range Maxwellian distribution at Tfi and Tgi. Net phase change rate is the difference of these two integrals. Knudsen assumed that all vapor molecules impacting the interface condense and all escaping liquid molecules remain in the vapor phase, the net phase change mass flux is then given by:

8

J 

M 2  Rg

 

Pgi

 Tgi



Pfi

  Tfi 

(1) In the literature this is usually referred to as the “Hertz-Knudsen” equation.

3. CV-2 Schematic Practically from the start it was recognized that when Eq. 1 is used to determine phase change rates of fluids other then monatomic Hg, the theoretical mass fluxes are considerably larger than experimentally measured values. This is especially true for water; in 1931, Alty (2) showed that for evaporation at low pressures the discrepancy can reach two orders of magnitude. Alty reasoned that Eq. 1 correctly predicts the number of molecules impacting or escaping from an interface, but not all of these molecules actually condense or evaporate. He proposed that a fraction of the vapor molecules impacting the interface are reflected. At equilibrium conditions the condensing and evaporating mass fluxes are equal; therefore, a similar mechanism must exist for the molecules moving across the interface from the liquid side. These assumptions were widely accepted, the correction factors required to force agreement between theory and experiment are usually referred to as “condensation and evaporation coefficients”, in some publications as “accommodation” coefficients. The CV-2 schematic is an example of the figures used to illustrate these processes (6, 8, 10, 11). It shows two distinct correction factors, the condensation coefficient αc and the evaporation coefficient αe. In a number of publications both coefficient are assumed to have the same value. Numerous experimental studies seeking to measure these correction factors were conducted; this included several approaches which explored more detailed reflection mechanisms. For example, momentum preserving and “specular”

9

collisions were considered, escaping molecules were assumed to be reflected back to the liquid with and without loss of kinetic energy. All of these studies share common features – the condensation/evaporation mechanisms are assumed to be symmetric and the correction factors are applied to the entire evaporation or condensation spectrum In 1953 Schrage (3) proposed a correction which accounts for the effect of bulk vapor velocity when net condensation or evaporation is present. The resulting widely used relationship for evaluating the phase-change mass flux is usually referred to as the “Hertz-Knudsen-Schrage” equation:

(2) Eq. 2 became the starting point for further development of analytical models and is still employed for evaluating condensation/evaporation mass fluxes. Given a circumscribed range of liquid/vapor conditions and accommodation coefficients which were determined for these conditions, this approach can yield useful results.

Pg Ty

(1-αc) Qc αcQc

αεQε (1-αε) Qε

Tf FIG 2. V-2 Condensation/Evaporation Rates with Corrections Factors 10

The problem is that though considerable effort was expended, no universal coefficients which would be valid for a wide range of conditions could be determined. In fact, as the number of studies increased, so did the range and variety of the adjustment coefficients published in the literature. By 2001 Marek & Straub (10) reviewed over 50 such studies and summed up their conclusions as follows: “The frequently applied assumption that evaporation and condensation coefficients are constant properties of a substance is not confirmed by the analysis carried out for water. Evaporation and condensation coefficients of water obtained theoretically or experimentally scatter in a range of more than two decades.” Two reasons are possible why the results of a computational model exhibit such a large variability: 1) The experimental measurements used to verify the model are subject to large uncertainties. 2) The computational model incorporates an assumption which does not correspond to physical reality. It will be shown, that when eq. 2 is used to evaluate condensation-evaporation mass fluxes for water, both of the listed reasons apply. The influence of measurement uncertainties is considered first.

4. CV-3 Schematic When net evaporation or condensation takes place, impurities present in the liquid or vapor tend to accumulate at the interface. This is illustrated in the CV-3 schematic; the left hand side shows that net evaporation carries non-volatile contaminants up to the interface. For non turbulent flow conditions, characteristic of low evaporation rates, contaminant concentration in the interface layer will increases with time, consequently the net evaporation rate should be time dependent. In fact, this has been observed in a number of experiments.

11

Net condensation

Net evaporation

Accum. of non-condensibles Accum. of impurities

FIG 3. CV-3 Impact of impurities on conden/evap rates

The right hand side represents net condensation of vapor in the presence of a non-condensable contaminant. The contaminant is carried to the interface, and accumulates there in a non-condensable rich layer. That is a well known phenomenon, in order to correct for it different empirical relationships have been developed for downward and upward facing condensing surfaces, as well as for non-condensables which are heavier or lighter then water vapor. The described phenomena are characteristic of all fluids; they have an especially large influence for water. As noted in the Introduction, water is a universal solvent of ionic compounds. In fact, it is such a good solvent that pure water does not exist in nature, furthermore, even purified water used in standard laboratory experiments contains sufficient impurities to influence evaporation rates. Simple distillation is not adequate. Marek and Stroub (10) recognized that the scatter of measured phase change rates can, to a significant part, be traced to this cause. In their conclusions they write: “Even small contamination of the surface significantly reduce the interfacial mass transfer. Evaporation and condensation in glass vessels can strongly be hindered by the accumulation of dissolved glass components on the interface.” 12

Steady state evaporation can amplify the concentration of impurities at the interface to a surprising degree. For example, if the liquid concentration of impurities is only 1 ppm, after evaporation of just a 1 cm column of water, their accumulation at the interface can rise up to ~15%. In effect the interface becomes a different fluid, with altered properties. A contributing factor is that Eq. 2 evaluates the net condensing or evaporating mass flux by subtracting two large numbers, this amplifies sensitivity to errors. Recognition of these experimental difficulties led to improved experimental procedures. An important aspect of these improvements was exceptionally careful purification of the water employed in the tests. Measured parameters of experiments using the carefully purified water showed some fundamental disagreements with accepted theory, the cause of these disagreements could not be ascribed to experimental uncertainty.

5. Theory vs. experiment Classical kinetic theory, on which Eq. 1-2 are based, requires that Tfi = Tgi at equilibrium, and that net evaporation can take place only when Tfi > Tgi. Though some earlier experimental studies indicated that this might not always be the case, they were not considered to be conclusive because of inadequate temperature measurement resolution. Only in the past decades a series of experiments, using specially designed test stands and extremely pure water, convincingly showed that net evaporation into a low pressure vapor can occur when Tfi < Tgi. Measurement results reported by Ward et. al. (9, 13) in Fig 4, and Badam et. al. (15) in Fig 5, show that steady state evaporation can proceed from liquid having a lower surface temperature into vapor at a higher temperature. This not only contradicted prevailing theory, it was counter intuitive. As the figures show, the “anomalous temperature jump”, as it was initially called, can be on the order of several degrees K.

13

FIG 4. Measurements of temp. profiles for evaporation and condensation from a surface (From Ward & Stanga Ref 8)

Both sets of experiments were conducted at reduced pressure, for the first experiments series (Fig 4) the interface is curved, a question thus remained whether the observed „temperature jump“ could be attributed to surface tension and the resulting pressure discontinuity. Subsequently the test stand geometry was redesigned so that temperature profiles could be measured along the centerline of a test channel (15). As shown in Fig 5, the results confirm the conclusion that net evaporation can proceed even when the vapor near the interface is at a significantly higher temperature then the evaporating liquid. Besides the exceptional purity of the water, an important aspect of these tests is the improved resolution of the measured temperature profiles. The TC bead used in the measurements presented in Fig 5 had a diameter of 25 μm, and could be positioned with an accuracy of 10 μm. The mean free path of the vapor molecules is ~44 μm at 245 Pa and ~18μm at 560 Pa. These two pressures span the range of the Badam tests. Measurement resolution near the interface is thus well within the average mean free path in the vapor region. 14

The publication of these results led to a renewed interest in the evaporation/ condensation issue. Several analytical approaches were proposed to deal with the “temperature gap” problem; they are customarily identified in the literature by three letter mnemonics. The most promising of these are: TIP (Thermodynamics of Irreversible Processes) (12) and SRT (Statistical Rate Theory) (12, 14) An overview of these methods is provided by Bond (7). The proposed methods are potentially promising in that they encompass both phases and broaden the range of thermodynamic variables (e.g. Cpf, Cpg and hfg are included). They analyze the thermodynamic equilibration of the phases in a small control volume around the interface; the direction of phase transfer is specified by the requirement of a positive change in system entropy. They can indeed show that under certain conditions net evaporation is possible even when Tfi < Tgi, however their practical use is limited. More detail regarding both methods is provided in Appendix A.

FIG 5. Measurements along the centerline of a rectangular test channel (From Badam etal Ref 14) 15

6. CV-4 Schematic The CV-4 schematic shown in Fig 6 represents a closer approximation of physical reality. In the earlier schematics the vapor and liquid regions are assumed to be homogenous and continuous; their interface is a zero thickness surface. Such representations have substantial analytical advantages, mass and energy balances and their transfer rates can be evaluated using closed form solutions of continuous functions. It is the traditional approach taught in textbooks and does not require extensive computational resources. When used appropriately these methods yield results which agree reasonably well with experimental data, for situations where significant deviations are present, empirical correction factors can be employed. The CV-4 schematic shows a simplified representation of the discrete molecular nature of water. The H2O molecules are assumed to be spheres which are free to move in the vapor phase, in the liquid they are held together by isotropic short range van der Waal type forces. The fact that the H2O molecules are dipoles

Knudsen layer Surface tension layer

FIG 6. CV-4 Molecules with isotropic attractions forces 16

and are able to form short lived attachments with other dipoles is not considered. This simplification makes it easier to use phase dependent macro properties, such as the latent energy of vaporization, heat capacities and the speed of sound to infer averaged molecular characteristics. An important difference from earlier representations is that in CV-4 the vapor and liquid regions are not homogeneous. The interface between the liquid and vapor is not a zero thickness surface which in mathematical models can be represented by a step function, but has a transverse dimension. This has important physical consequences some of which, like surface tension, are readily observed at the macro scale. The actual physical reality is more complex; however, the representation shown in CV-4 is adequate for a simple explanation of surface tension, in a number of basic textbooks on the subject it is presented as follows: In the bulk liquid each molecule is surrounded by 6 neighboring molecules, at the surface by only 5. A surface molecule thus experiences an unbalanced pull in the downward direction, the upper molecule layers are therefore more closely packed, an important consequence is that this generates a pressure component which is inversely proportional to the curvature of the surface. On the molecular level these forces can be quite large, in fact they are larger then the average translational kinetic energy of a vapor molecule impacting the interface. For example, at the conditions which generated Fig 5 the kinetic energy of the average vapor molecule is ~0.03 eV1, the average binding energy of a surface liquid molecule evaluated using surface tension is ~0.045 eV. It seems reasonable that the presence of such a barrier between the liquid and vapor regions will have an effect on evaporation/condensation processes, however in most of the studies referred to so far this has not been considered. Eq. 1 and 2 are based on the assumption that molecules escaping from the liquid region have a Maxwellian distribution characterized by Tfi, and that the impacting molecules have a distribution characterized by Tgi. The imposition of assumed phase change spectra leaves no room for the potentially disturbing effect of surface tension. 1  An eV (electron volt) is an energy unit used in neutronic calculations, it is equal to 1.602E-19 J

17

Analytical approaches based on kinetic theory transfer spatial non-homogeneity to the vapor region, this is shown in the CV-4 schematic as the “Knudsen” layer. The “Knudson” layer term is used in many technologies, it identifies a transition layer at the boundary between a solid and a fluid, or a fluid and a gas. In evaporation/condensation studies it is taken to be the transition region in which a Maxwellian distribution at Tfi merges with a Maxwellian distribution at Tgi. In some studies it is considered to have the thickness equal to a mean free path, in others an unspecified thickness of „several “mean free paths is quoted. The mathematical procedures analyzing what happens in the thus specified transition region can be quite complex, that is especially true if the condensation or evaporation generated bulk vapor motion is taken into account. Excellent examples of such studies are provided by Pao(5) in 1971 and Deng (16) in 2008.

7. Molecular characteristics of water at low pressures Kinetic theory provides the answer to the question of what happens when freely moving point particles interact until an equilibrium state is reached. The CV-4 schematic illustrates that the evaporation/condensation process raises additional questions, for example – how do the tightly bound liquid molecules manage to join the freely moving molecules in the vapor region? Precisely how that occurs remains to be answered, Table 1 lists some low pressure vapor and liquid properties which will help to put the characteristics of these interactions into context. Liquid properties are shown in the first line, analogous vapor properties in the second, the third line presents the liquid/vapor ratios. The vapor pressure is chosen so that it approximates the experimental pressures at which the temperature profiles shown in Fig 5 were measured. Table 1 illustrates that the liquid/vapor ratios range over many orders of magnitude. The principle cause of these variations is the very large molecular density difference.

18

Molecular density (Molec/mm^3) Vapor Liquid Ratio Liquid/ vapor

Time between Av distance Interactions between for an average molec molec. (in nano- m) (in nano-sec)

1.49 E14 3.34 E19

18.3 0.030

31. 0.013

2.24 E5

0.0017

2400

Interaction time Energy required of a surface to increase a molec with α 1 cm3 vol. in vapor or liquid vapor & liquid by molecules 1 K (in J) (in nano-sec) 450 1.12E-5 0.0026 4.2 180 000

380 000

Table 1 Comparison of conditions in the liquid and vapor adjacent to the interface for T = 273 K, P = 560 Pa

If the relative areas of the vapor and liquid regions shown in Fig 5 would be used to illustrate this density difference, less then 0.1% of a vapor molecule would appear in the CV-4 schematic. Such a large difference of the molecular densities implies that the vapor molecules have a small influence on the evaporation process. There are simply not enough of them. This conclusion is strengthened when the effect of other liquid/vapor properties are included. For example, in the liquid a perturbation is propagated at a speed of 1400 m/sec, in the vapor approximately 400 m/sec. The smaller distances between molecules and the increased speed implies that interaction time between liquid molecules is much shorter then between vapor molecules; in Table 1 it is listed as 0.013 nano-sec2. Column 4 shows the interaction time of a liquid molecule at the interface with a vapor molecule (in the first line), and with a neighboring liquid molecule (in the second line). The last line in column 4 shows that a liquid molecule at the surface will interact ~180 000 times with its neighbors for every time it interacts with an impacting vapor molecule. This strengthens the conclusion that for low vapor pressure conditions evaporation is determined almost entirely by the liquid phase. 2  Note that this estimate is obtained using the simple molecular model, quantum-mechanical approaches which considers the fluctuations of electrons, estimate an interaction time on the order of  0.001 nano-sec ~

19

Mass transfer is accompanied by energy transfer. The “energy carrying” capability of the two phases are very different, this is illustrated in column 5 of Table 1, which lists the amount of thermal energy that is required to increase the temperature of a unit volume of either phase by 1 K. The ratio of these values confirms that the energy required for evaporation must be transported predominantly through the liquid region. This conclusion is confirmed by the temperature measurements documented in ref 8 and 14 which show a negative temperature gradient in the liquid region. Fig 9 shows an example of this gradient measured during a “no heating”3 test (15). The temperature profile in the liquid shows that thermal energy is being conducted toward the interface from the liquid region, however, Fig 9 also shows that within the vapor region the gradient toward the interface is negative as well! Thus, thermal energy is diffusing towards the interface both from the liquid and vapor regions. In fact, for conditions when net evaporation occurs the interface becomes an “energy sink”! This observation is confirmed by the test results shown in Figs 4 and 5.

8. Agitational and total energy Considered in the classic kinetic theory context both the “temperature gap” and the “energy sink” phenomena are anomalies. Though it is universally recognized that the H2O molecule is not a “hard ball “as required by kinetic theory, it is widely assumed that deviations from this ideal can be compensated by using appropriate correction factors. In the majority of situations that is indeed the case, for gases and vapors the van der Waal constants can compensate small deviations from the PV = NRT relationship, and temperature differences can be used to evaluate energy transfer rates by employing heat capacities which include 3  A unique feature of these tests was that the vapor region could be heated by resistance wires placed just 3 mm above the interface. Fig 5 shows the measured temperature trace for a test with no-heating and three traces with increasing heating intensity.

20

Fig 7 Change of the temperature (K) and state of water with enthalp

an internal energy component. This works satisfactorily as long as the ratio of the internal and agitational (kinetic, in the case of gases) energy components remains almost constant. For processes involving phase change this is no longer true, for water the deviation from the ideal kinetic theory model is especially large. This is shown graphically in Fig 7 by plotting the change of phasic enthalpy as a function of temperature. The figure illustrates that the evaporation of water requires uniquely large inputs of thermal energy; the enthalpy difference between the liquid and vapor phases is significantly larger then between ice and liquid water. The uniqueness of water is illustrated further in fig 8 which depicts the agitational, internal and total heat capacities for all phases of water plus two selected materials. Ne is a monatomic gas (atomic weight 20), Mn a solid (atomic weight 55, ~ three times that of water). The ideal gas Ne has only agitational (translational) kinetic energy, its heat capacity can be evaluated using kinetic theory, the obtained result of 1039 J/kg K agrees very closely with measured values4. In comparison, the measured heat capacity of water vapor is ~1950 J/ kg K, and it includes besides the translational kinetic energy rotational and vibrational energy components. Ice and Mn are both solids; their heat capacities 4 Cp = (5/2) Rg/M

21

depend on the structure of their lattices and are not so easily compared. The Mn atom is almost 3 times heavier then an ice molecule, but since the multi-atomic H2O molecule has more degrees of freedom then the monatomic Mn, it has a ~ 4.4  times larger heat capacity. The vapor and solid phases of water thus have heat capacities which have a similar magnitude as the heat capacities of comparable materials, it is the liquid phase of water which displays unique characteristics. As shown in Fig 8, its heat capacity is 9 times as large as the heat capacity of Mn, more then twice as large as the heat capacities of ice or water vapor. The principal cause of this remarkably large heat capacity is a strong, fluctuating water-to-water molecule bond between an H atom and the O atom of a neighboring liquid molecule. Hydrogen and Oxygen atoms are present in other molecules besides water, consequently “hydrogen bonds” occur in those molecules as well, however, for the case of the H2O molecule the effect of the H-bond is unique. Water is the smallest molecule that can be formed between H and O atoms, and it has the highest H to O atom ratio, therefore multiple bonds can exist for each water molecule. The bond temporarily ties two molecules together, at a given temperature the formation and dissolution of these bonds in the 4.5 4 3.5

kJ/kg*K

3 2.5

Inter en Agit en

2 1.5 1 0.5 0

Ne – 20

H2O – ice

H2O – liqH 20 – vap

Fig 8. Agitational and internal molecular energies

22

Mn – 55

liquid phase is in equilibrium, therefore their fluctuation does not contribute to the agitational energy component of the liquid. However, when an unbalanced formation or dissolution of hydrogen bonds occurs, it can have a large effect. During evaporation H-bonds are dissolved, during condensation they are recreated, that is the cause of the exceptionally large latent heat energy of water shown in Fig 7. As the temperature (and thus the agitational energy component) increases, the rate at which the H-bonds are dissolved and reformed increases as well – that is the cause of the uniquely large heat capacity of liquid water shown in Fig 8.

Fig 9. “No heat” test from Badam et al (Ref 14) Vapor pressure 516 Pa

9. The hydrogen bond and its impact on the phase change of water The identification of the H-bond and its importance is not a new issue. The first mention of the unique role of H atoms in the formation of molecular bonds dates back to the beginning of the last century, its importance was recognized in the formation of organic compounds. In the past decades it has been studied

23

Fig 10. Schematic of the water molecule dipole and potential multi-molecular structures in the liquid phase

extensively in the field of biochemistry. The specialized scientific literature on the subject is very large. Text books provide the following generic definition: “A hydrogen bond is the attractive force between the hydrogen attached to an electronegative atom of one molecule and an electronegative atom of a different molecule. Usually the electronegative atom is oxygen, nitrogen, or fluorine, which has a partial negative charge. The hydrogen then has the partial positive charge.“ Some authors caution that it is not a true „bond“ since it fluctuates and in fact does so very rapidly. That is a semantic issue, there is no question that for molecules containing hydrogen atoms their tendency to interlink with neighboring molecules accords them a special importance. For water this importance is overwhelming, it is no exaggeration to state that the physical properties of water are determined predominantly by the H-bond. Fig 10 shows a computer graphics generated depiction of how the H-bond binds the H2O molecules in the liquid phase. It is colorful, and greatly simplified, but it does illustrates some basic aspects. Since the H2O molecule has two hydrogen atoms, it can form multiple hydrogen bonds. The O atom can be attracted to 2 neighboring H atoms, and the H atoms can form individual bonds. As noted, these bonds fluctuate rapidly, depending on temperature, on the average there are about 3 H-bonds per molecule.

24

The effect of the H-bond on the macro properties of water is illustrated in Figs 7 and 8; table 2 illustrates this importance on the molecular level. The table compares the strength of the H-bond first to a permanent covalent bond and the attractive van der Waal forces, second to the translational kinetic energy of vapor molecules. As seen, though the H-bond is weaker then the covalent bond holding the H2O molecule together, it is about 20 times stronger then the van der Waal attraction which determines the liquid phase for compounds lacking hydrogen. Its importance for evaporation/condensation processes can be assessed by comparing H-bond strength with the kinetic energy of vapor molecules. As seen, the internal energy stored in an H-bond, is about 10 times larger than k*Tg, the kinetic energy of a vapor molecule at the most probable speed of an equilibrium Maxwellian distribution. Even more telling, it is twice as large as the kinetic energy of the vapor molecule at the 1% cutoff of the Maxwellian distribution. In some condensation models it is suggested that impacting vapor molecules can „knock“ a liquid atom out of the surface, a comparison of the respective energies shows that this is improbable. Average bond strength (in eV) O-H covalent bond:

Relative strength (Compared to k*Tg)

5.1

217

Hydrogen bond:

0.240

10.2

Van der Waals forces:

0.013

0.55

Molecular agitational Energy k*Tg at Tg = 273K

0.023

1

0.13

5.5

Agitational energy limit for 1% of high energy molecules

Table 2. The hydrogen bond in water

The information provided in Tables 1 and 2 and Figs 7 to 10 show that computational models evaluating phase change processes of water must consider

25

the influence of the H-bond. A definitive resolution of this issue remains to be achieved; this study uses the assembled data to propose a probable molecular level evaporation mechanism. One of the basic conclusions on which the proposed model is based is that the evaporation is driven by the exchange of agitational energy between liquid molecules. The ratios of the liquid-liquid and liquid-vapor interaction probabilities shown in the last line of Table 1 are evaluated for low pressures, however, they are so large that even at significantly higher pressures interactions between liquid-liquid molecules will predominate. The second conclusion is that the breaking of an H-bond requires more energy then can be generated by a single liquid to liquid molecule interaction. That is evident from the bond strength values shown in Table 2. The table lists average bond values, the agitational energy component k*Tf is also an average, assuming it has a Gaussian distribution ~1% of the time this energy component can reach values which are equal to half the strength of an H-bond. This leads to the following important conclusion – a liquid molecule can accumulate sufficient agitational energy to break the restraining H-bonds only from multiple, nearly simultaneous interactions with neighboring liquid molecules. The corollary to this conclusion is that in comparison to the frequency of liquid-liquid molecule interactions the occurrence on an interaction sequence which leads to evaporation is rare. This conclusion can be verified experimentally, for example, the steady state evaporation rate during the „no heating” test shown in Fig 9 is 4.7 gm/s m2. The number of liquid-liquid molecular interactions in the upper 2 molecule layers is estimated to be ~8.*1029 interactions per second per m2. Dividing this by the number of molecules evaporated per sec yields a ratio of ~5.0*106.

26

10. The limits of kinetic theory Simplicity is a virtue; if a simple computational model generates satisfactory results, there is no need to complicate it. For many applications, models based on kinetic theory are entirely adequate. For numerous additional situations, the use of kinetic theory can be extended by using appropriate correction factors. However, when simplification becomes the cause of erroneous results and obscures some fundamental physical characteristics, it is no longer acceptable. For kinetic theory this limit is approached in the vicinity of an anisotropic surface and for conditions involving phase change. Kinetic theory can not account for the effect of H-bonds; therefore for water the resulting errors can be especially large. Evaporation or condensation processes are inherently complex, a definitive mechanism remains to be determined, however, using the H-bond based conclusions of the previous section, a probable scenario can be inferred. As has been noted and is shown in table 2, “by themselves” surface molecules are not able to break their constraining H-bonds. They can accomplish this only by using additional agitational energy which is transported toward the interface. For the conditions of the tests shown in Figs 4, 5, 9 and 14 flow regimes in both phases are laminar, agitational energy is thus transported by thermal diffusion along the negative temperature gradients shown in the figures. At the interface an evaporating molecule must accumulate sufficient energy from at least 2, possibly 3, neighboring molecules. This energy must be transferred nearly simultaneously and approximately in the same direction. Not surprisingly such an event has a low probability; this is confirmed by numerical examples. The agitational energy transported toward the interface is invested in the breaking of the H-bonds, in the accepted “energy book keeping” scheme it becomes part of the enthalpy of the evaporated molecule. However, it is in effect a “virtual energy” component, it can be recovered only when the molecule condenses and H-bonds are reformed. The large amount of agitational energy required for evaporation strongly cools the interface, as shown in the noted figures; surface temperature can decrease to such a degree, that agitational energy is transported toward the interface from the vapor region as well. 27

The details of this scenario remain to be verified, but the basic conclusion is clear – the evaporation event is not an interaction of two freely moving particles. In fact, it has no similarity to such an interaction; therefore the distribution of the agitational energy of evaporating molecules can not be Maxwellian5. In Sect 16 similar arguments are employed to reach the conclusion that the agitational energy of condensing molecules can not have a Maxwellian distribution. These conclusions have significant consequences. For example, Eq 1 and 2 evaluate net phase change rates by subtracting evaporating and impacting mass fluxes which are assumed to have Maxwellian distributions of agitational energy. A large data base has been accumulated facilitating the application of these equations under various conditions. That is useful information, however it should be recognized that Eq 2 is in fact an empirical relationship which should be used only under selected conditions and within a limited range. General theory based relationships remain to be developed, their development requires physically realistic distributions of evaporating and condensing fluxes. An additional consequence is that the “temperature gap” and “energy sink” issues, which are perceived as anomalies when it is assumed that the distributions of evaporating/condensing fluxes are Maxwellian, are in fact illustrations that the total energy must be conserved during phase change. Conservation of the total energy during phase transfer must include both the agitational and internal energies. As shown in Fig 8, the agitational energy component represents approximately only 20% of the total energy of a liquid molecule, after an evaporated liquid molecule has joined the vapor phase, this percentage increases to ~60%. Therefore, across an evaporating interface the total and agitational energy can not be continuous simultaneously. Temperature is a measure of agitational energy, on the vapor side of an evaporating interface two types of molecules intermingle: those which are evaporating from the liquid, and molecules diffusing towards the interface from the bulk 5  The emphasis of this conclusion does not imply that other investigators (1,2,3,4) did not recognize that kinetic theory is not readily applicable to the liquid phase. They certainly did. The use of the Maxwellian distribution was justified by reasoning that since at equilibrium conditions, when no net phase change occurs the impacting and evaporating mass fluxes are equal, the distributions of the impacting and evaporating molecules must be equal as well. If it is assumed that the impacting flux has a Maxwellian distribution, this must also hold for the evaporating flux.

28

region of the vapor. In the vicinity of an evaporating interface a temperature measurement will represent a mixture of these molecules; the magnitude of the “temperature jump” depends on the relative size of these two components. This is illustrated in Fig 5, which shows measured temperatures for tests conducted at a vapor pressure of ~560 Pa, but at different bulk vapor temperatures. As seen, though the vapor temperature changes significantly, the interface temperature of the liquid remains almost constant, the “temperature jump” becomes larger because on the vapor side the evaporating molecules mix with more energetic “bulk vapor” molecules. An evaporating surface is not an “energy sink”, but a “sink of agitational energy”. The agitational energy of an average liquid molecule is not sufficient for the dissolution of H-bonds; additional agitational energy must be supplied from neighboring molecules. Figs 5, 9 and 14 illustrate net evaporation conditions during which agitational energy is transported toward the interface through both the liquid and vapor media. Probably the most surprising consequence is one that has not been documented in the literature6. It can be summarized as follows: if the physically realistic distributions of agitational energy for evaporating and condensing molecules are different, then, for conditions where no net phase change occurs Tfi ≠ Tgi at the vapor-liquid interface. The inequality of vapor and liquid surface temperatures is required by the conservation of total energy. A steady state condition implies that condensation and evaporation rates are equal, thus the amount of energy required for the dissolution of H-bonds by evaporating molecules is matched exactly by the amount of energy generated by condensing molecules. If the kinetic energy that is carried away by evaporating molecules is not exactly equal to the kinetic energy of condensing molecules, the total energy transfer rate is not balanced. To maintain a steady state balance, additional agitational energy must be transferred across the interface. Diffusive transfer of agitational energy across an interface is possible only when Tfi ≠ Tgi.

6  If such information exists, it is not available in the readily accessible literature. Probably this can be explained by the circumstance that it is next to impossible to measure a precise liquid interface temperature.

29

11. Alternate computational approaches The final consequence is that for an adequate evaluation of the rate at which phase change proceeds alternative analytical and computational approaches are required. The potentially simplest improvement would be the substitution of physically realistic distributions of evaporating and condensing molecules into relationships similar to Eq 1 and 2. That goal defines one of the principal challenges in this area. Sections 16 and 17 take preliminary steps in that direction by proposing inferred distributions for evaporating and condensing molecules. Sections 13 and 14 analyze the suitability of the powerful finite difference computational approach for analyzing vapor molecule distributions in the transition region near an interface. This section asses the applicability of already developed thermal hydraulic code categories. The remarkable expansion of computer power in the past decades has led to the development of computational approaches which represent matter not as a continuum but at various levels of discretization. For the analysis of thermal hydraulic problems a broad range of computer codes have been developed, some of these codes include the capability to handle two phase flows. The question arises to what extent this computational capability could be applied in resolving the condensation/evaporation issue. A sizable category of codes were developed for specialized thermal hydraulic systems, for example, the RELAP code series (19) is designed to analyze thermal hydraulic behavior of nuclear power reactor cooling systems. The discrete control volumes, or “nodes” of these computational codes represent thermal-hydraulic components or segments of components. The linear dimensions of these nodes are meters or tens of centimeters, within these nodes averaged local equilibrium state properties are used. Their capability to handle a variety of two-phase flow conditions depends entirely on built in theoretical and empirical relationships. The general purpose fluid flow simulation codes (Computational Fluid Dynamics or CFD) are more promising. They are designed to model turbulent 30

flows; consequently discretization of space is much finer and extends down to the turbulent eddy level. Spatially integrated averages of pressure and temperature are employed. In some studies CFD has been adapted to analyze phase change related problems (22), this requires the imposition of empirical or heuristic models of condensation/evaporation processes. CFD methodology can not generate information beyond that which is incorporated in the imposed models, but it could be useful for analyzing the performance of proposed condensation/ evaporation models and benchmarking models against experimental data. The most demanding computational approach is MS (Molecular Simulation) which refines the discretization down to the molecular level. MS encompasses several sub branches, this includes the MC (Monte Carlo) and BD (Brownian Dynamics) methods. The potentially most promising approach, which has in fact already been applied to the condensation/evaporation issue, is MD (Molecular Dynamics). MD models intermolecular forces and their resulting interactions, including H-bond interactions. Like all of the approaches in this category, MD uses probabilistic methods to estimate where and when an interaction occurs, empirical or theoretical probability distributions to estimate how it proceeds. This requires very large computational resources, and even when these are available the time frames that can be analyzed for a statistically significant molecule population is very small, currently it is limited to the nano-sec range. Though the MD methodology comes closest to “simulating” the discrete, molecular nature of physical processes, that does not mean that it is a truly “fundamental” approach. In order to define a “virtual molecule” which can be evaluated mathematically a  number of assumptions and simplifications have to be made. No consensus exists regarding the best choice; in fact, well over 20 such models are currently in use. Their success in predicting fundamental properties of water, such as heat capacity, viscosity and self diffusivity, is mixed. At the current level of development MD methodology is best suited to analyze high probability states which are close to equilibrium; evaporation/condensation interactions do not fit these criteria. Nevertheless, some initial studies in this area have been reported. An example is a study of evaporation (17) from

31

a  ~6.5 nano-m diameter water droplet encompassing nearly 5000 molecules. During a simulation time of 0.5 nano-sec 70 evaporation events are observed. In this time frame the molecules near the surface experience about 10E7 interactions, MD simulation thus confirms that evaporation is a rare event. The interactions of evaporating molecules are presented in detail, they show that during a 1 p-sec time frame an evaporating liquid phase molecule can reform H-bonds up to 10 times, the last 1 or 2 H-bonds are broken in ~0.1 p-sec. In some cases the last H-bond is broken and reestablished several times before the molecule departs. The kinetic energy of the departing molecules is not provided, however, it is evident that the interaction sequence leading to evaporation is not a two particle event. MD simulation thus confirms that the distribution of agitational energy of evaporating molecules is not Maxwellian.

12. Insights from neutron transport theory The evaluation of neutron distributions in a medium which includes H atoms has some striking similarities to the evaluation of the distribution of H2O molecules at low vapor pressure conditions. Neutrons, like vapor molecules, travel over comparatively large distances between interactions, when interactions occur, for the case of neutrons, the dominant energy and direction altering mechanism is scattering from Hydrogen atoms, that is, from scatterers which have the same mass. For the case of H2O vapor molecules, scattering occurs exclusively between entities having the same mass. Scattering events between identical mass particles have unique characteristics: they are anisotropic with a pronounced “forward” scattering bias, depending on the approach angle the faster particle can transfer a large fraction of its kinetic energy including all of it7 to the target. The interaction probability of the two particles depends on the center-of-mass kinetic energy of the two colliding particles. 7 A billiard player does not have to study physics to know this.

32

Decades of effort in the neutron transport area has led to the development of computational techniques which are able to determine neutron distributions in space, energy and direction so successfully that experimental verification is no longer required. An extensive literature exists on the subject, including numerous monographs and text books (20, 21), this review is confined to the question – can these techniques be adapted for evaluating vapor molecule distribution in a region close to an anisotropic boundary? Besides the noted similarities between neutron and H2O vapor molecule interactions, there are significant differences as well; the adaptation of neutronic computational techniques will therefore require considerable effort. Additional information regarding this issue is provided in Appendix B, this section lists several insights which can be useful even if a fully developed computational code is not available. The “3 λ rule” Under appropriate conditions the Maxwellian distribution is used also in neutronic computations. Neutrons which have kinetic energies similar to the atoms from which they scatter are called “thermal” neutrons; they interact with atomic nuclei in a classical “two particle” manner, consequently the spectrum of their kinetic energy approaches a Maxwellian8. Averaged interaction rates can be evaluated easier if an analytical distribution can be employed; this was especially true in the early days of neutronic analysis when computer resources were limited. The Maxwellian distribution was then used; however, it was used only if conditions warranted its application. Analysis verified by experimental data showed that an equilibrium Maxwellian distribution is approached in regions which are at least three mean free path distance away from an anisotropic surface. For the case of vapor molecules an analogous transition layer is called the “Knudsen” layer. Its definition is less precise, in some publications it is considered to extend over a single mean free path. It is proposed that the “3 λ rule” should be used for the definition of non-equilibrium transition regions near interfaces. The 1/Vc dependence. No particle, not even a neutron is simply a “hard ball” object, it is surrounded by short range force fields which diminish as a function of 1/R², 8  Since a fraction of the neutrons are not just scattered but absorbed, the distribution is not a “true” Maxwellian, however, the similarity is sufficiently close to be useful.

33

where R is the distance to the center of the particle or molecule. The time that the particles spend in the vicinity of the 1/R² force field is proportional to 1/Vc, where Vc is the relative speed between the two interacting entities. As the Vc speed decreases, the time that the particle spend crossing the range of the force fields increases. With decreasing speed the momentum of the particles also decreases, this means that weaker force fields can deflect the particle from their original path. In effect, at lower Vc values the “effective size” of the particles increases, stated in probabalistic terms, their interaction probability increases. For simple particles like neutrons, the 1/V dependence is so exact that it is used for extrapolation of experimental data, this is illustrated in Fig 11 which superimposes the 1/Vc dependence of interaction probability (the “cross section”) over a Maxwellian distribution evaluated at 273K. Of course, H2O molecules are more complex then neutrons and nuclei; they are dipoles with a fluctuating polarity. The interaction of two H2O vapor molecules is consequently considerably more complex; nevertheless, some form of 1/Vc dependence of their interaction probability will be present. For example, the time that molecules which approach each other at a relative speed of

Fig 11.

34

Table 4 H2O Vapor Molecule Dimensional Characteristics

20 m/s spend in their mutual vicinity is 100 times longer then for molecules that do so at 2000 m/s. This longer time will influence their interaction probability. As illustrated in Fig 11 both relative speeds are possible. Table 4 lists dimensional characteristic of the H2O vapor molecule, Appendix B illustrates how these dimensions can be interpreted using definitions and methods employed in neutronic computations. Advantages of finite difference methods. Theoretical models of physical reality and the computational methods used to generate results based on these models should, in principle, be independent. Generally this is true, however in some cases the computational method does influence results. For example, because the evaluation of phase change rates relied on the Maxwellian distribution, solution methods employing continuous functions were preferred. This preference is influenced also by history, in the first half of the last century when the initial development took place computers were not available, there was thus no practical alternative. The major developments in the neutronic field occurred later; computers were available from the beginning and they became the dominant analytical tool. A basic computational technique employed in a large variety of ways is the “finite difference” method, some details regarding this widely used technique are presented in Appendix B, in this 35

Fig 12. Maxwellian distr. at 275 K Number density vs molecular speed. 4 energy group number density

section a numerical example is employed to show how it could be used in evaluating vapor molecule distributions in the vicinity of an anisotropic interface. The “50 C heated” test shown in Fig 5 is chosen as an example. In order to avoid meta-stable states of water the pressure is specified as 700 Pa9. The equilibrium Maxwellian distribution of the vapor ~3 mm away from the interface is shown in Fig 12. At this location Tg = 323K, ρg = 0.0047 kg/m3. Using Eq. 3 and the equilibrium Maxwellian distribution the pressure is evaluated to be 700 Pa.



3

4 Where F(v) is the Maxwellian distribution expressed as a function of vapor molecule velocity. 9 The triple point for water is at 611Pa, T = 273K

36

In the finite difference approach pressure and temperature are evaluated using equations 5 and 6 5

6 Where Γ(i) is the molecular density fraction, n the total number of energy groups. Fig 12 illustrates how the distribution could be represented by a simple four group segmentation of the velocity space, the respective Γi and the corresponding rms velocities are given in table 3. The pressure computed using eq. 5 is also 700 Pa. 4 energy grp. width in m/s

Eq. Maxw. @ Num. dens.

T = 323K Urms (m/s)

Non. eq. imposed distr. Num. dens. Urms (m/s)

0 to 400

0.216

300.

0.45

302.

400 to 800

0.553

603.

0.35

596.

800 to 1200

0.210

951.

0.15

938.

1200 to 1600

0.021

1320.

0.05

1302.

Table 3. Finite difference segmentation of eq. Maxwellian distribution and an imposednon-equilibrium spectrum in the 3λ transition region

This example illustrates that when properly employed the finite difference approach can duplicate results of continuous functions. The advantage of the discrete eq 5 over the continuous eq 3 is that it can accommodate all possible spectra, whereas eq 3 is useful only if the distribution F(v) is known. This is illustrated by fig 13 which shows a non-equilibrium distribution near an interface emitting a constant stream of low energy vapor molecules. The distribution close to such a boundary condition is strongly space dependent; at increasing distances from the interface it will approach the equilibrium Maxwellian distribution. Neither the distribution nor its special evolution can be represented by continuous functions however; it can be approximated using the finite difference method. This requires the segmentation of space as well of energy, and an iterative solution of vapor mass, energy and momentum balances. Assume that fig 13 represents the velocity

37

Fig 13. Equilibrium and non-eq. vapor molecule distributions. Pressure 700 Pa, Eq T = 275K

distribution of vapor molecules in a spatial node located in the non-equilibrium 3 λ region. Using eq 6 and the non-equilibrium Γ(i) and Urms values of table 3 the temperature in this node is computed to be 275 K, at a pressure of 700 Pa this yields a vapor density of 0.0055kg/m3. The finite difference method has the flexibility to accommodate a broad range of boundary conditions; another advantage is that in the calculation of mass/energy balances it can employ energy group dependent interaction probabilities. Consider again fig 13, the average velocity of the molecules in group 1 is 200 m/s, in group 4 it is 1400 m/s, this velocity difference influences molecular interaction probabilities and parameters as the diffusivity and heat conductivity. Using the finite difference approach this variation can be taken into account. For illustration purposes a simple 4 energy group segmentation is adequate, presently available computational resources make it possible to use almost arbitrarily small finite difference segments. Boundary conditions can then be represented more precisely and the evaluation of averages becomes simpler.

38

14. Experimental verification A limited amount of experimental data for the verification of vapor molecule distributions in the 3*λ transition region near an interface is already available. The traces presented in Fig 4, 5 and 9 have a resolution which shows the temperature variation within a mean free path of the vapor molecules. For some of the “heated” tests conducted at very low vapor pressures the 3*λ transition region extends for over 1 mm and includes a sufficient number of measurement points that the temperature gradient within the transition region is depicted in detail. Fig 14 is an example of a test in which the vapor at a distance of 3mm from the interface was heated to almost 80 C. Vapor pressure for this test is 218 Pa, the “temperature gap” ~14 C, liquid temperature at the interface ~-12 C10. A steady state evaporation rate of 1.23 gm/s m² is measured; this implies that bulk flow of vapor is present, however, the Re number in the vapor region has a magnitude of only ~13, thus the flow is definitely laminar.

Fig 14 . “Heated” test from Badam et al (14) P = 218 Pa 10 Though the conditions for this test are well bellow the triple point of water, according to the authors, the water remained liquid.

39

In addition to the measured temperature, fig 14 includes the extrapolation of the negative temperature gradient established far from the interface where an equilibrium Maxwellian distribution exists. The thermal energy flux diffusing from the vapor region toward the interface is equal to the product of the temperature gradient and the thermal conductivity. The evaporation rate is constant, there are no energy “sinks” in the transition region, therefore the rate at which thermal energy diffuses toward the interface must be constant. As shown in the figure, that does not seem to be the case, as the interface is approached, the measured temperature deviates toward smaller values. In fact, the negative temperature gradient increases as the interface is approached11, according to the classical -(dT/dz)ωg relationship, the thermal energy flux must increase. However, that can not be true, there are no energy sources in this part of the vapor region. The thermal energy flux is constant throughout; the only possible explanation is that in the transition region above an evaporating interface g, the conductivity of the vapor, becomes smaller. Like all other vapor properties, the thermal conductivity is measured at equilibrium conditions; it represents the conductivity of vapor molecules which have a Maxwellian distribution. The deviation shown in fig 14 illustrates that in the 3*λ transition region the molecular distribution is not at equilibrium The thermal conductivity of a vapor can be expressed as follows: ωg = ⅓ Ng Uav λ Cv 7 Ng appears also in the denominator of λ, Cv is constant, therefore ωg is directly proportional to Uav, the average velocity of the vapor molecules. This leads to the following important conclusion – compared to an equilibrium Maxwellian distribution an excess of low energy molecules is present in the transition region and this excess increases as the interface is approached. This confirms the conclusion reached in sect 10 that the distribution of agitational energies of evaporating molecules can not be Maxwellian. 11  Using hand book conductivity values it is estimated that ~45 % of the energy required to maintain the measured evaporation rate is conducted toward the interface from the liquid region, ~10% from the vapor. Fig 14 represents energy transfer only in the vertical direction. The balance of the required energy therefore must be supplied along the traverse dimensions. The surroundings are warmer, thermal energy from the transverse directions diffuses toward the interface, therefore, for ideal adiabatic 1-D test conditions the difference between the measured and extrapolated traces would be larger.

40

16. Inferred spectrum of evaporating molecules The energy spectra of evaporating and condensing H2O molecules remain to be determined, however some inferences about their probable characteristics can be made. The agitational energy of liquid molecules is smaller then the binding energy of H-bonds, this implies that interaction with a single neighboring liquid molecule is not sufficient to break it, therefore evaporation must be the result of multiple simultaneous (or nearly simultaneous) interactions. Such a process is indeed observed in MD simulations of evaporation (17). Most of the energy that an evaporating surface molecule accumulates from its neighbors is required to break the last H-bond, consequently the residual kinetic energy which propels the molecule from the interface is probably quite small. If this “extra” kinetic energy does not exceed the attractive van der Waal forces at the liquid interface, the molecule recondenses, the most probable energy of the evaporating molecules therefore should be a little larger then this “minimal escape” energy. Beyond that the kinetic energy of the evaporating molecules should decrease steeply. This inference is supported by the test results shown in Fig 14 the concentration of evaporating molecules is highest near the interface, their relatively small kinetic energy accounts for the decreased heat transfer coefficient near the interface. These observations provide the basis for the inferred agitational energy spectrum of evaporating molecules shown in Fig 16. The evaporating molecule should have some “escape velocity” to overcome the non-fluctuating van der Waal forces, after this threshold is reached the evaporation probability rises steeply to a peak and subsequently decays. The proposed energy spectrum shown in Fig 16 contradicts popular descriptions of the evaporation mechanism presented in some text books. For example, the “Wikipedia” entry on the subject includes the following observation: „As the faster-moving molecules escape, the remaining molecules have lower average kinetic energy, and the temperature of the liquid decreases. This phenomenon is also called evaporative cooling“.  The statement that evaporation leads to cooling is certainly correct, however the explanation of this phenomenon suggests that the evaporating molecules 41

Fig 16. Infered translational energy spectrum of evaporating molecules

have a high translational energy. If that were the case, there would be an excess of high energy molecules near the interface and the measured temperature trace in Fig 14 would deviate in the direction of higher temperatures.

17. Inferred spectrum of condensing molecules In the classic kinetic theory approach the evaporation and condensation mechanisms are assumed to be symmetric, this is illustrated in the CV-2 schematic shown in Fig 2. The arrows shown in the figure imply that the condensation mechanism is a mirror image of evaporation. For both processes a KT determined mass flux approaches the interface from opposite directions, at the interface some of the impinging molecules manage to cross it, a fraction is reflected. In Sect. 10 and 11 it is shown that the mechanisms leading to evaporation and condensation are completely different, therefore such a representation of the phase change processes is not correct. Evaporation of a liquid molecule depends almost entirely on

42

the liquid phase; sufficient thermal energy must be transferred to the evaporating liquid molecule to break the H-bonds. The evaporating molecule can be re-absorbed by the liquid phase if it does not have an adequate “escape velocity”. However , once it has passed this barrier “reflection” back into the liquid is a rare event. Condensation proceeds differently. The rate at which vapor molecules impinge upon the interface depends on conditions within the 3*λ transition region. The spectrum of these molecules is not known, however the magnitude of the impingement rate is probably reasonably close to the magnitude evaluated by kinetic theory. In their comprehensive review of published experimental results Marek & Straub (10) showed that when the impinging molecular flux is evaluated using a Maxwellian distribution, the experimentally verified fraction which is actually condensed varies by orders of magnitude. The average condensed fraction appears to be somewhere in the vicinity of ~0.1. Thus, according to experimental evidence, ~90% of the impinging vapor molecule will be reflected. What processes at the molecular level can account for a ~0.9 probability that an impinging molecule is reflected back into the vapor? The lateral cohesiveness of the surface tension layer probably plays an important role. The H2O molecule is a rotating dipole thus an approaching molecule is both attracted and repelled by electrostatic forces. Which tendency predominates will depend on the molecular velocity and angle of approach. Another hypothetical “reflection” mechanism could proceed as follows: 1. A condensing vapor molecule becomes a true liquid molecule only after it has established at least one H-bond. 2. The newly created H-bond transfers a large amount of thermal energy to neighboring liquid molecules, in effect creating a microscopic “hot-spot” of strongly agitated liquid molecules. 3. The strongly agitated liquid molecules increase the probability that a near by ­H-bond is broken and the same or a neighboring molecule is ejected into the vapor region. Several recent publications have analyzed the condensation rate of vapor molecules. Louden et al. (18) report the results of an MD simulation in which a best estimate value for the condensation coefficient of 0.77 was obtained. The authors 43

observe in their study that a condensation fraction as low as 0.2 can be obtained only if unrealistic parameters for the water molecule model are imposed. Tsuruta et al. (7) proposed a molecular condensation probability12 model based on heuristic arguments for the condensation rate of Argon. The authors reasoned that the condensation probability of an impinging Argon vapor molecule should depend on the vertically directed kinetic energy component and on the liquid temperature at the interface Q = ψ{ 1 – w exp(-Uc2 /(2*Rm*Tfi)} 8 Where Uc is the vertical component of molecular velocity, Q is the Uc dependent condensation probability of impinging vapor molecules, y and w are fitted constants. Bond & Struchtrup (12) extended the application of this model to water. Fig 17 shows a molecular kinetic energy dependent condensation evaluated using Eq. 8 (blue trace), and a probability proposed in this study (red trace). In addition to the “penetrability” argument employed by Tsuruta the proposed distribution takes into account the 1/v dependence of molecular interactions and the experimental evidence that most vapor molecules which impinge on the interface are reflected. To accommodate these considerations, the kinetic energy range of impinging molecules is divided into three segments: a “slow”, a “fast” and a broad “reflective” segment in the middle.

Fig 17. Infered molecular condensation probability 12  Not to be confused with the “condensation coefficient”, which is the product of the molecular “condensation probability” and the molecular flux impinging upon the interface.

44

From the viewpoint of the “slow” vapor molecules, the liquid molecules at the interface are held together by significantly stronger bonds than the kinetic energy of the approaching vapor molecule. There is thus no possibility to disturb these bonds, however, since the molecule moves slowly, the alternating attractive/repulsive forces average out allowing the vapor molecule to approach within the range of the constant van der Waal forces. At this point the 1/V trend discussed in Sect. 13 dominates, and the vapor molecule is attracted to the interface. At the opposite end of the spectrum the “fast” vapor molecules have kinetic energies which are within the range of the strength of the H-bond, in addition their approach times match the dipole fluctuation frequency. The vapor molecule interacts not with a layer of closely inter connected molecules, but rather with a specific molecule in that layer. It is assumed that this effect becomes noticeable when the kinetic energy of vapor molecules reaches ~1/3 of the H-bond strength; higher molecular kinetic energies increase the interaction probability. Finally, empirical evidence indicates that most of the impacting molecules are reflected, this means that for a broad range of kinetic energies the

Fig 18. Infered condensation spectrum

45

condensation probability must be small. Most of the impinging molecules are to be found between the “slow” and the “fast” kinetic energy segments. The condensation probability in this central segment is, therefore, assumed to be below 0.1. The consequences of the two proposed molecular condensation probabilities is illustrated in Fig 18 which shows the product of the condensation probability and the impinging molecular flux, the integral of this product yields the “condensation coefficient”. Fig 18 shows that the resulting condensation coefficients are quite different, the coefficient proposed in this study is biased toward low energy molecules, the condensation coefficient evaluated using Eq. 8 exceeds the ~0.1 best estimate value based on empirical data. At this point in time, the relevant question is not which of the two relationships is a “better” representation of reality, rather their disagreement points out a serious lack of knowledge. The evaluation of theory based phase change rates requires physically realistic condensation and evaporation spectra; these in turn require a better understanding of the condensation and evaporation mechanisms at the molecular lever.

18. Conclusions, proposals and recommendations Familiarity can hide complexity; unexamined complexity can lead to misunderstanding. Evaporation and condensation of water are examples of phenomena that are certainly familiar, however, because of their inherent complexity, they are still not completely understood. This study provides an overview of this complex issue. The development of analytical approaches is illustrated using a sequence of progressively more comprehensive “control volumes”. Quantitative estimates of liquid-liquid and vapor-liquid molecule are used to illustrate that at low pressures evaporation is dependent almost entirely on the liquid phase. The influence of the exceptionally strong hydrogen bonds in liquid phase water is evaluated.

46

In the course of the overview a number of conclusions are reached, inferred distributions of evaporating and condensing molecular spectra are proposed and recommendation regarding further development of analytical approaches are presented. This closing section presents a summary of the conclusions, proposals and recommendations. If they are considered from the view point of specific disciplines, some of the listed observations might seem obvious, they are included for completeness. Conclusions: *) The evaporation of a liquid molecule requires multiple, nearly simultaneous interactions with neighboring liquid molecules. Therefore, the agitational energy spectrum of evaporating molecules can not be Maxwellian. *) A 3 mean free path (3*λ) distance from an anisotropic boundary is required for the establishment of an equilibrium Maxwellian distribution. *) In the liquid phase H2O molecules are restrained by strong, short-lived hydrogen bonds which are the cause of the large internal energy component. It is shown that the “temperature jump” at an evaporating interface is not an anomaly, but the consequence of the greatly differing ratios of agitational to total energy for liquid and vapor molecules *) The interaction probability of potentially colliding vapor molecules depends on their relative (center-of-mass) kinetic energy. *) Closed form analytical approaches which represent the entire population of vapor molecules in a single relationship can not yield physically meaningful results in the 3*λ transition region. Proposals: *) Theoretical and experimental arguments suggest that the spectrum of evaporating molecules is definitely not Maxwellian, it is strongly biased toward low velocity molecules. An inferred evaporation spectrum is proposed. *) The agitational energy spectrum of condensing molecules is not known. The spectrum of the molecules impinging on the interface is influenced by the

47

distribution of molecules in the 3*λ transition region. Experimental evidence suggests that the condensation probability of molecules which do impinge upon the interface is quite low, on the order of ~0.1. An inferred condensation probability and the resulting spectrum of the condensing molecules is presented in the study. *) The agitational energy spectra of the evaporating and condensing molecules are different. Very probably the total agitational energy due to mass transfer across the interface is different as well. This implies that for steady state conditions when the condensing and evaporating streams are equal the liquid and vapor temperatures at the interface will differ, that is Tfi ≠ Tgi. Recommendations: *) Useful insights for evaluating vapor molecule distributions can be found in the neutronic field, of special interest are the well developed finite difference methods including the discrete ordinance computational techniques. Material for the adaptation of finite difference techniques for evaluating mass/energy distributions of vapor molecules is provided in Appendix B. *) The most promising analytical approach is MD (molecular dynamics), however, considering presently available computational resources MD is not suited for evaluating statistically rare molecular events such as evaporation. MD can be used for the evaluation of steady state molecular distributions in the 3*λ transition region for cases where the evaporation spectrum is provided as an input boundary condition.

48

Nomenclature Av Cp Cv D G hf hg hfg J M N Pf Pg Rg Rm Tf Tfi Tg Tgi q U Uc

Avogadros number heat capacity at constant pressure (J/kg K) heat capacity at constant volume (J/kg K) diffusion coefficient (cm) Gibbs free energy (J/kg) specific enthalpy of fluid (J/kg) specific enthalpy of vapor (J/kg) latent heat of evaporation (J/kg) mass flux (kg/s m2) molar weight of H2O (kg) molecular density (#/cm3) pressure in liquid (Pa) pressure in vapor (Pa) universal gas constant (J/mol K) gas constant for H2O (J/ kg K) liquid temperature (K) liquid temperature at interface (K) vapor temperature (K) vapor temperature at interface (K) heat flux (W/m2) molecular velocity (m/s) vertical molecular velocity component (m/s)

49

Greek symbols αe evaporation coefficient αc condensation coefficient Q condensation prob. of impinging vapor molecule ρg density vapor (kg/m3) ρf density liquid (kg/m3) λ mean free path (nm or μm) ω conductivity (W/m K) σ microscopic cross section (barns) Σ macroscopic cross section (1/cm) Γ molecular fraction ζj – m molecular scattering probability from group j to m

Subscripts f liquid property g vapor property gi vapor property at interface fi liquid property at interface rms root mean square average

References

with Irreversable Thermodynamics, Kinetic Theory and Statistical Rate Theory” Master thesis, Univ. of Victoria (2000) 9) Ward C. A., D. Stanga “Interfacial conditions during evaporation or condensation of water”, Physical Rev. E, Vol. 64, pp.  05159091– 05159099. (2001) 10) Marek R., J. Straub “Analysis of the evaporation coefficient and the condensation coefficient of water”, Int. J. Heat Mass Transfer 44 pp. 39–53 (2001) 11) Meland R., T. Ytrehus “Evaporation and condensation Knudsen layers for nonunity condensation coefficient” Phys. Fluids 15. 1348 (2003) 12) M. Bond and H. Struchtrup “Mean evaporation and condensation coefficients based on energy dependent condensation probability”, Physical Review E 70, 061605 (2004) 13) S. Popov, A. Melling, F. Durst. C.  A.  Ward “Apparatus for investigation of evaporation at free liquidvapor interfaces”, Int. J. of Heat and Mass Trf., 48, pp. 2299–2309, (2005) 14) Payam Rahimi, Ch. A. Ward “Kinetics of Evaporation: Statistical Rate Theory Approach” Int.  J. of Thermodynamics Vol. 8 pp. 1–14 (March 2005)

Selected chronological condensation/ evaporation publications 1) Knudsen M. “Die Maximale Verdampfungsgeschwindichkeit des Quecksilbers” Ann. Phys. Chem. 47, (1915) 2) Alty, T. “The Reflection of Vapor Molecules at a Liquid Surface”, Proc. Roy. Soc. A. Vol. 131, (Feb. 1931) 3) Schrage R. “A Theoretical Study of Interphase Mass Transfer” Columbia Univ. Press, New York, (1953) 4) Shankar P. N. “A kinetic theory of steady condensation” J. Fluid Mech. 40, pp 383–400, (1970) 5) Y.  P. Pao “Application of Kinetic Theory to the Problem of Evaporation and condensation” Phys. of Fluids, 14 (2) (1971) 6) C.  A Ward and G. Fang “Expressions for predicting liquid evaporation flux: Statistical rate theory approach”, Physical Review E, 59, pp. 429–440, (1999) 7) Tsuruta T.,  H. Tanaka and T. Masuoka“, Int. J. Heat Mass Transfer, 42, 4107, (1999) 8) Bond M. “Non-Equilibrium Evaporation and Condensation Modeled

51

Other references:

15) Badam V.  K, V.  Kumar, F. Durst, K.  Danov “Experimental and theoretical investigation of interfacial temperature jumps during evaporation”, Experimental Thermal and Fluid Science, 32, pp. 276–292 (2007) 16) Deng S., W. Wang, S. H. Yu “Pointwise Convergence to Knudsen Layers of the Botzman Equation” Comp. Math. Phys. 281 (2008) 17) P. Mason “Molecular dynamics study on the microscopic details of the evaporation of water” J Phys Chem A. 2011 Jun 16; Vil. 115 (23), pp. 6054–6058. doi: 10.1021/ jp1104517. Epub 2011 Feb 15. 18) P. Louden, R. Schoenborn, C. P. Lawrence “Molecular dynamics simulations of the condensation coefficient of water” Fluid Ph. Es. Vol. 349, (July 2013)

19) REALP5/MPD3 Code Manual, NUREG/CR-5535, Idaho Natl Eng. Lab. (1995) 20) R. L. Murray “Unified Neutron Transport Theory”. Springer Link, (1979) 21) E.  E. Lewis “Computational Methods of Neutron Transport”, American Nuclear Society, (1993) 22) K. Almenas & R. Lee “Nuclear Engineering. An Introduction”, Springer-Verlag (1992) 23) P. Vidal-Lopez, V. MartinezAlvarez, B. Gallego-Elvira, B. Martin-Gomez “Determination of synthetic wind functions for estimating open water evaporation with CFD “, Hydrological Processes, 26, pp. 3945–3952, (2012)

52

Appendix A

Limitations of SRT and TIP theories The SRT (Statistical Rate Theory), TIP (Thermodynamics of Irreversible Processes) and KT (Kinetic Theory) approaches share several common assumptions: 1) The kinetic energy distribution of the vapor molecules can be represented by a single analytical expression1. 2) The properties of the vapor (e.g viscosity, diffusivity, thermal conductivity etc.) can be obtained by evaluating appropriate averages using this expression. The consequences of these assumptions for the evaluation of phase change rates using KT are examined in the main part of this study; this Appendix considers the limitations that they impose on results obtained using SRT and TIP theory assumptions. The three approaches differ in the amount of information that is used for the evaluation of the physical state of liquid and vapor phases. KT considers only the agitational energy of the molecules, phase change rates are obtained by subtracting integrated half-range Maxwellian distributions at Tfi, and Tgi, therefore, net evaporation is possible only when Tfi > Tgi. The SRT and TIP approaches use additional physical properties, specifically – Cpg, Cpf and hfg. Thus more information regarding the physical state of the phases is available. The heat capacities incorporate contributions due to internal energy, the latent phase change energy hfg includes the energy required for the dissolution of H-bonds. Both approaches differ from KT in that they use the 2-nd law of thermodynamics to determine the direction of phase change. SRT relies entirely on the thermodynamic state of the liquid and vapor phases, the TIP approach employs additional empirical information regarding mass or energy fluxes. 1  In practice these properties are measured in a vapor region in which an assymptotic Maxwellian distribution prevails. Subsequently they are applied also in regions in which the agitational energy distribution of the vapor molecule population deviates from an assymptotic Maxwellian distribution.

53

One of the claims of the SRT approach as presented by Ward & Fang (6) is that it does not require any fitting parameters, that is, both the direction and rate of phase change is evaluated in a closed form expression based on the state of the liquid and vapor phases. Formally this is true; however the phase change rate is determined using an imposed assumption, and is consequently limited by the physical appropriateness of this assumption2. The phase change rate obtained using eq. 1 based on KT, has a similar shortcoming3, for both methods the absolute value of the calculated phase change rate is not reliable and can differ significantly from measured values. Bond and Struchtrup (12) have observed that though the initial assumptions of the KT and SRT approaches differ, the coefficients of the phase change rate equations are basically identical. This can be shown after some adjustments of the equations which convert them into a form consisting of a dimensional coefficient and a dimensionless relationship which is proportional to the evaporation/condensation driving force. The KT based phase change mass flux (eq. 1) can be rewritten as follows: A-1 Combining eq. 52 and 53 in ref 12, the SRT mass flux is given by: A-2 Since RM = Rg/M and Psat(Tf) is very close to Pfi the coefficients, of the JKT and JSRT mass fluxes are seen to be identical, the assumption used by SRT that entropy is maximized in a small volume at the interface leads to results which are equivalent to the KT assumption that the distribution of agitational energies of condensing and evaporating molecules are identical and that the evaporation and condensation fractions are equal to 1. 2  As presented in ref. 6, the SRT approach assumes that in a small volume surrounding the liquid/vapor interface thermodynamic equilibrium is reached when the entropy within the volume is maximized. The transition occurs in a small time interval δt > (2π h/( ∆Ef +∆Eg)). The limitation of this assumption is that the energy change in liquid and vapor phases are accorded the same weight. It is shown in this study that the evaporation rate at low pressures is determined almost entirely by the liquid phase. 3  This is discussed in sect 10 to 15 of the study. As noted, the major shortcoming is that evaporating and condensing molecular fluxes are assumed to have a Maxwellian distribution.

54

The difference between the KT and SRT approaches is represented by the dimensionless expression in the square brackets of A-1 and A-2. The evaporation rate becomes equal to the condensation rate when these expressions are equal to zero. The outstanding question can now be restated as follows: what is the relationship between the Tfi and Tgi temperatures when the terms in the square brackets sum to zero? For the case of the KT derived eq. A-1 this occurs when Tgi = Tfi, for eq. A-2 equilibrium is reached when ∆Ssrt = 0, where ∆Ssrt is the phase change entropy difference:

A-3

A-4 ∆Ssrt is identified by an srt subscript because it does not represent the complete phase change entropy difference, this issue will be considered later in the section dealing with the TIP approach, this section analyzes the dependence of ∆Ssrt on Cpf, Cpg and hfg As shown in Eq. A-3, ∆Ssrt includes two differences: the difference between the Gibbs free energy of the liquid and vapor phases, and the entropy difference due to the transfer of vapor enthalpy from Tfi to Tgi. When net evaporation prevails the first term is positive because though the Gibbs energy is negative, the entropy of vapor molecules is larger then the entropy of molecules bound in the liquid phase. On the other hand, the entropy of vapor molecules decreases as temperature rises, therefore if Tgi > Tfi, the second term of Eq. A-3 is negative. Consequently, depending on the relative magnitude of these two terms evaporation is possible (that is ∆Ssrt > 0) even when Tgi > Tfi. The range of temperatures for which, according to the SRT approach, this holds depends on the phasic properties Cpg , Cpf , hfg. This dependence is illustrated by evaluating the components of Eq. A-3 as a function of Tfi, Tgi temperature differences. Before proceeding, some comments are in order. The numerical evaluation of Eq A-4 is straightforward, however, the specification of the input parameters, and the

55

required phasic state properties hf, hg, Sf, Sg, Gg, Gg presents difficulties. The numerical evaluation of Eq. A-3 involves a sequence of subtractions. For cases which are of interest, the final result is obtained as the difference of large numbers having a similar magnitude. The precision of physical parameters usually does not extend beyond 3 decimal places; a result obtained by subtracting numbers which differ only in the third or fourth decimal place is not meaningful. Ward and Fang (6) present a sensitivity analysis which shows that the measurement precision which would be required to obtain meaningful condensation rates using Eq. A-2 is not attainable. Consequently Eq. A-2 is verified by comparing measured and experimental saturation pressures rather then condensation rates. SRT is an interesting alternative theoretical approach but of limited use for a practical evaluation of physical condensation related parameters. The theoretical limits of this alternate approach can be investigated by employing a purely numerical procedure. This can be achieved by using nominal, “steam table” magnitudes of h and S at measured experimental temperatures and evaluating their temperature dependence by employing constant phasic heat capacities. Cpf and Cpg are weakly dependent on temperature, therefore for limited temperature ranges this should yield meaningful results. The phasic state functions then can be determined as follows:



A-5

A-6

Where Tf = Tfi, Tg = Tgi, the reference temperature To is set at 273.15 K. For a range of temperatures Tg – To = ~20K, the above equations reproduce “steam table” values with a 4-th decimal place accuracy if the following phasic properties are employed: Cpf = 4198 J/kg K Cpg = 1855 J/kg K hfg = 2.501E6 J/kg A-7 If the definitions shown in A-5 and A-6 are substituted in Eq. A-3, the terms

56

making up ΔSSRT can be rearranged to obtain three relationships which depend only on one of the properties listed above. ∆Ssrt can then be expressed as follows: A-8 Where: A-9

A-10 A-11

As seen, Ff represents the entropy contribution determined only by Cpf, Fg only by Cpg, Ffg only by hfg. Experimental measurements obtained by Badam etal are used for a numerical comparison. Table 1 in ref 14 provides the following data: Tfi = 275.8K Tgi = 279.8K P = 736 Pa For this test inlet water temperature Is 30C, the measured steady state evaporation rate = 0.578 gm/m2s. Note that net evaporation proceeds though the vapor temperature immediately above the interface is 4 degrees higher then the liquid temperature. This implies that the vapor – liquid phase entropy difference ∆Sexp is positive. The vapor temperature at which ∆Sexp = 0 and no net phase change occurs, must be above Tgi = 279.8 K. Substituting the experimental Tgi temperatures into Eqs. A-9 to A-11 yields the following entropy change values: Ff(Tfi) = –0.19 J/kgK g(Tgi) = 45.26 J/kgK Ffg(Tgi) = –129.6 J/kgK Summing the three component parts: ΔSSRT(Tfi ,Tgi) = –84.6 J/kgK

A-12

The numerical example shows that at experimental conditions for which net evaporation is measured ΔSSRT is negative; the SRT approach thus predicts net condensation. The dependence of ΔSSRT on the vapor temperature above the interface is illustrated in Fig A-1 57

Fig A-1. Dependence of the STR vapor – liquid phase entropy difference on Tgv

(Evaluated using Eq. 8, at Tfi = 275.8K To = 273.15K ) In Fig A-1 ΔSSRT and the entropy differences dependent on Cpg and hfg are plotted against a variable vapor temperature. It shows that the equilibrium state, at which ΔSSRT = 0 and thus no net phase transfer occurs is predicted to be 276.5K, that is 0.7 K above the interface temperature of the liquid. The figure thus illustrates that the SRT approach can indeed predict net evaporation for conditions when Tgi > Tfi. That is an improvement over the capabilities of KT, however, this improvement is based on questionable assumptions and is of little practical value. As has been noted, the defining assumption employed by the SRT approach, is that thermodynamic equilibrium is established within a micro volume at the interface which encompasses both liquid and vapor molecules. The numerical example shown in A-12, illustrates the consequences of this assumption: since the energy of vapor phase molecules is higher, the vapor dependent Fg(Tgi) and Ffg(Tgi) terms are correspondingly larger then the liquid dependent Ff(Tfi) term.

58

As a result, the magnitude of ΔSSRT(Tfi ,Tgi) is dominated by the entropy of the vapor phase molecules. In sec. 8 to 10 of this study it is shown that for evaporation of water the opposite is true, that is, the rate of evaporation is determined predominantly by the liquid phase. That is an illustration of the fact that the second law of thermodynamics can predict only the direction of a process like phase change, but not its rate. It can not do so because it does not consider the actual mechanism by which a liquid molecule becomes a vapor molecule. An additional shortcoming is that, as defined by Eqs A-6, the latent heat of evaporation hfg is simply added to the vapor enthalpy. As a result phasic enthalpy (consequently also entropy) differences depend on the choice of the reference temperature To. Bond & Struchtrup (12) propose that the To dependence can be eliminated by redefining vapor entropy so that the contribution of the hfg enthalpy does not depend on Tg. The vapor entropy expression is redefined as follows:



A-13

A final shortcoming of the SRT approach is that ΔSSRT, as defined by Eq. A-3, does not represent the total phase change entropy difference. Consider again the Badam etal test used for a numerical comparison, steady state evaporation is observed even when Tgi – Tfi = 4K. Temperature differences lead to energy transfer, therefore, in addition to the evaporative mass transfer, energy transfer must occur as well. Since Tgi > Tfi, thermal energy is being transferred from the vapor to the liquid. Energy transfer over a finite temperature difference leads to an increase of entropy. The TIP approach (12) considers this additional entropy component, as shown in Eqs. A-14, the total entropy difference is given by: A-14 Where J is the mass flux, qg the heat flux, α and β are dimensional “accommodation coefficients” which are required for unit consistency. 59



A-15

Eq A-15 illustrates that ∆SSRT is not the total entropy difference, it represents the “mass flux” component. As defined by the above relationships, the TIP approach requires additional information, qg must be known before the entropy difference and thus the direction of the phase change process can be inferred. Fig A-2 presents all of the entropy difference components for the same 736 Pa test shown in fig A-1.

Fig A-2. Dependence of the TIP vapor – liquid phase entropy difference on Tgv

The heat flux qg from the vapor region has to be estimated, it probably does not exceed ~5% of the total energy required for evaporation, in Fig A-2 this value is increased to 15%. The positive Fq(t) term shifts the vapor temperature at which ∆STIP > 0 closer to the experimentally observed value but does not reach it. To summarize this brief overview – the SRT and TIP approaches are interesting attempts to deal with the “temperature jump” discrepancy by using the second law of thermodynamics. The 2-nd law is a powerful tool for determining which physical processes are not possible; however, it is not suited for predicting the rates at which thermodynamically possible processes occur. The examples provided in this appendix illustrate this characteristic. 60

Appendix B

Adaptation of finite difference methods for the evaluation of H2O vapor molecule distributions The adaptation of finite difference techniques for the determination of mass, energy and momentum balances in vapor regions would increase computational flexibility and accuracy. Finite difference methods are able to deal with boundary conditions at which spectral and directional distribution of vapor molecules reflect physical reality, and they can take into account the dependence of molecular interaction probability on the relative speed of the interacting molecules. The advantages are thus evident, however the development of this analytical capability presents a number of challenges. A major one is the generation of a database which specifies molecular interaction probabilities between molecules of the finite difference energy groups. For neutron-atom interactions a large database of experimentally verified energy dependent interaction probabilities is available, for vapor-vapor molecule interactions this type of information is sparse. This Appendix presents information for the generations of such a database and outlines the developmental steps needed to adapt finite difference techniques for analysis of molecular distributions in vapor regions.

B-1 Evaluation of vapor molecule interaction probabilities The development of a database of kinetic energy dependent interaction probabilities for H2O vapor molecules is potentially less and more challenging than for the case of neutron-nuclei interactions. “Less” cchallenging because the energy span that must be considered is much narrower, it extends approximately from 0 to 0.2 eV. For neutrons kinetic energies up to 10E7 eV must be considered. “More” 61

challenging, because the H2O molecule is a considerably more complex entity then a neutron. It is a vibrating, rotating dipole, because of this the polarity of two approaching molecules fluctuates. Their interaction and the consequences of such an interaction is a product of not one but several probability distributions. The relationship between the relative speed of two colliding particles and their interaction probability is more complex. The time that two approaching molecule spend in their mutual viscinity can differ by orders of magnitude, this “mutual influence time” significantly influences their interaction probability. For example, consider the kinetic energy segmentation shown in Flg 12: a “head –on” collision of two group 4 molecules can have relative speeds on the order of 3000 m/s, on the other hand, a faster molecule overtaking a slower one in group 1 can have a relative speed as low as ~3 m/s. In the first case the time interval during which their polarity reverses (due to rotation) is simmilar to the time of approach, in the other many fluctuations can occur while they are at a distance in which their force fields interact. Basic characteristics of the H2O molecule relevant for the evaluation of interaction probabilities are shown in table 4. The table lists geometric dimensional characteristics available in literature. The shortest dimension is the O to H atom distance, the largest is the distance at which the Lennard Jones potential used to model the H2O molecule for Molecular Dynamics computations decreased to 10% of its maximum value1. As a first apporximation it can be assumed that the “interaction cross section” of two colliding molecules is proportional to (2*R)2. This indicates that on the basis of the listed dimensions, the interaction probability of the molecules approaching each other at the lowest velocity is ~50 times higher than for those which have such a large relative speed that their kinetic energy exceeds the force fields extending past the molecular boundaries. This issue requires a thorough fundamental analysis. The results shown in Table 4 are evaluated using concepts and notation employed in neutronic computations and should be considered as preliminary estimates. 1  The last two lines iin table 4 are based on the Lennard-Jones potential employed in the TIP3P(23) H2O molecule model used in MD simulations. At ~0.31 n-m this potential has an attractive force maximum, at 0.61 n-m it decreases by ~90%. The strength of the Lenn-J potential at R = 0.61 n-m is ~0.003 eV, this matches the kinetic energy of a 200 m/s vapor molecule in velocity „group 1“ (Fig 11). In comparison, the kinetic energy of at mid point velocity for molecules of „group 4“ is 0.18 eV, such a molecule will not be affected by a residual potential of 0.003 eV.

62

In neutronic computations neutron-atom interaction probabilities are expressed in therms of “microscopic cross sections”. σ(E) represents the effective cross section of individual atoms, ∑(E) is the “macroscopic cross section” which is obtained by multiplying σ(E) by the respective atomic density (the result is usually expressed in 1/cm units). The volumetric neutron-atom interaction rate is obtained by multiplying ∑ by the “neutron flux”. φ (φ = neutron density * av. neutron speed). The evaluation of the mean free path λ takes into account the circumstance that for H2O vapor molecules all collisions are collisions between entities having the same mass. The adaptation of these definitions to vapor molecules requires a redefinition of units for σ(E). Nuclei have dimensions on the order of 10E-14 m, to avoid large negative exponents the neutronic microscopic cross-sections are expressed in “barn”. units. A “barn” has an area of 10E-28 m². Molecules have dimensions on the order of nano-m, analogous molecular unit should therefore be a on the order of a squared nano-m, (10E-18 m²). In the table this area unit is abbreviated as “sqn”2. Table 4 provides a first estimate overview of the differences that can be expected for interaction of molecules in different kinetic energy segments. A relevant issue for the interpretation of experimental data and analysis is that the mean free path of molecules is also strongly dependent on molecular velocity. The application of the 3*λ rule for the definition of the transition region near an anisotropic interface should consider the mean free path of the fastest molecules.

B-2 Finite difference issues The adaptation of finite difference techniques for the evaluation of vapor molecule distributions presents several challenges. Besides the noted similarities between neutronic and vapor molecule interactions, there are also significant differencesl. The major ones are: 1) Interactions between neutrons are so rare, that neutron-neutron scattering events can be neglected, consequently with respect to neutron density the mass/ 2  The name of the “barn” unit shows that the pioneering neutron physicist did not lack a sense of humor. The origin of the name comes from the phrase “You could not hit a barn door”. An imaginative physicist is needed to provide an analogous catchy name for the “square nano meter”

63

energy balance equations for neutrons are first order. For H2O vapor molecules, “self-scattering” is the dominant, in fact, aside of collisions with walls or interfaces, the only interaction, consequently, in the mass/energy balances the density of the H2O molecules appears as a square term, the equations are quadratic. 2) The interaction probability for both neutrons and molecules depends on the relative speed between the interacting particles. It is easier to apply this characteristic in neutronic computations because in the “lab” system one of the particles (the atom) is stationary. For vapor both particles are moving. Differences 1 and 2 complicate the solution of vapor mass/energy balance equations, however there are also differences which offer advantages. 3) The pressure and temperature of a neutron population can be defined, however in nuclear reactor cores or their vicinity it is so small that it can not be measured. For vapor molecules pressure and temperature are readily measured, this information can be used in the solution of vapor balance equations. 4) An important advantage is that for equilibrium conditions the solution of the vapor molecule mass/energy balances is provided by the Maxwellian distribution. The corollary to this observation is that in most cases for vapors equilibrium conditions prevail. The last point is probably the principal reason why finite difference techniques have not been developed for mass/balance caluclations of vapor molecules. Such a computational capability is required only within the transition region near boundaries which can not be characterized by a Maxwellian distribution. The ~3*λ wide region is usually relatively small, it becomes important when the boundary is an evaporating or condensing interface. Fig B1 shows the finite difference segmentation of a one-dimensional region near such an interface. The lower boundary is an evaporating interface, equilibrium conditions prevail at the upper boundary. The upward dimension (Z) is divided into N spatial segments represented by the index i, the kinetic energy of the molecules is divided into M energy “groups” designated by the index j. The molecular density within each segment is represented by Ngi*Γi ,j where Ng is

64

the total molecular density and Γi ,j is the fraction of the vapor molecules within the energy group j at vertical space segment i. If measurements of local pressure and temperature are available, Ng is determined by Eqs. B-1 to B-3 and does not have to appear in the energy/mass balance equations. B-1

B-2

B-3 Where Urms(i,j), is the “root-mean-square” velocity of the molecules in space/energy node (i,j). Within the transition region it is not known, as the number of energy group M increases, it can be approximated by the average velocity assuming the velocity distribution is linear.

Fig B-1 One dimensional space-energy matrix for vapor molecules in the transition region near an evaporating interface 65

Note that the segmentation shown in Fig B-1 does not exploit all of the capabilities of the finite difference method. The directional variable is not segmented, this implies that the distribution of molecular velocities is assumed to be isotropic. In neutronic computations this assumption separates “diffusion theory” from “transport theory” methods. The introduction of a directional variable leads to considerably more complicated algorithms. The mean free path of neutrons is orders of magnitude larger than for molecules, correspondingly the anisotropy of velocity distributions is greater. For neutronic computations the additional complexity is justified, it is probably not justified for vapor molecules except for cases when near vacuum conditions exist. In any case, adaptation of finite difference methods for vapor molecules should start with the simpler “diffusion theory” alternative. Consider the energy and mass balance of an internal node in the transition region represented by Γi, j. At steady state conditions the molecular density in node (i,j) is determined by the rate at which molecules are moving across the segmented spatial and energy interfaces. If a bulk velocity UB is present, then the in-out flow balance of molecules per second due to bulk motion can be represented by balance term A: A. Bulk flow component:



B-4

In addition to bulk flow, molecules in each energy group can be transfered by diffusion. In fact, since molecular velocities are much larger then UB, whenever group dependent concentration gradients exist, the diffusion term B can dominate. B. Diffusive transfer component: B-5 Where Vj and Dj are energy group j dependent velocity and diffusion coefficient. The most complicated balance Eq. term C represents the transfer of molecules accross the energy group interfaces ∆E:

66

C. Inter-group scat. component: B-6 Where ζm – j is the probability that a molecule from group m scatters into group j. Term C states that the change of molecular concentration in node (i,j) is equal to the sum of the molecules scattered out of node (i,j) minus the number of molecules scattered into energy group j form all the other energy groups. In the 3* λ transition region for steady state conditions: B-7 In the region sufficiently far from the interface an equilibrium condition is reached, the Maxwellian distribution can then be used to evaluate Ngi*Γi,j the molecular density in each energy group. In this region no net transfer between energy groups occurs, therefore the B and C terms are equal to zero. If the pressure and temperature does not change with Z, then this is true also for term A. Of special interesst is that C, the intergroup scattering term not only sums to zero, but in fact no net interchange rate of molecules between any two energy groups occurs. That is: B-8 Since the Γi ,j values are known, the inter-group scattering probabilities ζm – j can be evaluated by an itterative procedure. The generation of this matrix requires significant development effort. This can be justified by the circumstance that for the computation of phase change rates the use of continuous functions has reached its applicability limit, on the other hand computer simulation of phase change by methods like MD (Molecular Dynamics) is still not practical. Even after MD capabilities have increased by orders of magnitude, the method will require large computer resources and will be able to simulate only very small volumes over short time frames. A more generally applicable computational tool is required, Adaptation of the finite difference method can fill this need.

67

Contents 1. Introduction

5

2. CV-1 Schematic

7

3. CV-2 Schematic

9

4. CV-3 Schematic

11

5. Theory vs. experiment

13

6. CV-4 Schematic

16

7. Molecular characteristics of water at low pressures

18

8. Agitational and total energy

20

9. The hydrogen bond and its impact on the phase change of water

23

10. The limits of kinetic theory

27

11. Alternate computational approaches

30

12. Insights from neutron transport theory

32

14. Experimental verification

39

16. Inferred spectrum of evaporating molecules

41

17. Inferred spectrum of condensing molecules

42

18. Conclusions, proposals and recommendations

47

Appendix A

53

Appendix B

61

69

Al-146 Almenas, Kazys Evaporation/condensation of water. Unresolved issues. I. Phase change at low pressures, laminar conditions. Monograph / Kazys Almenas. Kaunas: Vytauto Didžiojo universitetas; Vilnius: Versus aureus, 2015. – 72 p. ISBN 978-609-467-098-5 (Print) ISBN 978-9955-34-531-2 (Print) ISBN 978-609-467-097-8 (Online) ISBN 978-9955-34-530-5 (Online) This study assesses empirical data and theoretical approaches relevant for condensation/evaporation processes at low vapor pressures. Classical theoretical and experimental studies based on kinetic theory are augmented by information from the fields of thermodynamics, physics, chemistry and potentially relevant computational approaches developed in neutron transport theory and molecular dynamics. In the course of the overview a number of conclusions are reached, inferred distributions of evaporating and condensing molecular spectra are proposed and recommendation regarding further development of analytical approaches are presented. UDK 556

Kazys Almenas Evaporation/condensation of water. Unresolved issues I. Phase change at low pressures, laminar conditions

Designer Kornelija BUOŽYTĖ 2014 12 30. Issuance by Order No. K15-003. Published by: Vytautas Magnus University K. Donelaičio g. 58, LT-44248 Kaunas www.vdu.lt | [email protected] “Versus aureus” Publishers Rūdninkų g. 10, LT-01135 Vilnius www.versus.lt | [email protected] Printed by: „BALTO print“ Utenos g. 41A, LT-08217 Vilnius www.baltoprint.lt