- Preface
- Introduction to Transport Phenomena
- Thermodynamics of Multiphase Systems
- Generalized Governing Equations in Mutliphase Systems: Local Instance Formulations
- Generalized Governing Equations for Multiphase Systems : Averaging Formulations
- Solid-Liquid-Vapor Phenomena and Interfacial Heat and Mass Transfer
- Melting and Solidification
- Sublimation and Vapor Deposition
- Condensation
- Evaporation
- Boiling
- Two-Phase Flow and Heat Filter
- Appendix A : Constants, Units and Conversion Factors
- Appendix B: Transport Properties
- Appendix C: Vectors and Tensors
- Index

9
EVAPORATION
9.1 Introduction
Evaporation processes generally occur from liquid films, drops, and jets. Films may flow on a heated or adiabatic surface as a result of gravity or vapor shear. Drops may evaporate from a heated substrate, or they may be suspended in a gas mixture or immiscible fluid. Jets may be cylindrical in shape or elongated (ribbon-like). In film evaporators, a liquid film is evaporated in order to: (1) cool the surface on which it flows, (2) cool the liquid itself, or (3) increase the concentration of some component(s) in the liquid. Heat may be transferred to the film surface by conduction or convection. The design of film evaporators is most frequently constrained by the extent of the liquid superheat, since nucleate boiling must be prevented. Nucleate boiling can lead to product deterioration because elevated temperatures increase chemical reaction rates. Also, the accumulation of deposits on walls (fouling) as a result of nucleate boiling hinders heat transfer. To avoid nucleate boiling, liquid superheating should not exceed 3 to 40 K, depending on the liquid. This chapter deals with evaporation with no nucleate boiling. In addition to the capability for low liquid superheating, film evaporators may be designed to achieve short residence time (that is, contact time between the liquid and the evaporator wall), which further reduces the tendency for chemical reactions to occur. The tendency for fouling is also reduced in more rapidly flowing films. These features are ideal for heat-sensitive fluids such as those in food processing and polymer devolatilization (Alhusseini et al., 1998). A thin liquid film in an evaporator can be produced by several arrangements, two of which are falling and climbing films (See Fig. 9.1). In falling film evaporators, the liquid is evenly distributed around the opening of a vertical channel – usually a tube – and falls under gravitational force as its surface evaporates to a vapor or gas mixture. If the vapor is drawn to the bottom, the flow is cocurrent as shown in Fig. 9.1(a) and the shear stress of the vapor on the liquid film increases its speed downward. When the vapor is drawn to the top, the flow is countercurrent as shown in Fig. 9.1(b). A climbing film evaporator Fig. 9.1(c) is designed so that liquid is fed into the bottom, and vapor shear pushes the liquid film along the wall to the top. In all of these cases, the lack of a hydrostatic
678
Transport Phenomena in Multiphase Systems
Figure 9.1 Examples of film evaporators.
pressure drop eliminates the corresponding temperature drop along the evaporator length, thereby increasing the uniformity in temperature. A more complicated example, shown in Fig. 9.1 (d), is a two-phase closed thermosyphon where evaporation takes place in the evaporator section and liquid and vapor flow in opposite directions. The residence time of a liquid in a falling film evaporator may be increased with a countercurrent configuration. If this configuration is chosen, however, the likelihood of entrainment, flooding, or liquid blocking should be considered. Unlike the cocurrent configuration, countercurrent flow is more prone to entrainment, which occurs when liquid is carried upward by the vapor (see Section 8.3.8). Liquid blocking – when liquid bridges the channel and blocks the vapor path – can also be more likely in a countercurrent configuration. In addition to the film evaporation described above, evaporations from liquid droplets attached to a heated wall or surrounded by hot gas can also find their application in irrigation of crops, firefighting, and combustion. Figure 9.2 shows evaporation from a liquid droplet attached to a heated wall; its application can be found in surface spray cooling, where liquid droplets are sprayed onto a hot surface. Heat is conducted through the droplet to the interface, where an abrupt temperature drop takes place due to evaporation. The mass fraction of the vapor component in the gas mixture, ωv , is greatest next to the interface, and by diffusion the mass fraction decreases to its bulk level with increasing distance from the wall. The size of the liquid drop and the temperature distribution at two different times ( t2 > t1 ) are shown in Fig. 9.2. As time goes on, the liquid droplet becomes smaller (as indicated by the dashed-line), while the temperature at the heated wall and interface remain unchanged.
Chapter 9 Evaporation
679
r rI
v
v,I v,
Tw t2 > t1 t2 TI T rI
Figure 9.2 Evaporation from a liquid droplet on a heated wall.
t1
When a liquid droplet is surrounded by a hot gas mixture, evaporation takes place on the surface of the droplet, as shown in Figure 9.3. Such processes are used in bulk cooling as well as spray fuel injection in diesel engines. The time required for such a droplet to completely evaporate as a member of a dispersed phase will determine the overall heat transfer (cooling effect). The partial pressure, pv, and the mass fraction, v, of the vapor in the gas mixture are the greatest adjacent to the surface of the droplet, and diffusion causes the partial pressure and mass fraction to decrease with increasing distance from the interface. The temperature profile in the drop is dominated by transient conduction, but convection takes place outside the drop. Since the evaporating liquid is directly supplied by the liquid droplet, evaporation on the drop surface is dominated by heat transfer, not by mass diffusion. The criteria for evaporation in pure vapor or gas mixtures are discussed in Section 9.2. Section 9.3 includes a discussion of evaporation from a liquid film on a horizontal adiabatic wall, as well as evaporation from a vertical falling film with waves on an adiabatic surface. The contributions of the waves are considered by introducing additional corrections on viscosity, thermal conductivity, and mass diffusivity. Evaporation from a vertical falling film on a heated surface is presented in Section 9.4; this includes classical Nusselt evaporation, i.e., evaporation from a laminar film with a smooth surface;
680
Transport Phenomena in Multiphase Systems
r rI
v
v,I v,
t2 > t1 T T t2 TI t1
rI
Figure 9.3 Evaporation from a liquid droplet suspended in vapor-gas mixture.
and laminar film with waves; as well as turbulent film. Section 9.5 discusses direct contact evaporation from liquid drops surrounded by hot gas and hot immiscible liquid, as well as direct evaporation of a liquid jet surrounded by hot gas. Evaporation in porous media and microchannels are discussed in Sections 9.6 and 9.7, respectively.
9.2 Classification and Criteria of Evaporation
The latent heat of vaporization acts as a heat sink during phase change. It is supplied either by migration through the liquid, as in heterogeneous evaporation, or directly to the interface, as in direct contact evaporation. In the former, the heat must migrate by conduction (and in some cases convection) through the liquid to the interface. Figure 9.4 illustrates the difference between heterogeneous and direct contact evaporation. When a beaker of water is placed on a hot plate, as shown in Fig. 9.4(a), only heterogeneous evaporation can occur, because heat must migrate from the bottom of the beaker to the surface of the water in order for evaporation to take place. The liquid’s temperature is highest at its point of contact with the bottom of the beaker. Adjacent to the bottom is a layer through
Chapter 9 Evaporation
681
Figure 9.4 Comparison of (a) heterogeneous and (b) direct contact evaporation.
which heat passes only by conduction. Convection takes place throughout most of the depth, causing the temperature to decrease slightly. Just below the surface, a thin layer exists in which the temperature drops abruptly as the interface is approached. The surface temperature of the liquid can be assumed to equal the saturation temperature corresponding to the partial pressure of the vapor at the surface. Evaporation from a liquid film or a droplet attached to a heated wall is heterogeneous evaporation. Figure 9.4(b) shows an example of direct contact evaporation, where a beaker of water sits on a table in a room of temperature T∞ with a relative humidity < 100%. Since the heat required for the evaporation is supplied at the surface, there is no temperature gradient through the depth of the liquid except just below the interface. Evaporation of a liquid droplet suspended in a hot gas mixture is another example of direct contact evaporation, because the latent heat of evaporation is supplied to the surface of the droplet by the surrounding hot gas mixture. Evaporation of a liquid film or droplet attached to an adiabatic wall is
682
Transport Phenomena in Multiphase Systems
direct evaporation because the latent heat is provided by the gas mixture, not the wall. Evaporation from a liquid surface is a mass transfer process, where the rate at which the molecules leave the liquid exceeds the rate at which the molecules move back into the liquid from the vapor. Therefore, the vapor pressure at the interface, which is the saturation pressure corresponding to the liquid surface temperature, psat (TI ) – if the effects of capillary and disjoining pressures are negligible – must exceed the vapor pressure, pv. The difference between the saturation pressure and the vapor pressure is negligibly small when a liquid is evaporating to its pure vapor. Consequently, the temperature at the interface is equal to the saturation temperature corresponding to the local vapor pressure. If evaporation is from a liquid to a mixture of vapor and noncondensable gas, the vapor generated at the interface must be diffused to the main stream of the mixture. In this case, evaporation will be further limited by the mass diffusion process. However, if the effect of vapor diffusion is negligible, the temperature at the interface equals the saturation temperature corresponding to the partial pressure of the vapor in the mixture. In this case, evaporation can occur only if saturation pressure corresponding to the interfacial temperature exceeds the partial pressure in the vapor. For example, a beaker of water is in thermal equilibrium with its surroundings at 20 °C and 50% relative humidity. The saturation pressure of water at 20 °C is 2.3 kPa, and the vapor partial pressure in air at 50% humidity is half of that, 1.15 kPa. Since the liquid saturation pressure exceeds the vapor partial pressure, evaporation occurs. Therefore, heterogeneous evaporation requires that the liquid is superheated to a temperature above the saturation temperature corresponding to the partial pressure of the vapor, TA > Tsat ( pv ) . For direct contact evaporation, the temperature of the vapor-gas mixture in contact with the liquid must exceed the saturation temperature corresponding to the partial pressure of the vapor. Let us reuse the example of a beaker of water in thermal equilibrium with its surroundings at 20 °C and 50% relative humidity. Since the gas temperature of 20 °C is greater than the saturation temperature of 8.8 °C, corresponding to the partial pressure of the vapor, evaporation occurs. The conditions under which evaporation to a mixture of vapor and noncondensable gas can occur are demonstrated in a T-s diagram, as shown in Fig. 9.5. The saturated liquid and vapor lines represent the temperature and entropy conditions necessary for equilibrium between liquid and vapor phases. The saturated vapor line and the isotherm T = TI intersect at a point I from which the interface saturation pressure isobar p = pI is drawn. For evaporation to occur, the interface saturation pressure pI must be greater than the vapor partial pressure pv. If the temperature of the vapor is above the interfacial temperature (see point a), the latent heat of evaporation is supplied by the gas; evaporation under this condition is referred to as hot evaporation. The liquid for hot evaporation can be either saturated or subcooled, as represented by curve a′b′ in Fig. 9.5. The above example with the beaker of water, shown in Fig. 9.4 (b), occurs at point a because the vapor temperature (20 °C) is greater than the interfacial temperature,
Chapter 9 Evaporation
683
s
Figure 9.5 Evaporation domains for a vapor-gas mixture.
which is at the saturation temperature corresponding to the partial pressure of the vapor (8.8 °C). At point b, which represents a case where the relative humidity is 100%, the gas temperature equals the saturation temperature and hot evaporation ceases. At point c, the temperature of the gas is below the saturation temperature, so the gas cannot provide the latent heat of vaporization. For evaporation occur at point c, heat must be supplied by a superheated liquid, represented by curve b′f ′ in Fig. 9.5. As discussed in Chapter 2, the superheated liquid is in a metastable state. If the temperature of the liquid is higher than that of point f ′ , vapor bubbles will appear in the liquid. Evaporation occurring under this condition is referred to as cool evaporation because the interface is warmer than the gas. At the intersection of isobar p = pv and the saturated liquid line (point d) the vapor is saturated. As a result, cool evaporation does not occur under equilibrium condition. Practically speaking, evaporation is still possible, but it will result in supersaturated vapor, which is in a metastable state as discussed in Chapter 2. If the vapor-gas mixture is already supersaturated, as represented by point e, evaporation is possible provided that sufficient superheat exists in the liquid. Evaporation to the vapor-gas mixture at point f will be impossible, because it represents the maximum possible supersaturation of the vapor-gas mixture.
684
Transport Phenomena in Multiphase Systems
9.3 Evaporation from an Adiabatic Wall
9.3.1 Evaporation from Horizontal Films
Evaporation from thin liquid films occurs in many industrial and natural processes, including drying, evaporative cooling, and sweating. Figure 9.6 shows the physical model of horizontal thin film evaporation under consideration (Carey, 1992). A film wets a horizontal surface over which flows a gas of ambient temperature T∞ and mass fraction of water ω∞ at a velocity of u∞ . Since the solid surface underneath the liquid is adiabatic, the latent heat of vaporization is provided by the vapor-gas mixture flowing above the liquid film. The heat source (vapor-gas mixture) and the evaporating liquid are in direct contact, making this scenario an example of direct contact evaporation. The liquid is evaporated and the latent heat is absorbed in the evaporation process. The resulting vapor injects into the boundary layer and is removed by the gas flow. The boundary layer becomes thicker, and a free stream of the gas flow is displaced from the surface being cooled. While phase change is the dominant mechanism of heat transfer at lower gas temperatures, vapor injection becomes more important at higher gas temperatures that correspond to higher evaporation rates. Due to diffusion and convection, the temperature rises and the concentration falls from their values at the film surface to their ambient values at the edge of the boundary layer. A nonslip boundary condition exists at the film surface, and the film itself is considered stationary with respect to the gas. Assuming a steady state, constant density, and incompressible flow, the continuity equation is
∂u ∂v + =0 ∂x ∂y
(9.1)
With no pressure gradients and constant viscosity assumed, the boundary layer momentum equation is written as ∂u ∂u ∂ 2u u +v =ν 2 (9.2) ∂x ∂y ∂y
u, T, y, v
m t
u u(x,y) T
T T(x,y) x, u
Liquid film
Figure 9.6 Evaporation from a thin liquid film on a horizontal surface.
Chapter 9 Evaporation
685
where ~ on top of viscosity signifies mass-averaged properties of the mixture. Neglecting viscous dissipation and assuming constant thermal diffusivity and specific heats, the boundary layer energy balance is ∂T ∂T ∂ 2T u +v =α 2 (9.3) ∂x ∂y ∂y where the axial conduction on the x-direction has been neglected because heat transfer occurred mainly in the y-direction. The mass fraction of water is also accounted for by assuming a constant mass diffusivity: ∂ω ∂ω ∂ 2ω +v =D 2 u (9.4) ∂x ∂y ∂y As y → ∞ the boundary conditions can be taken directly from Fig. 9.6:
u → u ∞ , T → T∞ , ω → ω∞ (9.5) The boundary conditions are now laid out at the film surface y = 0. There is a nonslip condition, and the mass flux m′′ evaporates normal to the surface. m′′ u = 0 , v = , ω = ωδ , T = Tδ = Tsat (ωδ , p∞ ) (9.6)
ρ
An energy balance at the interface accounts for the heat of vaporization and mass flux of the evaporating fluid. With the nonslip boundary condition, the only mode of heat transfer is conduction: § ∂T · m′′hAv = −q′′ = k ¨ (9.7) ¸ © ∂y ¹ y = 0 Substituting eq. (9.7) into eq. (9.6), the velocity component in the ydirection becomes k § ∂T · v y =0 = (9.8) ¨ ¸ ρ hAv © ∂y ¹ y =0 The mass flux at the interface is the result of both diffusion and convection. Therefore, the mass balance at the interface is written as § ∂ω · m′′ = − ρ D ¨ (9.9) ¸ + ρωδ v y = 0 © ∂y ¹ y = 0 Substituting eq. (9.7) into eq. (9.9), conservation of mass at the interface becomes § ∂ω · k § ∂T · −ρ D ¨ (9.10) ¸ + ρωδ v y = 0 = ¨ ¸ hAv © ∂y ¹ y = 0 © ∂y ¹ y = 0 Equations (9.1) – (9.2) are the same as the laminar boundary layer equations for forced convection over a flat plate with blowing on the liquid surface. By defining the stream function ψ as follows, ∂ψ ∂ψ u= v=− (9.11) ∂y ∂x
686
Transport Phenomena in Multiphase Systems
the continuity equation (9.1) is automatically satisfied. The momentum equation in terms of the stream function is then ∂ψ ∂ 2ψ ∂ψ ∂ 2ψ ∂ 3ψ − =ν 3 (9.12) ∂y ∂x∂y ∂x ∂y 2 ∂y Introducing the following similarity variable: T − Tδ ω − ωδ u ψ , θ= , ϕ= (9.13) η=y ∞, f = T∞ − Tδ νx ω∞ − ωδ ν u∞ x eq. (9.12) and eqs. (9.3) – (9.4) can be reduced to a set of ordinary differential equations: 1 f ′′′ + ff ′′ = 0 (9.14) 2 1 θ ′′ + Pr f θ ′ = 0 (9.15) 2 1 (9.16) ϕ ′′ + Scf ϕ ′ = 0 2 The velocity components in the x- and y-directions in terms of the dimensionless stream function f can be obtained by eqs. (9.11) and (9.13), i.e., 1 ν u∞ (9.17) u = u∞ f ′(η ), v = (η f ′ − f ) 2 x The boundary conditions at y→∞ represented by eq. (9.5) become f '(∞) = 1, θ (∞) = 1, ϕ (∞) = 1 (9.18) The boundary conditions represented by eq. (9.6), except the velocity component in the y-direction, can be rewritten to f '(0) = 0, θ (0) = 0, ϕ (0) = 0 (9.19) Conservation of mass and energy at the interface, eqs. (9.8) and (9.10), can be rewritten in terms of the similarity variables, i.e., 2Ja f (0) = − θ ′(0) (9.20) Pr ωδ − ω∞ Sc (9.21) φ ′(0) = Jaθ ′(0) 1 − ωδ Pr where Ja is the Jakob number and is defined as c p [T∞ − Tsat (ωδ ) ] (9.22) Ja = hAv and Sc is the Schmidt number defined as
(9.23) D The evaporation problem is now described by a set of ordinary differential equations (9.14) – (9.16) subjected to boundary conditions specified by eqs. (9.18) – (9.21). Since eq. (9.14) is a third-order ordinary differential equation,
Sc =
ν
Chapter 9 Evaporation
687
it requires three boundary conditions: f '(∞) = 1 in eq. (9.18), f '(0) = 0 in eq. (9.19), and eq. (9.20). Equations (9.15) and (9.16) are both second-order ordinary differential equations, each of them requiring two boundary conditions, which are specified in eqs. (9.18) and (9.19). Therefore, the problem is mathematically defined with boundary conditions stated by eqs. (9.18) – (9.20), which makes the boundary condition specified in eq. (9.21) an extra boundary condition – and makes the problem overstated. This happened because the heat and mass transfer are not independent of each other and the energy balance at the interface requires that eq. (9.21) be satisfied. The boundary value problem can be solved using the Runge-Kutta method in conjunction with a shooting method. The solution procedure begins with an assumed I and uses boundary conditions specified by eqs. (9.18) – (9.20) to solve eqs. (9.14) – (9.16). Once a solution is obtained, eq. (9.21) is employed to find ωδ . If ωδ obtained from eq. (9.21) agrees with the assumed value, the solution is complete. Otherwise, the assumed value of ωδ is corrected and the solution procedure is repeated until a converged solution is obtained. Once the converged solution is obtained, the local heat transfer coefficient can be evaluated by u∞ k § ∂T · θ ′(0) (9.24) hx = ¨ ¸ =k νx T∞ − Tδ © ∂y ¹ y = 0 The local Nusselt number is then 1 hx 2 Nu x = x = θ ′(0) Re x k The local mass transfer coefficient can be evaluated by u∞ D § ∂ω · hmx = ϕ ′(0) ¨ ¸ =D ω∞ − ωδ © ∂y ¹ y =0 νx
and the local Sherwood number is 1 hx 2 Shx = m = ϕ ′(0) Re x D where (9.25)
(9.26)
(9.27)
Re x =
u∞ x
ν
(9.28)
Although eqs. (9.25) and (9.27) have the same form as the case of singlephase forced convective heat or mass transfer over a flat plate, the solutions of θ ′(0) and ϕ ′(0) depend on Pr, Sc, Ja, and ∞. Consequently, simple empirical correlations similar to the case of single-phase forced convective heat or mass transfer over a flat plate are very difficult to obtain. The evaporation problem presented here is a coupled problem because the solutions of the momentum, energy, and species equations are coupled with eqs. (9.20) – (9.21). Since the blowing velocity at the surface of the liquid film, determined by eq. (9.8), is usually much lower than the incoming vapor-gas
688
Transport Phenomena in Multiphase Systems
mixture velocity, u∞ , one can assume that the blowing velocity at the liquid surface is negligible. Then eq. (9.20) can be replaced with f (0) = 0 (9.29) The heat and mass transfer problem described by eqs. (9.14) – (9.16), (9.18) – (9.19) and (9.29) then cease to be a conjugated problem, because each equation in eq. (9.14) – (9.16) can be solved independently. Numerical solutions of these ordinary differential equations yield the following results (Kays et al., 2004):
θ ′(0) = 0.332 Pr 3 ϕ ′(0) = 0.332Sc
1 3
1
(9.30)
(9.31) Substituting eqs. (9.30) and (9.31) into eqs. (9.25) and (9.27) yields 1 1 hx x 2 Nu x = = 0.332 Re x Pr 3 (9.32) k 1 1 hx 2 Shx = m = 0.332 Re x Sc 3 (9.33) k The average heat and mass transfer characteristics are often relevant for practical applications. The average Nusselt and Sherwood numbers based on average heat and mass transfer coefficients for a liquid film with a length of L are obtained by 1 1 hL 2 Nu = = 0.664 Re L Pr 3 (9.34) k 1 1 hm L 2 Sh = = 0.664 Re L Sc 3 (9.35) k where uL Re L = ∞ (9.36)
ν
The mass fraction of water at the liquid surface can be obtained by substituting eqs. (9.30) and (9.31) into eq. (9.21), i.e.,
ωδ =
ω∞ + (Sc / Pr)2 /3 Ja
1 + (Sc / Pr)2/3 Ja
(9.37)
The applications of the above empirical correlations can be demonstrated by the following example.
Example 9.1 Air at 1 atm with 30 ÛC and 55% humidity flows over a 0.8 m-long wet flat plate with a velocity of 1 m/s. Estimate the rates of heat transfer to and evaporation from the wet flat plate. Estimate the evaporationinduced blowing velocity on the surface and compare it with the incoming velocity of the humid air.
Chapter 9 Evaporation
689
Solution: Assuming the binary (air-water) mixture can be modeled as an ideal gas, the total density of the mixture is p M ( p − p1 ) M 2 ρ = ρ1 + ρ 2 = 1 1 + (9.38) RuT RuT where M1 and M2 are the molecular masses of water and air, respectively, which can be found from Table 1.1 as M 1 = 18.015kg/kmol and M 2 = 28.97kg/kmol . The mass fraction of the water in the mixture is then ρ p1 M 1 ω= 1 = (9.39) ρ p1 M 1 + ( p − p1 ) M 2 The saturation pressure corresponding to the air temperature of T∞ = 30 o C is p1, sat = 4.246 kPa . Since the relative humidity in the binary mixture is 55%, the partial pressure of the water in the mixture is p1 = 0.55 p1, sat = 2.335 kPa. Therefore, the mass fraction of water at the
locations far from the flat plate is 2.335 × 18.015 = 0.0144 ω∞ = 2.335 × 18.015 + (101.325 − 2.335) × 28.97 The properties of the moist air can be evaluated using those for pure air in Table B.1, because the mass fraction of the water is very low, i.e., ρ = 1.165 kg/m3 , c p = 1.006 kJ/kg-K, ν =0.16 × 10-4 m 2 /s,
k = 0.026W/m-K, α = 2.23 × 10−3 m 2 /s and Pr = 0.72 . The mass diffusivity of water in air can be found from Table B.71 as D12 =
0.26 × 10−4 m 2 /s , so the Schmidt number is Sc = ν / D = 0.62 . The surface temperature of the flat plate, T , is the saturation temperature corresponding to the partial pressure of the water. If the partial pressure at the surface can be treated as unchanged (an approximation which will be corrected later), the surface temperature is the saturation temperature corresponding to p1 = 2.335 kPa, i.e., Tδ = 20 o C . The latent heat of evaporation for water at 20 °C is hAv = 2463.8 kJ/kg. The Jakob number is
= 0.00408 2463.8 hAv The mass fraction of water at the surface of the flat plate must satisfy eq. (9.37), i.e., ω + (Sc / Pr)2 /3 Ja 0.0144 + (0.62 / 0.72)2/3 × 0.00408 ωδ = ∞ = = 0.0180 1 + (Sc / Pr)2/ 3 Ja 1 + (0.62 / 0.72)2/3 × 0.00408 which is based upon the assumption that the partial pressure of water is the same throughout the boundary layer. With the mass fraction obtained Ja = c p (T∞ − Tδ ) = 1.006 × ( 30 − 20 )
690
Transport Phenomena in Multiphase Systems
above, the partial pressure of the water at the surface can be corrected using eq. (9.39), i.e., ωδ M 2 p 0.0180 × 28.97 × 101.325 p1,δ = = = 2.901 kPa M1 − ωδ ( M1 − M 2 ) 18.015 − 0.018 × (18.015 − 28.97) The saturation temperature corresponding to the above partial pressure is Tδ = 23.49 o C . The mass fraction with this new surface temperature is found to be ωδ = 0.0168 . The Reynolds number is uL 1 × 0.8 Re L = ∞ = = 5 × 104 −4 ν 0.16 × 10 The average Nusselt and Sherwood numbers can be obtained with eqs. (9.34) – (9.35):
1 1 1 1 2 Nu = 0.664 Re L Pr 3 = 0.664 × (5 × 104 ) 2 × 0.72 3 = 133.08 1 1 1 1 2 Sh = 0.664 Re L Sc 3 = 0.664 × (5 × 104 ) 2 × 0.62 3 = 126.60 The average heat and mass transfer coefficients can be found as Nuk 133.08 × 0.026 = = 4.33 W/m 2 K h= 0.8 L ShD 126.6 × 0.26 × 10−4 hm = = = 0.00411 m/s 0.8 L The rates of heat transfer and evaporation are q′ = hL(T∞ − Tδ ) = 4.33 × 0.8 × (30 − 23.49) = 22.55 W/m
m′ = hm L ρ (ω∞ − ωδ )
= 0.00411 × 0.8 × 1.165 × (0.0180 − 0.0144) = 1.379 × 10−5 kg/s-m The average blowing velocity due to evaporation is m′ 1.379 × 10−5 = = 1.480 × 10−5 m/s vδ = ρ L 1.165 × 0.8 which is very low, compared with the velocity in the x-direction.
9.3.2 Evaporation from a Vertical Falling Film
In falling film evaporators, the flowing gas provides the latent heat of vaporization, and the gas cools as it passes through the evaporator. If the gas flows downward, it grows denser as it cools. The buoyancy force creates a chimney effect, by means of which the gas exits the bottom and more gas enters the top. An experimental study of this flow was conducted by Yan et al. (1991) and followed by a numerical simulation by Yan and Lin (1991), which is presented below. Figure 9.7 shows the physical model of evaporation from a vertical falling film on two parallel walls with a spacing of 2b. While the outer
Chapter 9 Evaporation
691
sides of the wall are insulated, the liquid falling films on the inner walls evaporate into a gas that flows through the channel. Assuming there is no pressure gradient and the inertial term is negligible, the momentum balance in the assumed laminar liquid film with negligible inertia is ∂uA º ∂ª ′ (9.40) «( μA + μA ) ∂y » + ρA g = 0 ∂y ¬ ¼ ′ where μA is a correction to the viscosity in response to presumed surface waves on the film. Considering conduction through the liquid and sensible heat loss, the energy balance is ∂T ∂T º ∂ª ρA c p ,A uA A = «( kA + kA′ ) A » (9.41) ∂x ∂y ¬ ∂y ¼
Figure 9.7 Evaporation of a thin liquid film on a vertical adiabatic surface.
692
Transport Phenomena in Multiphase Systems
where kA′ is a correction to the conductivity to account for the effect of the surface waves. Natural convection occurs in the gas region, and the continuity equation is ∂ ∂ (9.42) ( ρu )g + ( ρ v )g = 0 ∂x ∂y The momentum equation in the gas region is considered only in the axial direction, and constant density over the cross-section is assumed ∂u g ∂u g dpdyn ∂ ª ∂u g º ′ (ρu) g + ( ρ v) g =− + «( μ g + μ g ) (9.43) » + ( ρ g − ρ0 ) g ∂x ∂y dx ∂y ¬ ∂y ¼
where pdyn is dynamic pressure defined as
12 pdyn = p − ρ u0 (9.44) 2 The energy equation includes energy storage by the gas, conduction through the gas, and convection of the vapor component through the gas.
( ρ c pu) g
∂Tg ∂x
+ ( ρ c p v) g
∂Tg ∂y
(9.45)
∂T º ∂T ∂ω ∂ª ′ ′ = «( k g + k g ) g » + ρ g (α g + α g )( c p ,v − c p , g ) g ∂y ¬ ∂y ¼ ∂y ∂y
where the diffusivity is corrected to reflect the effect of surface waves. The last term on the right-hand side of eq. (9.45) represent the contribution of mass concentration gradient on the energy balance. A concentration balance is also required: ∂ω ∂ω ∂ ª ∂ω º (ρu) g (9.46) + ( ρ v) g = « ρ g ( D + D′ ) ∂x ∂y ∂y ¬ ∂y » ¼ The boundary conditions at the inlet are now specified with uniform temperature, velocity, and concentration: 2 x = 0: u = u0 T = T0 ω = ω0 pdyn = − ρ u0 / 2 (9.47) At the center of the channel, symmetric conditions must be satisfied, i.e., ∂u ∂T ∂ω = 0, = 0, =0 (9.48) y = 0: ∂y ∂y ∂y At the gas-liquid interface, the velocity, temperature, and concentration are at their interface values: y = b − δ : u = uδ , T = Tδ , ω = ωδ (9.49) With the nonslip condition, the shear stresses act equally and in opposite directions: ª ∂u º ª ∂u º τ δ = «( μ + μ ′ ) » = «( μ + μ ′ ) » (9.50) ∂y ¼ A ¬ ∂y ¼ g ¬
Chapter 9 Evaporation
693
The mass balance that accounts for the mass flux convected from the interface is ρ ( D + D′ ) § ∂ω · m′′ = (9.51) ¨ ¸ 1 − ωδ © ∂y ¹δ An energy balance that accounts for the heat required to release the mass from the liquid film is ª ∂T º ª ∂T º (9.52) «( k + k ′ ) ∂y » = «( k + k ′ ) ∂y » + m′′hAv ¬ ¼A ¬ ¼g The mass fraction at the interface is assumed to be that for thermodynamic equilibrium M v pδ ωδ = (9.53) M g ( p − pδ ) + M v pδ where Mv and Mg are the molecular masses of vapor and gas, respectively. At the channel wall there is also a nonslip condition with the liquid, and the wall is adiabatic –
y = b : uA = 0,
∂TA =0 ∂y
(9.54)
At the channel exit the dynamic pressure falls to zero:
x = 0:
pdyn = 0
(9.55)
Corrections to the viscosity due to presumed surface waves on the film are (Yan and Lin, 1991): ′ μA = cA ρAδ ′2 fd A (9.56)
′ μ g = cg ρ g δ ′2 fd g
(9.57)
where the values of cA and cg are 0.1 and 10, respectively. The amplitude, δ ′(m) , and frequency, f(Hz), of the wave are estimated by the following empirical correlations: −4 0.5505 We ≤ 0.038357 ° 4.1236 × 10 We δ′= ® (9.58) −4 0.2472 We > 0.038357 °1.5340 × 10 We ¯ We ≤ 0.0264443 86.70We0.371327 f =® 0.319413 We > 0.0264443 ¯71.79We where We is Weber number, defined as (9.59)
(9.60) Ω and the Reynolds number, Reδ , and surface tension parameter, Ω , are 4Γ Reδ = (9.61)
( Reδ / 4 ) We =
5/3
μA
694
Transport Phenomena in Multiphase Systems
Ω = 0.6736219
σ ρAν A4 / 3
(9.62)
The damping functions d A and d g in eqs. (9.56) and (9.57) are obtained by
§ b− y· (9.63) d A = 1 − exp ¨ − δ¸ © ¹ ª 5(b − δ − y ) º (9.64) d g = exp « − » δ ¬ ¼ The correction on thermal diffusivity, α ′ , can be estimated by introducing the wave Prandtl number, Pr ′ = ν ′ / α ′ . Similarly, the correction on mass diffusivity is estimated by introducing the wave Schmidt number, Sc′ = D′ / α ′ . It is recommended that the values of wave Prandtl and Schmidt numbers be unity and that the problem be solved numerically. Figure 9.8 shows predicted axial velocity, temperature and mass fraction for evaporation of falling ethanol liquid film. It can be seen from Fig. 9.8(a) that the inlet velocity of the gas stream is nearly uniform at the inlet and that the velocity profile becomes distorted with the maximum velocity shifted away from the central line. As can be seen from Fig. 9.8(b), the temperatures in the liquid and gas near the entrance (x/L=0.039 and 0.187) monotonically increase with y, which means that the direction of heat transfer is from channel wall to the liquidvapor interface, and then to the gas phase. At downstream (x/L = 0.468 and 0.93), the liquid and gas temperatures near the liquid-vapor interface is higher than the
(a)
(b)
(c)
Figure 9.8 Predicted axial (a) velocity, (b) temperature, and (c) mass fraction profile for evaporation of ethanol film ( TA ,i = 40 D C, Γ = 0.02 kg/m-s ; Yan and Lin, 1991; Reprinted with
permission from Elsevier).
Chapter 9 Evaporation
695
interfacial temperature, i.e., the latent heat of vaporization is supplied from both liquid and gas phases. The mass fraction of the ethanol in the central region of the channel increases with increasing x due to evaporation [see Fig. 9.8 (c)]. Experimental results using ethanol revealed that a larger inlet film temperature increases the liquid film cooling; this was to be expected because the temperature differential between the film and gas was increased. Lowering the liquid mass flow rate also increases the rates liquid film cooling; this is most likely due to longer residence times in the channel. It was also found that the exit temperature of the film could drop below the ambient temperature, which can be expected due to the volatility of ethanol. Yan et al. (1991) performed an experimental study of evaporative cooling of falling liquid film in an adiabatic vertical channel. Their experimental data yielded an empirical correlation for the mean heat flux, q′′, of water film evaporation: q′′ = 0.0113TA3.15b1 3 (9.65) ,0 where TA ,0 is the liquid film inlet temperature. It can be seen from eq. (9.65) that the heat flux is not a function of the mass flow rate of the liquid film.
Example 9.2 Falling water film evaporation occurs on the inner sides of two parallel plates separated by a distance of 2b = 3 cm. The length of the parallel plates is 1 m and the inlet temperature of the liquid film is 40 ÛC. The mass flow rate per unit width is Γ = 0.04 kg/s-m . Find the liquid film temperature at the outlet of the vertical channel. Solution: The specific heat of water should be evaluated at its mean temperature. Since the exit temperature is unknown, the inlet temperature is approximated as the mean temperature (which will be verified later). The specific heat at inlet temperature (40 ÛC) is c pA = 4.179 kJ/kg- o C .
The mean heat flux on the liquid film can be obtained from eq. (9.65), i.e., q′′ = 0.0113TA3.15b1 3 = 0.0113 × (40)3.15 (0.015)1 3 = 310.17 W/m 2 ,0 The mean heat flux is related to the liquid film temperature by q′′ = Γc pA (TA ,0 − TA , L ) / L (9.66) where the length of the parallel plates L = 1 m. The liquid film temperature at outlet is then q′′L 310.2 × 1 TA , L = TA ,0 − = 40 − = 38.1 DC Γc pA 0.04 × 4.179 × 103 Since the temperature change in the liquid film is very insignificant, the assumption that the inlet temperature was equivalent to the mean temperature was reasonable.
696
Transport Phenomena in Multiphase Systems
Example 9.3 Humid air with a temperature of T0 and mass fraction of water of ω0 enters at the top of a vertical channel composed of two parallel plates separated by a distance of b (Debbissi, et al., 2003; see Fig. 9.9). The inlet velocity and pressure of the humid air are respectively u0 and p0. There is an extremely thin liquid film on the inner side of the left wall and the outer side of the wall is insulated. The right wall is dry and maintained at a constant temperature, Tw. Since the humid air is cooled and its humidity increases as it flows through the vertical channel, the downward flow is accelerated by the chimney effect. Specify the governing equations and the corresponding boundary conditions of the problem. Solution: Since the liquid film is extremely thin, its thickness can be neglected and one can assume that the gas-vapor mixture flows in the entire channel (0<y<b). Since the flow is mainly in the x-direction, and the velocity component in the y-direction is several orders of magnitude smaller than that in the x-direction, one can assume that the flow in the channel is boundary layer type, i.e., diffusion in the x-direction can be neglected. Thus, the continuity equation is ∂ρ u ∂ρ v (9.67) + =0 ∂x ∂y
Figure 9.9 Evaporation into humid air between parallel plates.
Chapter 9 Evaporation
697
The momentum equation in the x-direction is ∂u ∂u 1 dpdyn ∂ 2u u +v =− − g β (T − T0 ) + β m g (ω − ω0 ) + ν 2 (9.68) ρ dx ∂x ∂y ∂y where and m are the coefficients of thermal expansion and concentration expansion, respectively. Neglecting viscous dissipation and assuming constant thermal diffusivity and specific heats, the energy balance is ∂T ∂T ∂ 2T c pv − c pg ∂T ∂ω u +v =α 2 + Dvg (9.69) cp ∂x ∂y ∂y ∂y ∂y where the second term on the right-hand side represents energy transport through interspecies diffusion, and Dvg is binary mass diffusivity from vapor to gas. The mass fraction of vapor in the mixture, , satisfies ∂ω ∂ω 1 ∂ § ∂ω · u +v = (9.70) ¨ ρ Dvg ¸ ∂x ∂y ρ ∂y © ∂y ¹ which is applicable in cases where the mass diffusivity is not a constant. The mass flow rate at any location x equals the mass flow rate at the inlet (x=0) plus accumulated evaporation from the left wall (0, x), i.e.,
³
b
0
ρ udy = ρ0u0b + ³ ρ v
0
x
y =0
dx
(9.71)
Since this is a boundary layer problem, we will need to specify boundary conditions at the inlet, and the left and right walls. The boundary conditions at the inlet are similar to those of eq. (9.47), i.e., 2 (9.72) x = 0: u = u0 T = T0 ω = ω0 pdyn = − ρ u0 / 2 The velocity components in the x-direction at both walls should be zero, as required by the non-slip condition: u = 0, y = 0, b (9.73) The velocity component in the y-direction at the left wall is determined by the evaporation rate, i.e., m′′ q′′ k ∂T =− = v= , y=0 (9.74) ρ ρ hAv ρ hAv ∂y y =0 while the velocity in the y-direction at the right wall is zero, i.e., v = 0, y = b (9.75) The latent heat of evaporation is supplied by the humid air flow through the vertical channel; therefore, the boundary condition of the energy equation at the left wall is ρ hAv Dvg ∂ω ∂T −k = , y=0 (9.76) ∂y y = 0 1 − ω y = 0 ∂y y = 0
698
Transport Phenomena in Multiphase Systems
The temperature of the right wall is constant, i.e., T = Tw , y = b The mass fraction of the vapor at the left wall is pv M v ω y =0 = pv M v + ( p − pv ) M g
(9.77) (9.78)
where pv is the partial pressure of the vapor at the left wall, which equals the saturation pressure corresponding to the temperature of the left wall. Since the right wall is dry, there is no mass transfer at the right wall, i.e., ∂C = 0, y = b (9.79) ∂y
Example 9.4 Determine the rate at which volatile liquid A evaporates into a stagnant gas B in a tube configuration as shown in Fig. 9.10, where the total pressure and temperature remain constant. Solution: This is the classic Stefan’s diffusion problem. The following assumptions are made to solve this problem:
1. 2. 3. 4. 5. 6.
Steady state 1D problem Constant total pressure and temperature Liquid A is permeable to gas B Ideal gas for gas mixture Flat equilibrium condition at the liquid/vapor interface Gas stream B
Gases A&B H y
Liquid A
Figure 9.10 Evaporation of liquid A into a stagnant gas B.
Chapter 9 Evaporation
699
Based on the above six assumptions, the continuity and species equations are reduced to the following form d ρv =0 (9.80) dy
§ dω A · d¨ρ ¸ dy ¹ dω ρ v A = DAB © (9.81) dy dy where A is the mass fraction of species A in the gas mixture and v is the velocity in the y-direction. It should be noted that energy and momentum equations are unnecessary since temperature and pressure are constant. However, mixture density is not constant. It is customary to solve this problem in terms of mole fraction, for simplicity, rather than mass fraction. The appropriate boundary conditions are At vapor/liquid interface y = 0 , ω A = ω A,δ v B = 0 (9.82) At the top of the tube y = H , ω A = ω A, H (9.83)
The velocity of the gas mixture, v, is related to the velocity of species, i.e., v = ω A v A + ω B vB (9.84) At the interface, vB = 0 , since liquid is impermeable to gas B. The mixture velocity at the interface is v = ω Av A (9.85) The rate of evaporation of species A at the interface is dω A m′′ = ρ A v A = ρ A v − ρ DAB (9.86) A dy Substituting eq. (9.85) into eq. (9.86) yields the fully related equation for evaporation at the interface dω A − ρ DAB dy m′′ = ρ Av A = (9.87) A 1 − ωA For a one-dimensional, steady-state problem dm′′ A (9.88) =0 dy Substituting eq. (9.87) into eq. (9.88) yields d § − ρ DAB · d ω A (9.89) =0 ¨ ¸ dx © 1 − ω A ¹ dy
− ρ DAB d ω A = m′′ = constant A 1 − ω A dy (9.90)
700
Transport Phenomena in Multiphase Systems
Since the gas mixture is assumed to be an ideal gas and isothermal, the total pressure is the sum of the partial pressures of components A and B p = p A + pB (9.91) where p = ρ RuT / M , p A = ρ A RuT / M A , and pB = ρ B RuT / M B Further simplification reduces eq. (9.91) to ρ1 1 ρA 1 = +B (9.92) M ρ MA ρ MB Using eq. (9.92) and definitions, with some simplification the mass fraction is pA ( M A / M B ) ωA = (9.93) p − p A (1 − M A / M B ) Similarly, the mixture density is given by § M A ·º pM M B ª = ρ= « p − p A ¨1 − ¸» RuT RuT « © M B ¹» ¬ ¼ Substituting eqs. (9.93) and (9.94) into eq. (9.90) yields § pM A ·§ DA→ B · dp A −¨ = m′′ = constant ¸¨ ¸ A © RuT ¹ © p − p A ¹ dy
(9.94)
(9.95)
Integrating eq. (9.95) and applying the boundary conditions yields § p − p A · y § p − p A, H · (9.96) = ln ln ¨ ¨ p− p ¸ H ¨ p− p ¸ ¸ ¨ ¸ A,δ ¹ A,δ ¹ © © where p A,δ and p A, H are the partial pressures of species A at the
liquid/vapor interface and at the top of the tube, respectively. Substituting eq. (9.96) into eq. (9.95) yields the rate of evaporation of species A pM A DAB § p − p A, H · (9.97) m′′ = ln ¨ A ¨ p− p ¸ ¸ RuTH A,δ ¹ © where temperature T and pressure p are assumed to be uniform within the gas mixture.
9.4 Falling Film Evaporation on a Heated Wall
As mentioned previously, falling film evaporators are used in both mass transfer and heat transfer applications. In general, the film’s mass flow rate decreases as evaporation occurs. If the film is being concentrated (as in fruit juice concentration), evaporation is also accompanied by an increase in viscosity. With the decrease in mass flow rate and the corresponding increase in viscosity, the Reynolds number decreases and the flow tends to progress from turbulent (if that ever existed) to wavy laminar, and then to laminar flow.
Chapter 9 Evaporation
701
9.4.1 Classical Nusselt Evaporation
The simplest evaporation heat transfer coefficient calculation uses the Nusselt treatment of falling film evaporation, as shown in Fig. 9.11. A flat plate with constant temperature Tw is placed vertically and a liquid film flows downward due to gravity. Evaporation takes place on the surface of the liquid film because the temperature of the vertical plate is above the saturation temperature of the vapor corresponding to the partial pressure of the vapor in the mixture. Since the heat must be transferred from the heat source (vertical wall) to the liquid-vapor interface, where evaporation takes place, this is an example of heterogeneous evaporation. The following assumptions are made: 1. The vapor is quiescent. 2. Inertia in the laminar liquid film is negligible. 3. Heat transfer through the liquid film is by conduction only, i.e., convective terms in the energy equation are negligible. 4. Constant wall temperature. 5. The process occurs at steady state. 6. A nonslip condition exists at the wall.
0,
Re0 y, v Vapor (v)
Liquid (A ) Tw Tv
g
L,
ReL
x, u
Figure 9.11 Schematic of Nusselt evaporation.
702
Transport Phenomena in Multiphase Systems
Proceeding from the assumptions of negligible inertia and constant fluid properties, the momentum balance is d 2u μA 2A = − ( ρA − ρv ) g (9.98) dy Integrating eq. (9.98) twice, with respect to y and using the boundary conditions of non-slip at the wall (u=0 at y=0) and zero shear at the liquid-vapor interface ( ∂u / ∂y = 0 at y = δ), the following velocity profile is obtained:
uA =
(ρ A − ρ v )g § ¨
μA
©
δ
¨ δy −
y2 · ¸ 2¸ ¹
(9.99)
The mass flow rate per unit width is defined by integrating the velocity profile from the wall to the interface.
Γ = ³ ρ A udy
0
(9.100)
The liquid velocity profile, eq. (9.99), is substituted into eq. (9.100) and the result is integrated from the wall to the interface. The result is then rearranged to obtain liquid film thickness:
ª º 3μ A Γ δ =« » ¬ ρ A ( ρ A − ρ v )g ¼
13
(9.101)
Since heat transfer through the film is by conduction only, the heat transfer coefficient is hx = kA / δ . Considering eq. (9.101), the local heat transfer coefficient becomes
ª ρ ( ρ − ρv ) g º (9.102) hx = kA « A A » 3μ A Γ ¬ ¼ Note that the film thickness δ drops out of the equation, and the heat transfer coefficient h is now a function of the mass flow rate per unit width, Γ . The local energy balance across the liquid film is dΓ (9.103) hAv = − hx (Tw − Tv ) dx Combining eqs. (9.102) and (9.103), one obtains,
13
k ª ρ ( ρ − ρv ) g º Γ dΓ = − A « A A » (Tw − Tv ) dx 3μ A hAv ¬ ¼ Integrating from the inlet to some arbitrary length L,
13 13 ΓL 13 13 L
(9.104)
k ª ρA ( ρA − ρv ) g º (9.105) » (Tw − Tv ) dx ³Γo Γ d Γ = − ³0 hAAv « 3μA ¬ ¼ A relation for at outlet is obtained in terms of the inlet mass flow rate, Γ0 , for the case of constant wall and vapor temperatures:
Chapter 9 Evaporation
703
4 kA ª ρ A ( ρ A − ρ v ) g º Γ =Γ − (9.106) « » (Tw − Tv ) L 3 hAv ¬ 3μA ¼ Now the mass flow rate per unit width is related to the Reynolds number by
13 43 L 43 o
(9.107) μA The local heat transfer coefficient from eq. (9.102) can be expressed in terms of Reynolds number as follows:
13 ª º μA2 1 §4· =¨ ¸ = 1.10 Re−1/ 3 (9.108) « » 13 © 3 ¹ Re « ρA ( ρA − ρ v ) g » ¬ ¼ Equation (9.108) for local heat transfer coefficient is valid for both laminar condensation and evaporation at both constant wall temperature and constant heat flux at wall. To calculate the average heat transfer coefficient, h , the energy flux over the entire length 0 < x < L is balanced with the mass flux from the interface: _ h ( Γ − Γ L ) μA hAv ( Re0 − Re L ) h = Av 0 = (9.109) L (Tw − Tv ) 4 L (Tw − Tv )
Re =
4Γ
hx kA
13
Equation (9.109) is valid for laminar, wavy laminar, and turbulent with constant wall temperature. Equation (9.103), which reflects the local energy balance across the liquid film, can be written in terms of Reynolds number as 4 (Tw − Tv ) d Re =− dx (9.110) hx μA hAv Integrating eq. (9.110) in the interval of (0, L), one obtains Re L d Re 4 (Tw − Tv ) L =− (9.111) Re0 μA hAv hx Combining eqs. (9.109) and (9.111), the average heat transfer coefficient is expressed as _ ( Re − Re L ) h = − Re 0 (9.112) L1 d Re Re0 h x Equation (9.112) for h is valid for laminar, wavy laminar, and turbulent flow. Substituting eq. (9.108) into eq. (9.112), an expression for the average heat transfer coefficient for laminar flow in terms of Reynolds number is obtained:
³
³
_
h kA
ª º μ A2 « » ¬ ρ A ( ρ A − ρ v )g ¼
13
§ 4· =¨ ¸ © 3¹
43
(Re
(Re o − Re L )
43 o
− Re 4 3 L
)
(9.113)
In order to use eq. (9.113) to obtain the average heat transfer coefficient, it is necessary to know Reynolds number at x=L, Re L . It can be obtained by
704
Transport Phenomena in Multiphase Systems
rewriting the relation for follows:
43 43 L 43 o
, eq. (9.106), in terms of Reynolds number, as
13
§ 4 · kA (Tw − Tv ) L ª ρ A ( ρ A − ρ v ) g º (9.114) Re = Re − 4 ¨ ¸ « » μA4 3 hAv 3μA ©3¹ ¬ ¼ The above analysis is valid for falling film evaporation on a vertical wall maintained at a constant temperature. Falling film evaporation from a vertical wall with constant heat flux is discussed in Example 9.5.
Example 9.5 A liquid film flows downward and evaporates on a vertical plate subjected to a constant heat flux. Determine the average heat transfer coefficient. Solution: The physical model of the problem under consideration is similar to that shown in Fig. 9.11 except that the thermal boundary condition is changed. The local heat transfer coefficient for falling film evaporation can still be calculated by eq. (9.102) or (9.108), because no assumption about boundary conditions were used to obtain them. For the case of constant heat flux, the energy balance accounts for the heat addition required for evaporation, i.e., dΓ ′′ = − qw (9.115) hAv dx ′′ Since qw is not a function of x, eq. (9.115) can be integrated in the interval of (0, L) to yield q′′ Γ = Γ0 − w L (9.116) hAv Equation (9.116) can be rewritten in terms of Reynolds number as ′′ 4 qw Re L = Re0 − L (9.117) μA hAv which can be used to obtain the Reynolds number at x=L. The average heat transfer coefficient can be obtained by _ h ( Γ − Γ L ) μA hAv ( Re0 − Re L ) h = Av 0 = (9.118) L (Tw − Tv ) 4 L (Tw − Tv )
which is similar to eq. (9.109) except that the wall temperature is the average value along the wall (0<x<L), i.e.,
Tw = ³ Tw dx
0
L
(9.119)
For the case with constant heat flux, the energy balance across the liquid film in terms of Reynolds number can still be represented by eq. (9.110). Integrating eq. (9.110) in the interval of (0, L) and considering the fact that Tw is not a constant, one obtains
Chapter 9 Evaporation
705
³
Re L
Re0
4 (Tw − Tv ) L d Re =− μA hAv hx
(9.120)
Combining eqs. (9.118) and (9.120) to eliminate (Tw − Tv ) yields the same equation as eq. (9.112), which eventually leads to eq. (9.113). Therefore, eq. (9.113) is also valid for the case with constant heat flux heating, while the Reynolds number at x = L must be calculated by eq. (9.117).
9.4.2 Laminar Falling Film with Surface Waves
Wavy flows of thin liquid films have higher heat transfer coefficients than smooth thin films. This effect is due to the former’s greater interfacial surface area and mixing action. Faghri and Seban (1985) analyzed a system in which a liquid film at an initial nondimensional temperature of zero flowed down a vertical wall and evaporated at the free liquid-vapor interface. The laminar liquid film thickness varied sinusoidally and the Reynolds number ranged between 35 and 472. The energy equation for a two-dimensional situation using the same vertical wall shown in Fig. 9.11 is as follows:
§ ∂ 2T ∂ 2T · ∂T ∂T ∂T +u +v = α¨ 2 + 2 ¸ ¨ ∂x ∂t ∂x ∂y ∂y ¸ ¹ ©
The continuity equation for this system is
(9.121)
∂u ∂v =0 + ∂x ∂y
§ η2 · u = 3u ( x, t )¨η − ¸ ¨ 2¸ © ¹ where the nondimensional coordinate η is defined as y η=
(9.122)
Faghri and Seban (1985) assumed a quasi-parabolic velocity profile in the xdirection as follows: (9.123)
and u is the local mean velocity. Equation (9.122) can be rearranged and integrated with respect to y to obtain the velocity component v: y ª ∂u § η 2 η 3 · ∂u ∂δ § η 2 η 3 · º (9.125) v = − ³ dy = −3 « δ ¨ − ¸ − u ¨ − ¸» 6¹ 3 ¹¼ ∂x © 2 0 ∂x ¬ ∂x © 2 It can be seen that eq. (9.123) has been substituted in order to find an expression for the v component of the velocity where the local mean velocity and film thickness δ are still unknown as functions of the wavy surface. To determine these two variables, the following mass balance can be written as a function of time and space:
δ
(9.124)
706
Transport Phenomena in Multiphase Systems
∂δ ∂δ ∂ = − ³ udy = − (u δ ) 0 ∂t ∂x ∂x
(9.126)
The film thickness can be written in terms of its average value, δ , and a local amplitude φ. δ = δ (1 + φ ) (9.127) Assuming that the wave formation has a periodic characteristic with a velocity c, the following relations can be written for the fluctuation of local mean velocity and film thickness:
∂u ∂u = −c ∂t ∂x ∂δ ∂δ = −c ∂t ∂x
(9.128) (9.129)
Combining eqs. (9.126) – (9.129) and integrating, the following expression is obtained: (c − u )(1 + φ ) = [c − u0 ] (9.130) where u0 is the average velocity for the average film thickness δ . For very small local amplitudes, φ < 1, u may be approximated by expanding eq. (9.130) and neglecting the third-order terms. u = u0 + ( c − u0 ) φ − ( c − u0 )φ 2 (9.131) Substituting eqs. (9.127) and (9.131) into eq. (9.125), one obtains the velocity profile in terms of the wavy surface parameters: 2 3 ∂φ °§ c − 1· (1 − 2φ ) δ § η − η · v = −3u0δ ¸ ¨ ¸ ®¨ δ©2 6¹ ∂x ¯© u0 ¹ ° (9.132) ½ ª §c · §c · 2 º § η 2 η 3 ·° − «1 + ¨ − 1¸ φ − ¨ − 1¸ φ » ¨ − ¸ ¾ 3 ¹° « © u0 ¹ © u0 ¹ »© 2 ¬ ¼ ¿
Assuming the wave is sinusoidal,
ª§ 2π · º (9.133) ¸( x − ct )» ¬© λ ¹ ¼ and introducing the nondimensional variable ξ defined as § 2π · ξ = ¨ ¸( x − ct ) (9.134) ©λ¹ where A is magnitude of the wave, and λ is the wavelength, eq. (9.121) is transformed from the (t, x, y) to the (ξ, η) coordinate system,
φ = A sin «¨
C1
∂T ∂T ∂ 2T ∂ 2T ∂ 2T + C2 = C3 2 + C4 + C5 2 ∂ξ ∂η ∂η ∂η∂ξ ∂ξ
(9.135)
Chapter 9 Evaporation
707
where the constants are as follows:
C1 =
C2 = 2πη δ
2π
λ
(u − c )
2
(9.136)
§ 2π · δ (c − u ) A cos ξ − ¨ ¸ α η A sin ξ λδ ©λ¹δ
2
v § 2π · § δ · 2 2 − 2¨ ¸ α ¨ ¸η A cos ξ + δ © λ ¹ ©δ ¹
(9.137)
α § 2π · C3 = 2 + α ¨ ¸ δ ©λ¹
2
§δ ¨ ¨δ ©
2
· 22 ¸ η A cos 2 ξ ¸ ¹
2
(9.138)
§ 2π · δ C 4 = −2α ¨ ¸ η A cos ξ ©λ¹ δ § 2π · C5 = α ¨ ¸ ©λ¹
2
(9.139) (9.140)
These coefficients are evaluated with u from eqs. (9.123) and (9.131), and v from eq. (9.132). In the evaporating film with a fluid with a high latent heat of vaporization, the average film thickness over a wave length, δ , will not vary substantially with distance, x. Furthermore, the temporal average temperature profile will not vary with x far from where heating begins. Therefore, the boundary condition for eq. (9.135) that describes evaporation on a wavy film is approximated as below. The equality of temperatures at corresponding points in the period gives T (0,η ) = T ( 2π ,η ) (9.141) Boundary conditions for dimensionless temperature at the heated vertical wall and wavy surface are, respectively, T (ξ ,0) = 1 (9.142) T (ξ ,1) = 0 (9.143) Therefore, the energy equation with the above boundary conditions can be solved to obtain the wavy condensate film’s temperature profile. With this profile, the heat flux at the wall can be found from
′′ qw = −
k ∂T δ ∂η η =0
(9.144)
and the heat flux normal to the wavy surface (η = 1) is as follows:
· q′′ ª 2π § δ 1 º ∂T δ = « ¨ A cos ξ ¸ i − j» (9.145) k ¬ λ ©δ ¹ δ ¼ ∂η η =1 where i and j are unit vectors in the x- and y-directions, respectively. The average heat flow over a period 0 < ξ < 2π is found by integrating eqs. (9.144)
and (9.145). For a sinusoidal wavy layer, such as the one considered here, this
708
Transport Phenomena in Multiphase Systems
average is identical to the average obtained for the actual surface length. The average heat flux at the wall (η = 0) is approximated by
′′ qw 1 = k 2π
³
2π
0
1 ∂T dξ δ ∂η 0
ª § ∂δ · 2 º «1 + ¨ ¸ » dξ « © ∂x ¹ » ¬ ¼
1/ 2
(9.146)
and the average normal heat flux at the wavy surface is as follows: ′′ qδ 1 = k 2π
³
2π
0
1 ∂T δ ∂η η =1
(9.147)
Table 9.1 Comparison of wavy film analysis with experimental data (Faghri and Seban, 1985).
References 1 2 3 T(°C) 4Γ / μ
Rogovan et al. (1969) 20 35 4.132
δ (10-4ft)
w0 (ft/s)
2π / λ (ft )
-1
Kapitza and Kapitza (1975) 25 92 5.96
Rogovan et al. (1969) 20 268 9.41
Rogovan et al. (1969) 20 472 11.5
4 5 6 7 8 9 10 11 12 13
*
0.279 191 2.12 0.30 7.2 1.7 1.057 1.056
0.383 241 1.93 0.52 6.2 1.309 1.314 1.17 1.171 1.171 1.7 1.296 1.300 1.169 1.169
0.836 159 1.65 0.55 7.2 1.567 1.611 1.19 1.199 1.197 1.7 1.543 1.580 1.199 1.197 0.40 7.2 1.393 1.433 1.09 1.092 1.091 1.7
1.164 136 1.55 (0.55) 7.0 1.724 1.816 1.19 1.199 1.197 1.7 1.694 1.772 1.199 1.197
c / u0
A
ν /α
hδ / k hδ / k
w
1.059 1.058 1.05 1.049 1.049
1.372 1.405 1.092 1.091
δ
1/ 1 − A2
hδ / k hδ / k
*
w
1.048 1.048
*
δ
Conduction. Reprinted with permission from Elsevier.
The calculations were made using 40 increments in both ξ and η for the value of 2π / λ , A, and c / u0 . Table 9.1 lists these experimental determinations made by the cited authors. Rows 1 through 7 give the experimental conditions and measurements. Row 8 gives the Prandtl numbers for which the calculations were made. One, on the order of 7, corresponds to the experimental conditions, and the other, 1.7, was selected as a comparison to show the effect of Prandtl number. Row 9 gives the calculated average Nusselt number at the wall, and Row 10 gives the average Nusselt number at the outer edge of the layer. These values should be the same, and the difference indicates the failure of the calculation to satisfy the energy balance. This difference is small; it increases with the Reynolds number, reflecting some increase in truncation error. Row 11 gives the
Chapter 9 Evaporation
709
average Nusselt number as evaluated for unidimensional conduction. On this basis, the local Nusselt number is hδ / k = 1 , and then
h 1 = k 2π
For the sinusoidal wave, h 1 2π dξ = ³0 1 + A sin ξ k 2πδ
³
2π
dξ
0
δ
(9.148)
(9.149) π º 1 ªπ dξ dξ 1 = + = 2πδ « ³0 1 + A sin ξ ³0 1 − A sin ξ » δ 1 − A2 ¬ ¼ The difference between the calculated average Nusselt numbers is due to the contribution of convection and to the two-dimensional nature of the conduction that exists because of the variation in layer thickness. The result for w = v = c = 0 corresponds to the two-dimensional conduction solution for the wave, and the Nusselt numbers for such a calculation are shown in Rows 12 and 13. These Rows should be identical, and the truncation error in the calculation creates a discrepancy between them, as it does in the cases in which there is fluid motion. These values are essentially the same as Row 11, for unidimensional conduction, and this correspondence shows that two-dimensional conduction effects are negligible. Therefore, the difference between the Nusselt numbers in Rows 9 and 11 reflect only the effect of convection. Figure 9.12 shows the values of the Nusselt number (h δ / k ) as given by Row 9 of Table 9.1 for a Prandtl number of about 7. The value of 472 for (4Γ / μ ) is also shown for an arbitrarily-higher value of A = 0.55, because the
Figure 9.12 The average Nusselt number for evaporation of a wavy falling film as a function of the Reynolds number: Curve A, Kutateladze (1982); Curve B, Kutateladze (1963); Curve C, Hirschburg and Florschuetz (1982); and for f+ = 0.65 (Faghri and Seban, 1985; Reprinted with permission from Elsevier).
710
Transport Phenomena in Multiphase Systems
hw / h hδ / h
hw / h
hδ / h
Figure 9.13 The local Nusselt number for evaporation of wavy falling film (Faghri and Seban, 1985; Reprinted with permission from Elsevier).
value of A probably should not decrease as the Reynolds number increases. Figure 9.11 also contains lines A – to show the specification of Zazuli given by Kutateladze (1963) – and B – to show that of Kutateladze (1982). Curve C is that of Hirschburg and Florschuetz (1982) for the intermediate wave solution designated by them as f+ = 0.65, which best fits the average of the data. Three data points from Chun and Seban (1971) for evaporation of water with a Prandtl number of about 5.6 are shown by plus (+) symbols, and four numerical results of Faghri and Seban (1985) are shown by circles ( ) in the same figure. Figure 9.13 shows the variation of the local heat transfer coefficient as a function of ξ , normalized with respect to the average value over a period. The results are for the Prandtl number of 8.2 and the Reynolds numbers of 35 and 472 (the same portrayal for the Prandtl number of 1.7 is not much different.) The ratio hw / h for the wall varies considerably for the low Reynolds number, but not very much for the high Reynolds number. For the surface, the ratio hδ / h varies substantially for both cases. The numerical results for the Nusselt numbers of Faghri and Seban (1985) are on the order of, but tend to be higher than, the correlation equations of Kutateladze (1963, 1982). Comparison of Rows 9 and 11 in Table 9.1 reveals the considerable effects of convection and two-dimensional conduction. For practical purposes, Chun and Seban (1971) suggested that the following correlation based on wavy laminar film condensation is valid for wavy laminar falling film evaporation:
§ Re · hx = 0.876 ¨ ¸ ©4¹
0.11
hNusselt
(9.150)
Chapter 9 Evaporation
711
An empirical correlation for the local heat transfer coefficient for laminarwavy flow is obtained by combining eq. (9.150) with eq. (9.108),
hx kA ª º μA2 « » « ρA ( ρA − ρv ) g » ¬ ¼
13
= 0.828Re −0.22 ,
Re > Re wavy
(9.151)
Chun and Seban (1971) suggested the following empirical correlation to predict the onset of wavy laminar flow: Re ≥ Re wavy = 2.43Ka −1/11 (9.152) where Ka is Kapitza number, a ratio of surface tension to viscous force: μA4 g (9.153) Ka = ( ρA − ρv ) σ 3 The average heat transfer coefficient is often relevant to the practical application and can be obtained by using eq. (9.112), which is also valid for laminar flow with waves. Substituting eq. (9.151) into eq. (9.112), one obtains an empirical correlation of the average heat transfer coefficient as follows:
_
h kA
ª º μA2 « » « ρA ( ρA − ρv ) g » ¬ ¼
13
=
( Re
( Reo − Re L )
1.22 o
− Re1.22 ) L
(9.154)
Substituting eq. (9.109) into eq. (9.154), an equation correlating Re L , the Reynolds number at x = L, is obtained:
Re
1.22 L
= Re
1.22 0
k L (Tw − Tv ) ª ρA ( ρA − ρv ) g º −4 A « » μA hAv μA2 ¬ ¼
13
(9.155)
Example 9.6 A uniformly-distributed liquid water film at 0.01 kg/s is introduced at the top of a tube with an inner diameter of 4 cm and a height of 2 m. The inner surface temperature of the tube is uniformly 105 °C and the saturation temperature of the steam is 100 °C. Find the evaporation rate and the average heat transfer coefficient. Solution: The properties of the liquid water and the steam can be approximated as those at the saturation temperature: kA = 0.68 W/m- o C,
ρA = 958.77 kg/m3 , μA = 2.79 × 10−4 N-s/m 2 , σ = 58.9 × 10−3 N/m, ρv =
0.596 kg/m3 , hAv = 2251.2 kJ/kg. It is assumed that the liquid film thickness is much smaller than the radius of the tube. Therefore, the falling film analysis for evaporation over a vertical flat plate can be used. The Kapitza number is obtained from eq. (9.153): μA4 g (2.79 × 10−4 ) 4 × 9.8 = = 3.032 × 10−13 Ka = −3 3 3 ( ρA − ρv )σ (958.77 − 0.596) × (58.9 × 10 )
712
Transport Phenomena in Multiphase Systems
Therefore, the Reynolds number for wavy film is Re wavy = 2.43Ka −1/11 = 2.43 × (3.03 × 10−13 ) −1/11 = 33.4 The inlet mass flow rate per unit perimeter is m 0.01 Γ0 = 0 = = 7.95 × 10−2 kg/s-m π D π × 0.04 The Reynolds number at the inlet is 4Γ 4 × 7.95 × 10−2 Re0 = 0 = = 1139.8 > Re wavy μA 2.79 × 10−4 which is also under the transition Reynolds number for turbulent, Returb = 5840 Pr −1.05 = 1862 [see eq. (9.156)]. This means that the liquid film is in the laminar wavy regime. The Reynolds number at outlet can be obtained from eq. (9.155), i.e., Re1.22 L kL = Re1.22 − 4 A 0
1.22
μA hfg
(Tw − Tv ) ª ρA ( ρA − ρv ) g º
« ¬
13
μA2
» ¼
13
= 1139.8
ª 958.77 × ( 958.77 − 0.596) × 9.8 º 0.68 × 2 × (105 −100) − 4× × » −4 3« 2.79 ×10 × 2251.2 ×10 ¬ (2.79 ×10−4 )2 ¼
= 3257.5 Thus
Re L = 757.6 which is still in the wavy regime so the entire surface is wavy. The mass flow rate of the liquid at the outlet is Re 757.6 mL = m0 L = 0.01 × = 6.64 × 10−3 Re0 1139.8 The evaporation rate is therefore mv = m0 − mL = 0.01 − 6.64 × 10−3 = 3.36 × 10−3 kg / m The heat transfer coefficient can be obtained from eq. (9.154), i.e.,
ªρ (ρ − ρ ) g º h = kA « A A 2 v » μA ¬ ¼
_ 13
( Re
( Reo − Re L )
1.22 o
− Re1.22 ) L
13
ª 958.77 × ( 958.77 − 0.596 ) × 9.8 º 1139.8 − 757.6 =0.68 × « »× −4 2 (2.79 × 10 ) 1139.81.22 − 757.61.22 ¬ ¼ =6017 W/m 2 - o C
Chapter 9 Evaporation
713
9.4.3 Turbulent Falling Film
Stephan (1992) presented an empirical relationship between Reynolds number and Prandtl number that could be used to determine whether a falling film is completely turbulent: Re ≥ Returb = 5840 Pr −1.05 (9.156) In order to model turbulent flow, the following development from Seban and Faghri (1976) proceeds with assumptions similar to those made in the classical Nusselt laminar analysis. The momentum equation is
d (vA + ε M ) du + g = 0 dy dy
where ε M is the turbulent eddy diffusivity. Defining the following dimensionless variables,
(9.157)
uτ = τ w ρ A , u + = u uτ , y + = yuτ vA
the momentum equation becomes
(9.158)
d dy +
· du + gv §εM ¨ + 1¸ + + 3A = 0 ¸ dy ¨v uτ ¹ ©A
(9.159)
which is subject to the following boundary conditions:
u+ = 0 ,
du = 0, dy +
+
y+ = 0
y+ = δ
+
(9.160) (9.161)
Integrating eq. (9.159) twice and using eqs. (9.160) and (9.161) to determine the integral constant, the dimensionless velocity becomes y+ 1 − y + δ + + (9.162) u= dy + 0 (1 + ε M vA )
³
(
)
The energy equation is
u
∂T ∂ ∂T = (ε H + α ) ∂x ∂y ∂y
(9.163)
where ε H is turbulent eddy thermal diffusivity. Equation (9.163) is subject to the following boundary conditions: T = Tsat x = 0 (9.164)
q′′ ∂T =− w, y=0 kA ∂y T = Tsat , y = δ
By defining dimensionless temperature as T ρ A cuτ T+ = ′′ qw
(9.165) (9.166) (9.167)
714
Transport Phenomena in Multiphase Systems
and using the dimensionless variables from eq. (9.158), the energy equation becomes
u+
∂T ∂ =+ + ∂x ∂y
§ α ε H · ∂T + ¨+ ¸+ ¨v ¸ © A v A ¹ ∂y
(9.168)
The solution requires choosing a turbulence model for the specification of ε M /ν and of ε H / ε M . The detailed turbulence model and solution procedure can be found in Seban and Faghri (1976). The predicted result agrees very well with the experimental result of Chun and Seban (1971) for constant heat flux at wall. To predict local heat transfer in turbulent falling film flow on a wall heated under constant heat flux, the following empirical correlation, recommended by Chun and Seban (1971), can be used: ª º μA2 0.4 0.65 −3 (9.169) « » = 3.8 × 10 Re Pr « ρA ( ρA − ρv ) g » ¬ ¼ The average heat transfer coefficient can be obtained by substituting eq. (9.169) into eq. (9.112), i.e., h kA
_
13
h kA
ª º μA2 « » ¬ ρA ( ρA − ρv ) g ¼
13
= 2.28 × 10−3
( Reo0.6 − Re0.6 ) L
( Reo − Re L )
Pr 0.65
(9.170)
As was demonstrated in Example 9.4, eq. (9.170) is valid for cases with either constant wall temperature or constant heat flux. However, the methods for determining the Reynolds number at x = L will vary for different thermal boundary conditions on the wall. For cases with constant wall temperature, the Reynolds number at x = L, Re L , is obtained by substituting eq. (9.109) into eq. (9.170), i.e.,
Re
0.6 L
k L (Tw − Tv ) ª ρA ( ρA − ρv ) g º = Re − 9.12 × 10 A « » μA hAv μA2 ¬ ¼
0.6 0 −3
13
(9.171)
For the case with constant heat flux, eq. (9.117) can be used to determine the Reynolds number at x = L, Re L . While gravitational force drives the falling film evaporation discussed here, liquid flow in a liquid film can also be driven by centrifugal force. Rahman and Faghri (1992) analyzed the processes of heating and evaporation in a thin liquid film adjacent to a horizontal disk rotating about a vertical axis at a constant angular velocity. The fluid emanates axisymmetrically from a source at the center of the disk and then is carried downstream by inertial and centrifugal forces. Closed-form analytical solutions were derived for fully developed flow and heat transfer. Simplified analyses were also presented for developing heat transfer in a fully-developed flow field. Rice et al. (2005) provided a detailed analysis for evaporation from a thin liquid film on a rotating disk including conjugate effects.
Chapter 9 Evaporation
715
9.4.4 Surface Spray Cooling
When liquid drops are sprayed onto a hot surface, heat conduction causes the liquid to superheat and evaporate into a gas (see Figure 9.2). The cooling capacity of the process depends on the time required for an average-sized drop to cool. If the shape of the droplet can be assumed to be hemispherical during the evaporation process, the energy balance at the interface can be written as π / 2 ∂T ( r ,θ ) dr AI ρA hAv 2π rI2 I = kA ³ 2π rI2 cos θ dθ − h 2π rI2 (T∞ − Tsat ) (9.172) 0 dt ∂r where rI is the radius of the droplet. Equation (9.172) can be simplified as π / 2 ∂T ( r ,θ ) dr AI (9.173) ρA hAv I = kA cos θ dθ − h (T∞ − Tsat ) 0 dt ∂r Analytical solution of eq. (9.173) is very difficult because the temperature distribution in the liquid droplet is two-dimensional in nature. A scale analysis similar to that in Lock (1994) is presented here to estimate the time required for droplet with an initial radius of Ri to evaporate completely. The scale of the radius of the droplet, rI, is Ri. Thus the scale of the temperature gradient at the interface is
³
∂T (rI , θ ) Tw − Tsat ∂r Ri
(9.174)
The scale analysis of eq. (9.173) yields R T −T ρA hAv i kA w sat + h(T∞ − Tsat ) tf Ri which can be rearranged as ρA hAv Ri2 hR (Tw − Tsat ) + i (T∞ − Tsat ) kA t f kA
(9.175)
(9.176)
Since the liquid droplet is very small, hRi / kA is very small, so the second term on the right-hand side of eq. (9.175) is much smaller than the first. Equation (9.175) can be simplified by neglecting the second term on the righthand side, i.e., R2 ρA hAv i (Tw − Tsat ) (9.177) kA t f The physical significance of eq. (9.177) is that the latent heat of evaporation is balanced primarily by conduction in the liquid droplet, while the effect of convection on the surface of the droplet is negligible. The time that it takes to completely evaporate the droplet with an initial radius of Ri is therefore estimated by ρA hAv Ri2 (9.178) tf kA (Tw − Tsat )
716
Transport Phenomena in Multiphase Systems
′′ For the cases where the substrate heat flux beneath the droplet, qw , is known, eq. (9.177) should be rewritten in terms of heat flux. The scale of the heat flux at the heating surface is k (T − T ) ′′ qw A w sat (9.179) Ri Combining eqs. (9.178) and (9.179) yields ρh R (9.180) t f A Av i ′′ qw which indicates that a smaller drop will provide a larger cooling effect. While eqs. (9.178) and (9.180) provide the order of magnitude of the time required to completely evaporate the droplet, quantitative estimation of the evaporation time is often also desirable. A simple approximate analysis will be done below by estimating the conduction in the liquid droplet by the following correlation: T −T (9.181) qd = kA A w sat δ where A and δ are the average cross-sectional area of heat conduction and the average path length of the conduction, respectively. For a hemispherical droplet, the contact area between the droplet and the heated wall is Aw = π rI2
and the interfacial area of the droplet is AI = 2π rI2 . Thus, we can take the average conduction area as 1 3 A = ( Aw + AI ) = π rI2 (9.182) 2 2 The average path length for conduction is V (2 / 3)π rI3 4 = rI (9.183) δ= = A (3/ 2)π rI2 9 Therefore, the conduction in the liquid droplet becomes 27 qd = π kA rI (Tw − Tsat ) (9.184) 8 Replacing the first term on the right-hand side of eq. (9.173) with eq. (9.184), and dropping the second term on the right-hand side of eq. (9.173) (see the above scale analysis), the energy balance for the hemispherical droplet becomes dr 27 ρA hAv 2π rI2 I = π kA rI (Tw − Tsat ) (9.185) 8 dt which can be rearranged as dr 27 kA (Tw − Tsat ) (9.186) rI I = dt 16 ρA hAv which is subject to the following initial conditions: rI = Ri , t = 0 (9.187)
Chapter 9 Evaporation
717
Integrating eq. (9.186) and considering its initial condition, eq. (9.187), one obtains the transient radius of the liquid droplet: 27 kA (Tw − Tsat )t (9.188) rI = Ri2 − ρ A hAv 8 The time required for the droplet to evaporate completely can be obtained by letting rI in eq. (9.188) equal zero, i.e., 8ρA hAv Ri2 (9.189) 27 kA (Tw − Tsat ) which agrees with the results obtained by scale analysis, eq. (9.178). The above analysis is otherwise simplistic because it assumes that the shape of the droplet is hemispherical. In fact, the shape of the droplet depends on the velocity with which it impacts the wall as well as the wettability of the liquid droplet on the wall. If the contact area between the droplet and the wall is larger (due to higher impacting velocity or good wettability), the time for conduction through the liquid is shorter and the life of the drop will be shorter as well. tf =
Example 9.7 A hemispherical liquid water droplet with an initial radius of 20 m is attached to a heated wall and exposed to pure water vapor at 1 atm. The temperature of the heated wall is 180 °C. Estimate the time required for the liquid droplet to evaporate completely. Solution: The saturation temperature corresponding to 1 atm is Tsat = 100 °C. The properties of the liquid water can be approximately evaluated at the saturation temperature: kA = 0.68 W/m- o C, ρA = 958.77 kg/m3 , hAv = 2251.2 kJ/kg. The time required for the liquid droplet to evaporate completely can be obtained from eq. (9.189), i.e., 8ρ A hAv Ri2 tf = 27kA (Tw − Tsat )
=
8 × 958.77 × 2251.2 × 103 × (20 × 10−6 ) 2 = 0.0047 s = 4.7 ms 27 × 0.68 × (180 − 100)
9.4.5 Evaporation from a Wedge or Cone Embedded in a Porous Medium
Heat and mass transport to a liquid film falling along a vertical or inclined surface in a porous medium occur in many types of chemical and energy production processes. Some examples are geothermal energy production, packedbed reactors and heat exchangers, oil recovery from petroleum reservoirs, heat transfer from buried vertical pipe lines, and gravity-assisted heat pipes. Flow and
718
Transport Phenomena in Multiphase Systems
y Free surface
2 x r g
Figure 9.14 Film evaporation on surface of wedge or cone embedded in a porous media.
transport of heat and mass in a thin liquid film on the surface of a wedge or cone embedded in a porous medium during evaporation, gas absorption at the free surface, and heating and solid dissolution at the wall have been studied by Rahman and Faghri (1993a,b) and will be presented here. The coordinate system used for the following analysis and numerical computation is shown in Fig. 9.14. The x-axis is directed along the length of the wedge or cone and the y-axis is perpendicular to the solid surface. The axis of symmetry of the wedge is oriented vertically along the direction of the gravitational force, while the solid wall is inclined by an angle . The continuity, momentum, and energy equations are ∂ (r n u ) ∂ (r n v) + =0 (9.190) ∂x ∂y ª ∂ § 1 ∂ (r nu ) · ∂ 2u º ν u g ( ρA − ρv ) ∂u ∂u 1 ∂p +v =− +ν « ¨ n + cos θ »− ¸+ ρA ∂x ρA ∂x ∂y ∂x ¹ ∂y 2 ¼ K ¬ ∂x © r (9.191) ª ∂ § 1 ∂ (r n v) · ∂ 2u º ν v g ( ρA − ρ v ) ∂v ∂v 1 ∂p u +v =− sin θ +ν « ¨ n »− + ¸+ ρA ∂y ρA ∂x ∂y ∂x ¹ ∂y 2 ¼ K ¬ ∂x © r u (9.192) ª ∂ § 1 ∂(r φ ) · ∂ φ º (9.193) « ¨n » ¸+ ∂x ¹ ∂y 2 ¼ ¬ ∂x © r Here, n = 0 or 1 denotes equations applicable to a wedge or a cone, respectively. The properties of the fluid and the medium are assumed to be u ∂φ ∂φ +v = Γ eff ∂x ∂y
n
2
Chapter 9 Evaporation
719
constant. The terms containing permeability quantify the resistance to the flow due to the presence of the porous matrix; as K → ∞ , the flow equations approach those for a liquid medium with no solid matrix. The variable φ represents temperature or concentration, and Γ eff is the corresponding diffusion coefficient for the different heat and mass transfer processes considered here. For heat transfer, Γ eff is the effective thermal diffusivity, which is defined as
Γ eff =
ε kA + (1 − ε )ks ερA c pA
(9.194)
where ε is the porosity of the medium. For mass transfer processes, mass diffusivity is ν Γ eff = A (9.195)
Sc
The porous matrix diffuses heat and therefore participates in the heat transfer process, but so far as mass diffusion is concerned, the fluid is the only active medium. An experimental or analytical correlation (if available) can be substituted for the linear averaging of thermal conductivity – as used in eq. (9.194) – without changing the analysis presented here. Due to the small rate of evaporation, the flow is not affected by the loss of fluid at the free surface. In that situation, the latent cooling due to evaporation has a negligible influence on heat transfer from the wall. Mass transfer during absorption or dissolution neither affects the flow or heat transfer, since the rate of mass diffusion is very small compared to the fluid flow rate, and the heat of absorption or dissolution is negligible compared to the sensible heat transfer. Therefore, the momentum and transport equations are decoupled, and the flow field can be solved independently without considering any associated transport processes. For a thin film flow, the thickness of the film is usually much smaller than the physical dimension of the object over which flow is taking place. Similarly, the velocity normal to the wall is much smaller than the velocity component parallel to the wall. A scaling argument can be used to show that convection in the flow direction is dominant over diffusion in that direction, and therefore one can neglect diffusion in the flow direction. Moreover, in a gravity-driven flow through a porous matrix, small effective Reynolds numbers cause the momentum boundary layer to develop quickly and encompass the entire thickness after a short distance. Therefore, inertia can be neglected in most regions. The convection of heat or mass can still be significant, particularly if the Prandtl or Schmidt number was large, and thus was retained here. Keeping only the most significant terms, the governing equations (9.190) – (9.193) were transformed to ∂ 2u μu (9.196) μ 2− + g ( ρA − ρv ) cos θ = 0 ∂y K
720
Transport Phenomena in Multiphase Systems
p = p0 + ( ρA − ρ v ) g (δ − y )sin θ
2
(9.197)
∂φ ∂φ ∂φ +v = Γ eff 2 (9.198) ∂x ∂y ∂y where p0 is the ambient pressure. Equation (9.196) is a second-order nonhomogeneous ordinary differential equation. Adding the homogeneous and particular solutions and applying boundary conditions at the wall (u = 0) and on the free surface (du/dy = 0), the following velocity distribution was obtained. ª cosh δ * (1 − y + ) º (9.199) u = uD «1 − » cosh δ * ¬ ¼ u where uD = K
μA
( ρA − ρv ) g cosθ
(9.200)
is the velocity obtained by direct application of Darcy’s law. The quantity δ * is the film height nondimensionalized with the permeability of the medium ( δ * = δ / K ). It can be noticed that u increases with y+ and becomes maximum on the free surface. The average film velocity was obtained by integrating the local velocity across the film thickness. § tan δ * · u = u D ¨1 − (9.201) ¸ δ* ¹ © It may be noted that eqs. (9.199) and (9.201) are applicable to both the wedge (n = 0) and the cone (n = 1), provided that for a cone, δ r (i.e., away from the vertex). The thickness of the film can be related to the velocity by satisfying the conservation of mass at any x-location. The volumetric flow rates for a wedge and a cone are respectively Q = uδ (9.202) Q = 2π xu δ sin θ (9.203) Substituting eq. (9.201) into eqs. (9.202) and (9.203), the following two equations for the film thickness are obtained. δ * − tan δ * = A for wedge (9.204) B δ * (δ * − tanh δ * ) 2 = + for cone (9.205) x where Q A= (9.206) uD K QΓ eff B= (9.207) 2 2π sin θ K 3 / 2u D
Chapter 9 Evaporation
721
(9.208) uδ 2 The solutions of eqs. (9.204) and (9.205) give the film thickness distribution for the wedge and cone, respectively. It may be noted that δ * remains constant for a wedge, whereas it decreases with x for the case of a cone because more area is available to the flow as the fluid moves downstream. The transport equation (9.198) can be nondimensionalized as ∂φ ∂ 2φ For wedge: u + + = +2 ∂x ∂y
º ∂φ y+ ∂ 2φ 1 · ∂φ u + ª 1 § = +2 For cone: u + ¨ 3 − * ¸ + + + « * − » δ ¹ ∂x x ¬ δ 1 − cosh( y +δ * ) + sinh( y +δ * ) ¼ ∂y + ∂y © (9.209) In addition to the simplified analytical solution for δ * , Rahman and Faghri (1993b) obtained a complete numerical solution of these equations using a twodimensional finite-difference method. To appropriately handle the variation of film height with flow, a curvilinear body-fitted coordinate system was used. The free surface formed one of the boundaries of the computation domain, and grid cells between the wall and free surface were generated by an algebraic interpolation between boundary values. The cell faces, in general, were nonorthogonal to each other. The equations were solved in the physical domain using the co-variant components of velocity and force vectors. The distribution of dimensionless film thickness along the wedge is shown in Fig. 9.15. Ethyl alcohol at 20 ºC was used for fluid properties of the film. The porous matrix was made of glass beads 0.98 mm in diameter, which provide a
x+ =
Γ eff x
Figure 9.15 Variation of film thickness along a wedge.
722
Transport Phenomena in Multiphase Systems
porosity of ε = 0.35 (for uniformly-packed spheres) and a permeability of K = 3.3 × 10−9 m 2 . The wedge angle was θ = 15D . A flow rate of Q = 9.5×10-5 m2/s was used, which produced laminar flow through the pore structure as well as a laminar fluid film. These specifications resulted in a value of A = 250 for the system. The range of parameter A chosen for the analytical solution was established by preserving the validity of the volume-averaging technique used in the formulation as one limit and by keeping the flow as laminar as the other limit. It can be seen that the film thickness increased with the parameter A, which was dependent on the flow rate and the permeability of the medium. Comparison of the numerical and analytical solutions revealed that the film height distributions were practically coincident. Therefore, neglecting terms for a simplified analysis introduced a very small error in the film thickness.
Figure 9.16 Variation of film thickness along a cone.
Chapter 9 Evaporation
723
For the different transport problems considered here, the film thicknesses corresponding to the thin film flow adjacent to a cone are shown in Fig. 9.16. Since the vertex of the cone is a singular point with zero area, the flow is expected to change very rapidly in the vicinity of the vertex. Moreover, the flow at that location will be determined by how the fluid is introduced to the cone (in the form of a liquid jet, etc.). Therefore, the analysis and computations presented here considered flow at some distance from that region so that a thin film could be attained. The most significant parameter for flow over a cone was B . Unlike the parameter A used in relation to the wedge, this parameter contains Γeff, which varied with the type of transport (heat transfer, gas absorption, or solid dissolution). For numerical computation, the cone angle was θ = 15D . The same porous medium and fluid combination was used as that for the computation of flow over a wedge. A flow rate of Q = 5.9×10-6 m3/s was used for the numerical computation. The analytical and numerical methods had negligible difference in their predictions of film thickness. The results and heat and mass transfer can be found in Rahman and Faghri (1993b).
9.5 Direct Contact Evaporation
Direct contact evaporation processes are the most effective cooling strategies when the coolant is mixed with the medium it will cool. Since direct contact evaporation requires no container to separate the coolant from the medium to be cooled, there are no corresponding thermal resistances.
9.5.1 Evaporation of a Liquid Droplet in a Hot Gas
Figure 9.3 describes the situation in which a liquid drop evaporates to a surrounding hot gas. Since the liquid droplet is surrounded by the hot gas, it is assumed that the liquid droplet is spherical in shape and the evaporation is a spherically symmetric problem. Neglecting radiation effect, the energy equation for a liquid droplet is ∂TA α A ∂ § 2 ∂TA · =2 (9.210) ¨r ¸ ∂t r ∂r © ∂r ¹ where the thermophysical properties of the liquid are assumed to be temperatureindependent. The initial temperature of the droplet is TA (r , t ) = Ti , t = 0 (9.211) The symmetric condition at the center of the liquid droplet requires that ∂TA (9.212) = 0, r = 0 ∂r
724
Transport Phenomena in Multiphase Systems
Since the evaporation takes place on the surface of the liquid droplet, its surface temperature equals the saturation temperature corresponding to the partial pressure of the vapor at the interface, i.e., T (r , t ) = Tsat , r = rI (9.213) An energy balance at the surface of the liquid droplet accounts for the latent heat exchanged through conduction inside the drop and convection away from the droplet; that energy balance is dr ∂T ρ A hAv I = kA A − h(T∞ − Tsat ), r = rI (9.214) dt ∂r where h is the heat transfer coefficient between the gas and the liquid droplet. If the liquid droplet is evaporating to the vapor-gas mixture, the mass balance at the interface can be written as dr ρ A I = − ρ g hm (ω I − ω∞ ) (9.215) dt where hm is the mass transfer coefficient, ρ g is the density of the mixture, and ω is the mass fraction of the evaporating component in the mixture. Analytical solution of the problem defined by eq. (9.210) – (9.215) will yield a very complex expression of the final results. Since the liquid droplet is usually very small, one can assume that the heat convected to the droplet surface is balanced by the latent heat of evaporation, so heat conduction to the liquid droplet can be neglected. This assumption is valid because the bulk temperature of the droplet will quickly approach the saturation temperature (Lock, 1994), i.e., h(T∞ − Tsat ) drI =− , r = rI (9.216) ρ A hAv dt which is subject to the following initial condition: rI = Ri , t = 0 (9.217) Integrating eq. (9.216) and considering its initial condition, eq. (9.217), the transient radius of the liquid droplet is obtained. h(T∞ − Tsat )t rI = Ri − (9.218) ρA hAv In arriving at eq. (9.218), the heat transfer coefficient was assumed to be constant throughout the entire evaporation process. The time required to evaporate the droplet completely, tf, can be obtained by setting rI equal to zero in eq. (9.218), i.e., ρA hAv Ri tf = (9.219) h(T∞ − Tsat ) The time required to evaporate the droplet completely can also be estimated from the mass balance equation (9.215). By following a procedure similar to that which obtained eq. (9.219), one obtains ρA Ri (9.220) tf = ρ g hm (ωI − ω∞ )
Chapter 9 Evaporation
725
Combining eqs. (9.219) – (9.220), one obtains h Ja g ω I − ω∞ = hm ρ g c pg
(9.221)
where Jag is the Jakob number, defined as c pg (T∞ − Tsat ) (9.222) Ja g = hAv The above analysis is based on the assumption that conduction to the liquid droplet is negligible. Evaporation of a liquid droplet with consideration of the effect of conduction heat loss will be investigated in the following example 9.8.
Example 9.8 Heat conduction in a liquid droplet can be considered using an integral approximate solution with a parabolic temperature distribution in the liquid droplet. Find the transient radius of the droplet, with conduction in the liquid droplet considered. Solution: Integrating eq. (9.210) in the interval of (0, rI) with respect to r, one obtains rI ∂T ∂T r 2 A dr = α A rI2 A (9.223) 0 ∂t ∂r r = rI
³
The left-hand side of eq. (9.223) can be rewritten using Leibnitz’s rule, i.e., dr ∂T d rI 2 (9.224) r TA dr − rI2Tsat I = α A rI2 A 0 dt dt ∂r r = rI
³
which is the integral equation of the heat conduction problem. Assuming the temperature distribution in the liquid droplet is parabolic profile; and using eqs. (9.212) – (9.213) to eliminate two of the three constants, the temperature distribution in the liquid droplet becomes (Dombrovsky and Sazhin, 2003)
§r· (9.225) TA (r , t ) = Tc (t ) + [Tsat − Tc (t )] ¨ ¸ © rI ¹ where Tc(t) represents the temperature in the center of the liquid droplet, which is still unknown at this point. Substituting eq. (9.225) into eq. (9.224), one obtains § α A 1 drI · d ( ΔT ) (9.226) = −15ΔT ¨ 2 + ¸ dt rI dt ¹ © rI
2
where ΔT = Tsat − Tc (t ) . Substituting eq. (9.225) into eq. (9.214), the energy balance at interface becomes dr 2k (T − T ) ρA hAv I = A sat c − h(T∞ − Tsat ), r = rI dt rI
726
Transport Phenomena in Multiphase Systems
which can be rearranged to º drI 2 kA ª 1 kg Nu (T∞ − Tsat ) » = « ΔT − dt ρA hAv rI ¬ 4 kA ¼ where Nu = 2hrI / k g is the known Nusselt Number.
(9.227)
Equations (9.226) – (9.227) are a set of ordinary differential equations with rI and ΔT as unknowns. They can be solved numerically with the following initial conditions at t = 0 : rI = Ri ΔT = Tsat − Ti Elperin et al. (2005) numerically investigated simultaneous heat and mass transfer during evaporation and condensation on the surface of a stagnant droplet in the presence of admixtures that contain noncondensable solvable gas. The equations for Stefan velocity and the transient radius of the droplet were obtained by using conservation of species mass at droplet surface with absorption and evaporation at the surface accounted for. Their results showed that droplet evaporation rate, droplet temperature, interfacial absorbate concentration, and rate of mass transfer are highly interdependent.
9.5.2 Evaporation of a Liquid Jet in a Pure Vapor
The physical model of the problem under consideration is shown in Fig. 9.17 (Lock, 1994). A liquid jet flows out from a nozzle with a radius of r0 and is surrounded by pure vapor at saturation temperature. At the exit of the nozzle, the velocity and temperature are uniformly T0 and u0, respectively. It is further assumed that the velocity in the jet remains uniformly equal to u0 as it continues flowing; then the jet can be treated as slug flow. This is possible when the friction between the liquid jet and the surrounding vapor is negligible. The temperature of the jet will be affected by the hot gas as soon as it exits the nozzle. Because evaporation occurs on the surface of the liquid jet, its surface temperature is equal to the saturation temperature corresponding to the vapor pressure. The energy equation for the liquid jet is ∂T α ∂ § ∂T · (9.228) u0 A = A ¨ r A ¸ ∂x r ∂r © ∂r ¹ where the thermophysical properties have been assumed to be constants. The initial temperature of the jet is TA (r , x) = T0 , t = 0 (9.229)
Chapter 9 Evaporation
727
u0, T0
2Ri
r x
δt 2 rI
Free-surface Superheated liquid
Vapor
Figure 9.17 Evaporation from a jet surrounded by a hot gas (Lock, 1994; Reproduced by permission from Oxford University Press).
The boundary conditions at the center and surface of the liquid jet are ∂TA (9.230) = 0, r = 0 ∂r T (r , x) = Tsat , r = rI ( x) (9.231) where rI(x) is the radius of the liquid jet, which equals r0 at x =0 but decreases with increasing x. The energy balance at the surface of the liquid jet is dr ∂T ρA hAv u0 I = kA A , r = rI ( x) (9.232) dx ∂r which demonstrates that superheated liquid supplies the latent heat of vaporization. Since the vapor temperature equals the saturation temperature, there is no convection on the surface of the liquid jet. As the liquid jet moves downstream, a thermal boundary layer, r0 − δ t < r < r0 , develops in the liquid jet. The evaporation of the jet can be divided into two stages: the early stage when δ t r0 and a later stage when δ t = r0 . These two stages are similar to the thermal entrance problem and fullydeveloped flow for forced convection in a tube. During the first stage, the effect of the curvature of the cylindrical liquid jet can be neglected because the thermal boundary layer is much thinner than the
728
Transport Phenomena in Multiphase Systems
radius of the tube. The solution of the temperature distribution in the thermal boundary layer can be obtained by using the solution of heat conduction in a semi-infinite solid, i.e.,
§ y T ( y, x) = Tsat + (T0 − Tsat )erf ¨ ¨ 2 α x/u 0 A ©
· ¸ ¸ ¹
(9.233)
where the y-coordinate originates at the jet’s surface and points toward its center. It is related to r by y = rI − r . Substituting eq. (9.233) into the energy balance equation at the interface, eq. (9.232), one obtains k (T − T ) drI (9.234) = − A 0 sat dx ρA hAv πα A u0 x The variation in the radius of the liquid jet during the first stage can be obtained by integrating eq. (9.234), i.e., 2k (T − T ) x (9.235) rI = Ri − A 0 sat ρ A hAv πα A u0 Evaporation during the second stage ( δ t = rI ), which occurs at the lower portion of the jet, is much slower than the first stage and ceases after the superheat in the liquid jet vanishes. The final radius of the liquid jet, rI , f , can be obtained by a simple energy balance: ρA c pAπ Ri2 (T0 − Tsat ) = ρA hAvπ ( Ri2 − rI2, f ) (9.236)
Rearranging eq. (9.236), one obtains the final radius of the liquid jet as rI , f = Ri 1 − Ja (9.237) where Ja = is the Jakob number. c pA (T0 − Tsat ) hAv (9.238)
9.6 Evaporation inside Pores and Slots/Microchannels
Heat transfer during the evaporation of liquid from a liquid/vapor meniscus heated from solid wall plays an important role in predicting the performance characteristics of various two-phase devices with capillary structures such as heat pipes, evaporators, heat sinks, and separators. Heat transfer during the evaporation of liquid from cylindrical and/or hemispherical liquid-vapor menisci in pores or slots is shown in Fig. 9.18. During evaporation from a pore with small wall-liquid temperature drops, the vapor flow does not affect the morphology of the liquid-vapor interface. It can be considered to have nearly-constant curvature except for the microfilm region [see Fig. 9.18(a)]. For extremely high heat flux, on the other hand, the vapor flow can significantly affect the morphology of the
Chapter 9 Evaporation
729
Figure 9.18 Evaporation of liquid from cylindrical pores: (a) low or moderate evaporation rates, and (b) high evaporation rates (Khrustalev and Faghri, 1997).
liquid-vapor interface, as shown in Fig. 9.18(b). Analytical and numerical solution of evaporations inside cylindrical pores under moderate heat fluxes will be considered in Sections 9.6.1 and 9.6.2. Evaporation inside microchannels under high heat flux will be discussed in Section 9.6.3. Section 9.6.4 discusses evaporation in an inclined micro channel.
9.6.1 Evaporation from Cylindrical Pore under Low/Moderate Heat Flux
Evaporation in a microchannel under low/moderate heat flux was studied by Khrustalev and Faghri (1995a). In the scenario under consideration, liquid evaporates from the surface of the menisci situated at the liquid-vapor interface. Figure 9.19 shows the schematic of the cylindrical pore (microchannel) and liquid meniscus. A description of the heat transfer during evaporation from a pore is given here with the two following main assumptions. 1. The temperature of the solid-liquid interface Ts can be considered constant along the s-coordinate for small s. 2. The curvature of the axisymmetrical liquid-vapor surface of the meniscus is defined by the main radius of curvature K = 2 / Rmen and is independent of s. The validity of the second assumption has been proven numerically by Khrustalev and Faghri (1995b); it was shown that this assumption could give an error of less than 5% when calculating the overall heat transfer coefficient during evaporation from a capillary groove. Since heat transfer during evaporation from thin films in a pore is similar to that in a capillary groove, this assumption can be justified for the present analysis.
730
Transport Phenomena in Multiphase Systems
Figure 9.19 Schematic of evaporation from a cylindrical pore
The local heat flux through the liquid film due to heat conduction is T −T ′′ qA = kA s δ (9.239)
δA
where the local thickness of the liquid layer δ A and the temperature of the free liquid film surface Tδ are functions of the s-coordinate. Tδ is affected by the disjoining and capillary pressures. It also depends on the value of the interfacial resistance, which is defined for the case of comparatively small heat flux at the ′′ interface, qδ , by the following relation given by the kinetic theory:
§ 2α · hAv ª pvδ ( psat )δ º ′′ qδ = − ¨ − (9.240) « » ¸ Tδ » © 2 − α ¹ 2π Rg « Tv ¬ ¼ where pvδ and ( psat )δ are the saturation pressures corresponding to Tv and at the thin liquid film interface, respectively. The relation between the saturation vapor pressure over the thin evaporating film, ( psat )δ , which is affected by the disjoining pressure, and the normal saturation pressure corresponding to Tδ , psat (Tδ ), is given by the following relation (see Chapter 5): ª ( p ) − psat (Tδ ) + pd − σ K º ( psat )δ = psat (Tδ ) exp « sat δ (9.241) » ρ A Rg Tδ « » ¬ ¼ which reflects the fact that under the influence of the disjoining and capillary pressures, the liquid free surface saturation pressure ( psat )δ is different from normal saturation pressure psat (Tδ ) . It varies along the thin film (or s-
Chapter 9 Evaporation
731
coordinate) while pvδ and Tv are the same for any value of s. This variation is also due to the fact that Tδ changes along s. While the evaporating film thins as it approaches the point s = 0, the difference between ( psat )δ , given by eq. (9.241), and the pressure obtained for a given Tδ using the saturation table becomes larger. This difference is the reason for the existence of the thin nonevaporating superheated film, which is in equilibrium state in spite of the fact that Tδ > Tsat . ′′ ′′ Under steady-state conditions, qA = qδ , and it follows from eqs. (9.239) and (9.240) that δ § 2α · hAv ª pvδ ( psat )δ º − (9.242) Tδ = Ts + ¨ « » ¸ kA © 2 − α ¹ 2π Rg « Tv Tδ » ¬ ¼ Equations (9.241) and (9.242) determine the interfacial temperature Tδ and pressure ( psat )δ for a given vapor pressure, pvδ , the solid-liquid interface temperature Ts , and the liquid film thickness δ A ( s ) . As the liquid film thins, the disjoining pressure, pd, and the interfacial temperature, Tδ, increase. Under specific conditions, a non-evaporating film thickness is present which gives the equality of the interfacial and solid surface temperatures, Tδ = Ts . This is the thickness of the equilibrium nonevaporating film δ 0 , For the disjoining pressure in water, the following equation for was used in the present analysis (Holm and Goplen, 1979):
ª § δ ·b º pd = ρ A Rg Tδ ln « a¨ ¸» « © 3 .3 ¹ » ¬ ¼
1 ª p (T ) − pv Tw / Tv + σ K § pv Tw · º ½ ° ° δ 0 = 3.3 ® exp « sat w + ln ¨ ¨ p (T ) T ¸ » ¾ ¸ ρA Rg Tw « v ¹» ° °a © sat w ¬ ¼¿ ¯ The total heat flow through a single pore is defined as ª Rp Rp T − T Rp T − T q p = ³ s δ 2π rds = ³ s δ 2π Rmen sin «arctan 0 0 2 2 δ A / kA δ A / kA « s + Rmen − R p ¬
1/ b
(9.243)
where a = 1.5336 and b = 0.0243. From equations (9.241) – (9.243), the following expression for the thickness of the equilibrium film is given:
(9.244)
º » ds » ¼ (9.245) The surface of the pore wall is totally covered with microroughnesses, from which the characteristic size varies, for example, Rr = 10-8 to 10-6 m. Apparently, the thin liquid film formation can be affected by some of these microroughnesses. The following approximation for the liquid film thickness was given by Khrustalev and Faghri (1995b): for δ 0 ≤ δ A ≤ δ 0 + Rr and Rr δ 0 ,
δ A = δ 0 + Rr − Rr2 − s 2 − Rmen + ( Rmen 2 + s 2 + 2 Rmen s sin θ f
)
1/ 2
(9.246)
732
Transport Phenomena in Multiphase Systems
where for a surface with microroughness θ f = 0 , and the liquid film thickness in the interval, δ A ≥ δ 0 + Rr , is
δ = δ 0 + Rr − Rmen + ( Rmen 2 + s 2 + 2 Rmen s sin θ men )
1/ 2
(9.247)
For the smooth-surface model, θ f is the angle between the solid-liquid and liquid-vapor interfaces at the point on s where the disjoining pressure and dK/ds become zero. Note that for small values of the accommodation coefficient (for example α = 0.05), the value of Rr has not affected the total heat flow rate through the liquid film (Khrustalev and Faghri, 1995b). The interfacial radius of curvature is related to the pressure difference between the liquid and the vapor by the extended Laplace-Young equation: ρ 2v 2 § 1 2σ 1· + v 2vδ ¨ − pvδ − pAδ = (9.248) ¸ ε © ρA ρv ¹ Rmen where vvδ is the vapor mean blowing velocity specified for a given meniscus, and ε is the porosity. The latter is needed in this equation because the evaporation takes place into the dry region of the porous structure. The temperature of the saturated vapor near the interface, Tv, is related to its pressure by the saturation conditions: Tv = Tsat ( pvδ ) (9.249) Then the heat transfer coefficient during evaporation from the porous surface is defined as εsqp hep = (9.250) 2 π R p (Ts − Tv ) where ε s = Ap / At is the surface porosity defined as the ratio of the surface of the pores to the total surface of the porous structure for a given cross-section (it is assumed that ε s = ε ; Khrustalev and Faghri, 1995a). Figure 9.20 shows that the heat transfer coefficient during evaporation from the porous surface, he,p,
Figure 9.20 Heat transfer coefficient during evaporation from the porous surface vs. liquid meniscus contact angle (Khrustalev and Faghri, 1995a).
Chapter 9 Evaporation
733
depends significantly on the curvature of the liquid meniscus. For smaller pore sizes, the values of he,p are larger because a larger relative surface is occupied by the thin films.
9.6.2 Fluid Flow Effect in Pore/Slots during Evaporation
Khrustalev and Faghri (1996) investigated the influence of the liquid and vapor flows on the heat transfer in evaporation from a liquid-vapor meniscus, which is presented in this subsection. Evaporation of liquids from capillary grooves, pores, and slots takes place as a result of the solid wall superheat. The temperature of the solid-liquid interface, Tw, is higher than that of the vaporliquid interface. Since the solid wall thermal conductivity is significantly higher than that of the liquid in many technical applications, constant wall temperature is assumed in the present analysis. The solid wall in the vicinity of the interline, which is an apparent liquid-vapor-solid boundary, can be dry or covered with nonevaporating liquid film; in both cases, the temperature of the solid wall or of the nonevaporating liquid free surface is Tw. Since the temperature of the vapor over the evaporating meniscus is lower than Tw, the vapor is affected by the superheated wall and/or liquid such that it also becomes superheated. This means that the temperature of the vapor some distance from the liquid-vapor interface can be higher than Tsat ( pv ) . Note that for comparatively thick liquid films ( δ ≥ 1 m ), the liquid-vapor interface temperature is very close to the vapor temperature. This is true because the interfacial thermal resistance is small compared to that of the liquid films. Based on this consideration, the momentum conservation and energy equations for the evaporating liquid and vapor over the liquid-vapor interface are included in the present model of the evaporating extended meniscus. The evaporation of liquid from the liquid-vapor meniscus interface in a narrow slot as considered here is illustrated by Fig. 9.21(a). Gravitational body forces are negligible in comparison to surface tension forces, due to the small size of the meniscus. It is assumed that the liquid-vapor interface curvature is constant: K = 1/Rmen. This assumption can be justified by the results of Swanson and Peterson (1994), where one-sided formulation of the heat transfer and flow in liquid was used. Actually, the curvature of the liquid-vapor interface of the evaporating non-isothermal meniscus is not exactly constant, especially in the microfilm region. However, the assumption of constant curvature is sufficient to model the heat transfer and flow in the vapor phase over the evaporating meniscus. The evaporating microfilm is the thin-film part of the meniscus where the disjoining pressure is significant and where the temperature of the liquidvapor interface differs from that of the vapor phase over the given point of the interface. This microfilm is considered using the approach in Section 5.5.4, while the fluid flow and heat transfer in the major part of the domain shown in Fig. 9.21(a) are treated by two-dimensional analysis.
734
Transport Phenomena in Multiphase Systems
(a) Computational domain (b) Control volume that contains interface
Figure 9.21 Evaporating liquid-vapor meniscus in a slot (Khrustalev and Faghri, 1996).
For an incompressible fluid with constant vapor and liquid viscosities and temperature dependent vapor density, the continuity and momentum equations for both the vapor and liquid are: ∂u ∂v + =0 (9.251) ∂x ∂y
ρ ¨u
§ ∂ 2u ∂ 2u · § ∂u ∂u · ∂p +v ¸=− + μ¨ 2 + 2 ¸ ∂y ¹ ∂x ∂y ¹ © ∂x © ∂x
(9.252) (9.253)
§ ∂ 2v ∂ 2v · § ∂v ∂v · ∂p +v ¸=− + μ¨ 2 + 2 ¸ ∂y ¹ ∂y ∂y ¹ © ∂x © ∂x The energy equation for both the vapor and liquid is § ∂ 2T ∂ 2T · § ∂T ∂T · +v = k¨ 2 + 2 ¸ ρcp ¨ u ¸ ∂y ¹ ∂y ¹ © ∂x © ∂x
ρ ¨u
(9.254)
The wall temperature is assumed constant (except for the control volume enclosing the microfilm). T x = 0 = Tw (9.255) The inlet flow temperature is known. T y = 0 = Tw Nonslip conditions at the wall are assumed. u x=0 = v x=0 = 0 (9.256) (9.257)
Chapter 9 Evaporation
735
At the inlet, the velocity profile is assumed to be fully developed: u y =0 = 0
(9.258)
3mA ,in x (2 X L − x ) (9.259) 3 2 ρA X L At the symmetry line, x = X L , u = 0, dv / dx = 0, and dT / dx = 0 . Heat transfer at the interface was considered through analysis of the mass and energy balance for a rectangular control volume containing a liquid-vapor interface [see Fig. 9.21 (b)]. The microfilm region, where the temperature of the liquid-vapor interface changes from Tw at the interline, which is the apparent liquid-vapor-solid boundary, to Tsat at the end of the microfilm, has been enclosed into a special control volume with the dimensions Δx × Δy . This is the vapor control volume since the microfilm thickness, δ , is much smaller than Δx . This means that in the two-dimensional analysis, the liquid microfilm is represented only by the vapor blowing velocity at the solid wall, Vv , mic , and by the temperature of the wall, Tmic, (this is the mean temperature of the microfilm free surface) in the special control volume. Note that for the liquid control volume located under the microfilm control volume [see Fig. 9.21 (b)], the velocity of the vapor at the upper liquid-vapor face is increased in the numerical procedure by vA = ( ρv / ρA ) V v ,mic Δy / Δx in order to satisfy conservation of v y =0 =
mass. The values Vv , mic and Tmic are the microfilm characteristics averaged over y along its length, which does not exceed y.
′′ Tmic = ³ qδ ( y )Tδ ( y ) dy
0
Δy
³
Δy
0
′′ qδ ( y ) dy
(9.260) (9.261)
′′ where qδ is the heat flux at the interface, which is equal to that due to heat conduction through the liquid film; eq. (9.240) defines that for the case of a comparatively small heat flux. Note that the interfacial thermal resistance is very important for the microfilm and results in the difference between the interface temperature, Tδ , and that of the vapor over the interface, Tsat. For the rest of the meniscus, the interfacial thermal resistance is negligible in comparison to that of the liquid film, and Tδ ~ Tsat, which is justified by calculations made at the end of the microfilm. The effective evaporative heat transfer coefficient characterizing the intensity of the evaporative heat transfer is defined as follows: mA ,in hAv he = (9.262) X L (Tw − Tsat ) where Tsat is the vapor temperature over the liquid-vapor interface, which is assumed constant due to saturation conditions and because the curvature of the interface was constant.
′′ V v , mic = ³ qδ ( y ) dy
0
Δy
( ΔyhAv ρ v )
736
Transport Phenomena in Multiphase Systems
(a) ΔT = 5 K
(b) ΔT = 10 K
(a) ΔT = 5 K
(b) ΔT = 10 K
Figure 9.22 Vapor flow over evaporating meniscus (Khrustalev and Faghri, 1996).
Figure 9.23 Liquid flow in the evaporating meniscus (Khrustalev and Faghri, 1996).
The problem was solved using the SIMPLE algorithm (Patankar, 1980) employing the power law discretization scheme and a line-by-line TDMA method. All control volumes in the computational domain were divided into two categories: “liquid” and “vapor,” with corresponding thermophysical properties. The following nonequality was used to identify the vapor control volumes: 2 x > X L − [ Rmen − ( y − Ymen − Rmen ) 2 ] (9.263) Numerical results were obtained for water at the saturation temperature of 373 K for X L = 50 m, YL = 400 m , YB = 100 m , and θ men = 15D . The accommodation coefficient, α , was set equal to unity. Typical vapor velocity fields over the evaporating meniscus are shown in Fig. 9.22 for comparatively small (a) and large (b) temperature drops ( ΔT = Tw − Tsat ) between the wall and the saturation temperature. The vapor velocity in the microfilm control volume can be as high as 100 m/s for ΔT = 15 K . Such a high-blowing velocity results in a recirculation zone in the vapor near the wall over the interline, as shown in Fig. 9.22(b). Since the velocities in the liquid are much smaller than those in the vapor, they are not shown in Fig. 9.22. For small temperature drops, the vapor flows along the x-coordinate near the heated wall and over the interline. However, for larger temperature drops, the vapor flow changes its direction in the above-mentioned region due to appearance of the recirculation zone. This effect can be important for micro heat pipes and other applications. The liquid flow is shown in Fig. 9.23. The liquid-vapor interface in Fig. 9.23 is represented by the velocity vectors at the vapor side of the interface. The liquid flow pattern did not change with temperature drop, as can be seen from comparison of Fig. 9.23 (a) and 9.23 (b).
Chapter 9 Evaporation
737
Figure 9.24 Effective heat transfer coefficient versus superheat of water at atmospheric pressure (Khrustalev and Faghri, 1996).
The numerical results for the effective heat transfer coefficient he were compared to those obtained using the one-dimensional model (see Section 9.6.1) by Khrustalev and Faghri (1995a). For small ΔT the agreement was very good: the one-dimensional model underestimated the heat transfer coefficients by only about 3%, as shown in Fig. 9.24. For larger temperature drops, however, the difference between the present and the one-dimensional model became 30 % (for T = 10 °C) due to the convective heat transfer.
9.6.3 Evaporation under High Heat Flux
At high heat flux, high velocity vapor flow may create an extended thick liquid film attached to the evaporating meniscus in spite of the capillary pressure drop between the hemispherical meniscus and the nearly flat thick film. Khrustalev and Faghri (1997) analyzed the evaporation of a pure liquid in the vicinity of a hemispherical liquid-vapor meniscus formed within a circular micropore, as shown in Fig. 9.18(b). Since the radius of the pore and the liquid film thickness are of the same order of magnitude, it is necessary to describe the problem in a cylindrical coordinate system, as shown in Fig. 9.25. Both the vapor and the liquid flow along the z-direction. The origin of the coordinate system is fixed with respect to the minimum film thickness. It is assumed that the vapor is at saturation condition. Heat transfer in the liquid film is in the radial direction by
738
Transport Phenomena in Multiphase Systems
Figure 9.25 Schematic of extended thick liquid film in a circular microchannel (Khrustalev and Faghri, 1997).
conduction only, and convective effects are neglected. The velocity profile in the fluid region is assumed to be fully developed and to behave as a modified laminar-type flow. This accounts for surface tension and mass transfer effects through the boundary conditions, but otherwise a single, axial direction momentum equation completely describes the flow field in the fluid. Combining expressions for the conservation of mass and the conservation of energy in the liquid film, one obtains R 1§ q· (9.264) ³R −δ rwA ( r ) dr = 2πρ A ¨ mA ,in − hAv ¸ © ¹ where mA ,in is the incoming liquid mass flow rate at z = 0 in the liquid film layer and q(z) is the heat flow rate through a given cross-section due to evaporation in the film from z = 0 to z , i.e.,
′′ q = 2π R ³ qR dz
0
z
(9.265)
For conduction heat transfer through the liquid film, Tw − Tδ ′′ qR = kA R ln [ R /( R − δ )]
(9.266)
Drawn from the definition of total heat flow rate, an expression for the gradient of heat flow rate in the z-direction is Tw − Tδ dq = 2π kA (9.267) ln [ R /( R − δ ) ] dz
Chapter 9 Evaporation
739
The differential equation of the velocity profile from the conservation of momentum in the film is 1 ∂ § ∂wA · 1 dpA (9.268) ¨r ¸= r ∂r © ∂r ¹ μA dz The boundary conditions for this equation are the nonslip condition at the wall and equality of shear stress at the liquid-vapor interface: (9.269) r = R , wA = 0
∂wA dσ dTδ · 1§ f = ¨ − v ρ v wv2 − =E (9.270) ∂r dT dz ¸ μA © 2 ¹ where f v is the vapor fraction factor. The solution to eq. (9.268) with eqs. (9.269) and (9.270) as boundary conditions, gives 2 ( R − δ ) ln r º + E R − δ ln r 1 dpA ª 1 2 2 (9.271) wA = − « (R − r ) + » ( ) 2 R» R μ A dz « 4 ¬ ¼ Substituting eq. (9.271) into eq. (9.264), the following equation for liquid pressure in terms of the mass flow rate and heat flow rate is obtained: ª 1 §q º · − mA ,in ¸ + E ( R − δ ) F » μA « ¨ dpA ¹ ¬ 2πρ A © hAv ¼ = (9.272) 2 2 2 ·º dz ª R 4 ( R − δ ) § ( R − δ ) − R ¸» « + ¨F + ¨ 2 8 4 ¸» « 16 © ¹¼ ¬ where r = R −δ ,
F=
(R − δ )2 § ln
2 ¨ ©
1 · R2 R + ¸− R −δ 2 ¹ 4
(9.273)
The pressure difference in the vapor and liquid phases due to capillary and disjoining pressure is 3 − 2 ½ ª § dδ ·2 º 2 1 dδ · ° °d δ § pv − pA = σ ® 2 «1 + ¨ cos ¨ arctan (9.274) ¸» + ¸ ¾ − pd R −δ dz ¹ ° © ° dz « © dz ¹ » ¬ ¼ ¯ ¿ where the two terms in the braces represent the two radii of curvature in the microchannel. Equation (9.274) can be rewritten as the following two ordinary differential equations: dδ =δ′ (9.275) dz 3/ 2 ª p − p + p cos(arctan δ ′) º dδ ′ A d = (1 + δ ′2 ) « v − (9.276) » σ dz R −δ ¬ ¼
740
Transport Phenomena in Multiphase Systems
The vapor flow is assumed to be compressible and quasi-one-dimensional (Faghri, 1995): dpv 1 ªd = ( − β v ρv wv2 Av ) − f v ρ wv2 ( R − δ ) dz Av « dz ¬ (9.277) dδ ·º § 2 −2π ( R − δ ) ρ v vv ,δ sin ¨ arctan ¸» dz ¹ ¼ ©
where β v = 1.33 for small Reynolds numbers. evaporation is dq 1 vv ,δ = dz 2π ( R − δ ) ρ v hAv The blowing velocity due to (9.278)
The vapor is assumed to behave as an ideal gas.
pv R g Tv The density gradient in the z-direction is then dρ v 1 § dp v 1 pv dTv · ¸ ¨ = − dz R g ¨ dz Tv Tv2 dz ¸ ¹ ©
ρv =
(9.279)
(9.280)
The temperature and pressure at saturation are related by the ClausiusClapeyron equation – 2 dTv dpv Rg Tv = (9.281) dz dz pv hAv At this point, seven expressions – eqs. (9.267), (9.272), (9.276), (9.277), (9.278), (9.280), and (9.281) – involving seven unknown variables – q , δ ,
δ ′ , pA , pv , ρv , and Tv – completely describe the stated assumptions.
Boundary conditions are required for each of the seven first-order differential equations at z = 0: q=0 (9.282) 2σ (9.283) pA = pv ,in − + pd R − δ in
pv = pv ,in = psat (Tv ,in )
δ = δ in δ′ = 0
(9.284) (9.285) (9.286) (9.287) (9.288)
ρ v ,in =
pv ,in Rg Tv ,in
Tv = Tv ,in
Chapter 9 Evaporation
741
As can be seen from Fig. 9.25, the liquid film ends with a microfilm; therefore, its thickness at z = Lδ equals the nonevaporating film thickness 0, determined by eq. (9.244). The interfacial temperature Tδ is required to obtain the gradient of the heat transfer rate along the z-direction, as indicated by eq. (9.266). It can be determined by solving eqs. (9.241) and (9.242) simultaneously. The mass flow rate of the liquid at z = 0, mA ,in , can be determined by the following boundary condition at the end of the microfilm: q( Lδ ) (9.289) mA ,in = hAv The problem described by seven ordinary differential equations – eqs. (9.267), (9.272), (9.276), (9.277), (9.278), (9.280) and (9.281) – with corresponding boundary conditions – eqs. (9.282) – (9.288) – can be solved using a Runge-Kutta method in conjunction with the shooting method to satisfy the constitutive equation, eq. (9.286). More detailed information about the numerical procedure can be found in Khrustalev and Faghri (1997).
Figure 9.26 Variation of characteristics of the liquid-vapor interface along the cylindrical channel: (a) thickness of the liquid film, and (b) pressure in the vapor and liquid phases (Khrustalev and Faghri, 1997).
742
Transport Phenomena in Multiphase Systems
Figure 9.27 Variation of characteristics of the liquid-vapor interface along the cylindrical channel (a) heat flux at the wall, and (b) temperature of the vapor and liquid-vapor interface (Khrustalev and Faghri, 1997).
Numerical solution is performed for water evaporation from a cylindrical pore with an inner diameter of 20 m and wall temperature of 388 K. The accommodation coefficient is set equal to unity. Fig. 9.26(a) shows the variation in liquid film thickness along the z-direction. The liquid film thickness initially increases slowly along the z-coordinate, and then rapidly decreases approaching the microfilm region at the end of the film. The film is concave at z=0, where it attaches to the hemispherical meniscus; however, it is convex over most of its length. This is due to the pressure gradient in the liquid and vapor along the film. Figure 9.26(b) demonstrates that the vapor pressure is decreasing and the liquid pressure is increasing along the z-coordinate. The liquid flows in the film due mainly to frictional vapor-liquid interaction at the interface. Figure 9.27(a) demonstrates that the rate of evaporation decreases towards the thickest part of the liquid film and reaches its maximum in the microfilm region. However, the microfilm region at the end of the liquid film does not make any significant contribution to the total mass flow rate at the pore outlet, because
Chapter 9 Evaporation
743
(a) the length of the microfilm region is much smaller than in the case where a hemispherical meniscus ends directly with a microfilm; and (b) disjoining pressure causes the liquid-vapor interfacial temperature Tδ to sharply increase at the end of the film, as shown in Fig. 9.27 (b). Khrustalev and Faghri’s (1997) work demonstrated, for the first time, that high-vapor velocities during the evaporation of pure liquids in micropores could allow a thick liquid film to attach to a hemispherical meniscus. In addition to the evaporation of pure liquid in a circular channel discussed in this section, Coquard et al. (2005) investigated evaporation in a capillary tube with a square cross-section. The evaporation rate is much higher than in a circular tube because of the liquid flow along the corner induced by capillary force. They focused on slow evaporation controlled by mass transfer and identified three regimes: the capillary regime, capillary-viscous regime, and capillary-gravity regime.
9.6.4 Evaporation in an Inclined Microchannel
The analysis in the preceding subsection dealt with evaporation of a liquid to pure vapor when the evaporation was driven by heat transfer. Chakraborty and Som (2005) analytically investigated heat transfer in an evaporating thin liquid film moving slowly along the walls of an inclined microchannel, and this work will be introduced in this subsection.
Quiescent air, T∞, ω∞
Channel centerline
Meniscus Liquid
Figure 9.28 Evaporation in an inclined microchannel (Chakraborty and Som, 2005; Reprinted with permission from Elsevier).
744
Transport Phenomena in Multiphase Systems
Figure 9.28 shows the lower half of the microchannel, in which quiescent air with a temperature T∞ and mass fraction ω∞ is in the core of the channel while liquid film forms on the microchannel wall. It is assumed that the transport phenomena are two-dimensional and that the diffusion in the axial direction is negligible. The temperature distribution in the liquid film is assumed to be linear along the y-direction. The momentum equation that governs the liquid flow is ∂ 2u ∂p μA 2 = A − ρA g sin θ = − F ( x) (9.290) ∂y ∂x and is subject to the following boundary conditions u = 0 at y = 0 (9.291) ∂u (9.292) = 0 at y = δ ∂y Integrating eq. (9.290) twice and considering eqs. (9.291) and (9.292) yields 1 (9.293) u= F ( x)(2δ y − y 2 ) 2μA The energy equation in the vapor phase is ∂T ∂ 2T Vδ v = α v 2v (9.294) ∂y ∂y where Vs is the Stefan velocity due to evaporation at the liquid surface. Equation (9.294) is subject to the following boundary conditions Tv = Tδ at y = δ (9.295) Tv = T∞ at y = h (9.296) The solution of eq. (9.294) with boundary conditions eqs. (9.295) and (9.296) is Tv − Tδ exp(Vδ y / α v ) − exp(Vδ δ / α v ) = (9.297) Tδ − T∞ exp(Vδ δ / α v ) − exp(Vδ h / α v ) Similarly, the mass fraction of the vapor is ω − ωδ exp(Vδ y / Dv ) − exp(Vδ δ / Dv ) = (9.298) ωδ − ω∞ exp(Vδ δ / Dv ) − exp(Vδ h / Dv ) where Dv is the mass diffusivity of the vapor in the air. The mass balance at the interface is dΓ = − ρvVδ (9.299) dx where Γ = ³ ρ A udy is the mass flow rate per unit width. Considering the
0
δ
velocity profile in eq. (9.293), the mass balance becomes dδ = −3ν A ρ vVδ δ 3 F ′( x) + 3δ 2 F ( x) dx
(9.300)
Chapter 9 Evaporation
745
Assuming the interface is impermeable to nonevaporating species, the mass balance at the interface can be written as ∂ω − Dv ∂y y =δ = −Vδ (9.301) 1 − ωδ Substituting eq. (9.298) into eq. (9.301) yields § 1 − ω∞ · D (9.302) Vδ = v ln ¨ ¸ h − δ © 1 − ωδ ¹
The energy balance at the interface is ∂T ∂T − kv v + kA A = − ρv hAvVδ ∂y y =δ ∂y y =δ
(9.303)
Substituting eq. (9.297) into eq. (9.303) and considering the linear temperature profile in the liquid film, the interfacial temperature becomes − ρ h V + kATw / δ − T∞ f Tδ = v Av δ (9.304) kA / δ − f where V exp(Vδ δ / α v ) f = kv δ (9.305) α v exp(Vδ δ / α v ) − exp(Vδ h / α v ) The pressure in the liquid film is A σ + (9.306) pA = p∞ − R ( x) δ 3 where p∞ is the vapor pressure, R( x) is the radius of curvature at x, and the last term on the right-hand side accounts for the effect of disjoining pressure.
Figure 9.29 Local Nusselt number (Chakraborty and Som, 2005; Reprinted with permission from Elsevier).
746
Transport Phenomena in Multiphase Systems
Therefore, the function F ( x) in eq. (9.290) becomes σ dR 3 A d δ + + ρ A g sin θ (9.307) F ( x) = − 2 R dx δ 4 dx Substituting eqs. (9.302) and (9.304) into eq. (9.300) and considering the variation of ωδ , Tδ , and hAv , a fourth order ordinary differential equation of δ can be obtained. The boundary conditions of the ordinary differential equation include δ = δ 0 at x = 0 and δ = d δ / dx = d 2δ / dx 2 = 0 at x → ∞ . The liquid film thickness can be obtained by a fourth order Runge-Kutta method. The local Nusselt number is then obtained by hx x Nu x = x = (9.308) kA δ Figure 9.29 shows the variation of the local Nusselt number along the axial direction for different δ 0 / h . It can be seen that the Nusselt number increases with increasing x and eventually reaches a very high value. This represents the location where dryout occurs and the liquid film becomes very small. The Nusselt number increases with decreasing initial film thickness δ 0 / h , but its effect becomes insignificant for δ 0 / h < 0.01 .
9.7 Evaporation from Inverted Meniscus in Porous Media
Evaporators that are capable of withstanding high-heat fluxes, for example larger than 100 W/cm2, are of great interest for electronic component cooling systems. The most promising evaporator design is the so-called “inverted meniscus” type. To predict the critical heat flux and effective heat transfer coefficients in the evaporator, a mathematical model has been developed by Khrustalev and Faghri (1995a). The model includes the following interconnected problems, which are treated simultaneously in the frames of the numerical analysis: (1) heat transfer during evaporation from a pore, (2) heat transfer and vapor flow in the dry region of a porous structure with a stable side boundary, the location of which depends on the operational conditions, and (3) heat conduction in a solid fin (or wall) with a non-uniform heat sink on the side surfaces. Heat transfer during evaporation from a pore has been presented in Section 9.6.1. Heat transfer and vapor flow in the dry region, as well as heat conduction in the solid fin, will be presented in this subsection. Schematics of two configurations for the characteristic elements of inverted meniscus evaporators are shown in Fig. 9.30. In the first configuration [Figs. 9.30 (a) and (b)], a heated triangular fin is inserted into the porous plate and sintered with it in order to provide good thermal contact. In the second configuration [Fig. 9.30 (c)], the heated wall is flat. With small heat flux, evaporation of the liquid, which saturates the porous element, can take place exclusively from the surface of the porous body into the vapor channel as shown in Fig. 9.30(a). Extremely high heat fluxes are significantly more interesting for industrial applications, and
Chapter 9 Evaporation
747
Figure 9.30 Schematic of the modeled elements of the inverted meniscus evaporators: (a) with triangular fin for low heat fluxes, (b) with triangular fin for high heat flux, and (c) with flat heated wall for high heat fluxes (A is shown in Fig. 9.18; Khrustalev and Faghri, 1995a; Reprinted with permission from Elsevier).
the existence of a stable vapor blanket along the heated solid surface inside the uniform porous structure was anticipated, as shown in Figs. 9.30(b) and (c). Evaporation takes place into the dry region of the porous structure at the liquidvapor interface, the location of which shifts depending on the operational conditions. Heat from the heated surface is conducted through the porous element’s dry region to the interface. The varying curvatures of the menisci along
748
Transport Phenomena in Multiphase Systems
the liquid-vapor interface create a capillary pressure gradient, which in turn drives the vapor flow. While the vapor flow takes place in a comparatively narrow porous passage, the liquid with the same total mass and flow rate (steady state) is filtered perpendicularly through the entire porous element to the liquidvapor interface. The pressure gradient in the liquid along this interface is negligible compared to that in vapor. This assumption can be justified by the case where the maximum pressure drop in liquid over the wetted region of the characteristic element under consideration is negligible compared to the pressure drop in vapor in the dry region. This assumption allows for one-dimensional approximation description of the heat transfer in the vapor blanket. The thickness of the vapor blanket, δ v , increases with the heat flux, which can lead to increasing thermal resistance in the element. In the situation where v|x=0 = 0 is of the same order of magnitude as the minimum thickness of the porous element [see Fig. 9.30 (b)], the vapor can penetrate into the liquid channels and obstruct the liquid supply of the evaporator, eventually resulting in dryout. The value of the heat flux at which the dryout takes place can be considered as critical. The operating parameters of the evaporator depend upon the heat and mass circulation in the entire system (for example in a heat pipe), with the evaporator in consideration. In the present model, the physical situation for the characteristic element is determined by three parameters: (1) the pressure in the liquid near the interface, pAδ , (2) the temperature of the solid surface, T0 , at x = 0, and (3) the liquid-vapor meniscus radius at the end of the vapor blanket ( x = Lvδ ), Rmen,o. Note that the superheat of the fin exists at the following condition: (9.309) T0 > Tsat ( pAδ + 2σ / Rmen ,min ) where the subscript “sat” denotes the normal saturation temperature corresponding to a pressure ( pAδ + 2σ / Rmen ,min ) and Rmen, min = R p / cos θ men, min (Rp is pore radius). Inequality (9.309) characterizes how the value of the solidliquid superheat, Rmen,o (curvature radius of meniscus at outlet), is related to the fluid circulation in the entire device. For the case of the evaporator with the forced liquid supply it can be set Rmen,o Rmen, min , because the pressure drop in the liquid is not due to the capillary pressure in this case. Since it can be anticipated that the temperature drops in metallic fin or wall are much smaller than those across the dry zone of the porous structure, because kw keff, the heat conduction in the solid fin or wall is considered using a l-D approach. For the case of the flat wall [Fig. 9.30(c)] it means that dTw/dy is not included in the consideration. The heat conduction in the triangular metallic fin is described by the following equation [Fig. 9.30(b)], which was obtained as a result of energy balance on a differential element consideration: keff cos γ d 2Tw dTw 1 + + (Ts − Tw ) =0 (9.310) 2 dx dx x xδ v ( x)kw sin γ
Chapter 9 Evaporation
749
where Ts is the local temperature of the porous structure at the liquid-vapor interface location. Similarly, the heat conduction equation for the wall in Fig. 9.30 (c) is keff d 2Tw q′′ + (Ts − Tw ) + 0 =0 (9.311) 2 dx twδ v ( x)kw tw kw The boundary conditions for equations (9.310) and (9.311) are Tw = T0 , x = 0 (9.312) dTw = 0, x = 0 (9.313) dx ′′ For the second configuration, q0 is the uniform heat flux at the outer surface of ′′ the heated part of the flat wall. The value of q0 and the functions δ v ( x) and Ts ( x) should be given by the results of the vapor flow and heat transfer in the dry region solution considered below. The local heat flux resulting from heat conduction across the dry region of the porous structure from the solid surface to the liquid-vapor interface where evaporation takes place is T −T ′′ qloc = keff w s (9.314) δ v ( x) which is valid for the case kv keff and c p ,v (Tw − Ts ) hAv . Hence, the mean velocity of the vapor flow for a given x along the solid surface is (the mass and energy conservation balances) keff x xT −T 1 w s ′′ uv ( x ) = qloc ( x) dx = (9.315) ³0 ³0 δ v ( x) dx δ v ( x)hAv ρv δ v ( x) hAv ρ v where uv ( x) is the mean vapor velocity along the x-coordinate. The modified Darcy’s equations for the vapor flow in both directions through a porous structure where the value of 0.55 is used for a dimensionless form-drag constant are ∂pv μ 0.55 = − v uv ( x ) − (9.316) ρv uv2 ( x) ∂x K K ∂pv μ 0.55 = − v vv ( y ) − (9.317) ρv vv2 ( y ) ∂y K K where uv and vv are the area-averaged vapor velocities. The corresponding continuity equation is ∂uv ∂vv + =0 (9.318) ∂x ∂y The Darcy’s equation is semi-empirical and describes the flow with the uniform velocity profile; therefore, it is assumed in the present analysis that uv does not depend on y. Taking the definitions of the mean vapor pressure δ δ ( pv = ³0 v pv dy δ v ) and axial velocity ( uv = ³0 v uv dy δ v ) for a given x into consideration and integrating eq. (9.316) over y, the following equation can be
750
Transport Phenomena in Multiphase Systems
obtained for the gradient of the mean vapor pressure along the x-coordinate dpv μ 0.55 = − v uv ( x ) − ρv uv2 ( x) (9.319) dx K K Since the situation when δ v Lvb is considered, the vapor pressure drop across the vapor blanket is much smaller than that along the x-coordinate. At the solid fin (or wall) surface vv y = 0 = 0 , and at the liquid-vapor interface
vv
y =δ v
= vvδ ζ where ζ = cos[arctan(d δ v dx)] is the cosine of the angle between
the y coordinate and the normal to the liquid-vapor interface and vvδ is the blowing velocity (normal to the liquid-vapor interface) : T −T (9.320) vvδ = keff w s δ v hAv ρ v Equation (9.320) implies that the total amount of energy transferred from the heated solid surface to the liquid-vapor interface by the heat conduction across the dry porous zone is spent on vaporization of the liquid. Since the axial velocity profile is nearly uniform, it follows from equation (9.318) that vv = vvδ ζ y / δ v . Integrating eq. (9.317) twice over y for a given x and implementing the definition of pv the difference between the vapor pressure near the liquid-vapor interface, pvδ , and the mean vapor pressure of the vapor flow, pv , for a given x can be estimated as follows 0.55 § v ζμ · ρ v vv2δ ζ 2 ¸ (9.321) pvδ − pv = δ v ¨ vδ v + 4K © 3K ¹ Combining eqs. (9.315) and (9.319), we have for the vapor filtration flow pressure gradient along the x-coordinate: 0.55 ª keff x Tw − Ts º dx » (9.322) « 0 δv ρv K ¬ δ v hAv 0 δ v ¼ The boundary condition for eq. (9.322) follows from eqs. (9.248), (9.320) and (9.321) ρv2 vv2δ x =0 § 1 1 · 2σ pv x = 0 = pAδ + + ¨−¸ Rmen x = 0 ε2 © ρA ρv ¹ (9.323) keff μvζ (Tw − Ts ) x = 0 0.55 2 2 (ζ vvδ δ v ) − − x =0 3KhAv ρv 4K The local heat flux at the liquid-vapor interface due to the evaporation of the liquid is: ′′ qloc = [Ts ( x) − Tv ( x)]he , p (9.324)
vv keff dpv =− dx δ v hAv K
³
x
Tw − Ts
2
dx −
³
Combining eqs. (9.314) and (9.324) because of the steady state situation, the expression for the local temperature of the porous structure at the liquid-vapor interface location is:
Chapter 9 Evaporation
751
Ts ( x) =
Tw ( x) + he , p ( x)Tv ( x)δ v ( x) / keff 1 + he , p ( x)δ v ( x) / keff
(9.325)
Substituting eqs. (9.320) and (9.321) into eq. (9.248) and differentiating it, the following equation for the radius of the meniscus curvature can be obtained d § 2σ · dpv 2 ρ v vvδ § 1 1 · keff μv ª§ dTw dTs · − − δv ¨ ¸= ¨−¸ 2 dx © Rmen ¹ dx ε © ρA ρv ¹ hAvδ v2 «¨ dx dx ¸ ¹ ¬© −(Tw − Ts ) keff μv dδ v º » + 3Kh ρ dx ¼ Av v ª § dTw dTs «ζ ¨ dx − dx ¬© dζ º · ¸ + (Tw − Ts ) dx » ¹ ¼
(9.326) 2 0.55 ρv § keff · 2 ª dTw dTs · dδ v º § + − δ v − (Tw − Ts ) 2 ¨ ¸ ®ζ « 2(Tw − Ts ) ¨ dx dx ¸ dx » 4 K © hAv ρvδ v ¹ ¯ ¬ © ¹ ¼ dζ ½ +2ζδ v (Tw − Ts ) 2 ¾ dx ¿ with the boundary condition Rmen x = 0 = C0 (9.327) where C0 should be chosen from the constitutive condition for the minimum value of the meniscus radius along the liquid-vapor interface min { Rmen ( x)} = R p / cosθ men ,min (9.328)
Now, the condition of the liquid-vapor interface mechanical equilibrium should be considered in order to find its location or δ v ( x) . In the analysis by Solov’ev and Kovalev (1987) it was assumed that δ v ( x) = const ⋅ x 0.33 , which is not quite satisfactory for several reasons. For example, in a hypothetical situation in which the definite start point along the x-coordinate is x1, and there is no evaporation from the liquid-vapor interface, δ v x > x = const , a condition which is
1
not satisfied by the discussed expression. Wulz and Embacher (1990) have modeled the vapor flow in the uniform zone of the dry porous structure. The ′′ thickness of the vapor zone was determined as 0.1 mm at qmax = 17500 W/m2 by comparing the calculated temperature difference between the fin top and the phase boundary with the value determined by the experiment using the simpler vapor zone model. Chung and Catton (1993) have considered the problem of steam injection into a slow water flow through porous media, where the interface location was also unknown. They have found “…that the interface can be idealized as a stream line as far as the momentum equations are concerned.” Here the concept of a streamline is used indirectly as explained below. For the asymptotic case, K → ∞, μv → 0, integrating Euler's equation along a streamline gives ρv uv2 ρv vv2 + = pv x = 0 (9.329) pv + 2 2 where the terms containing uv2 and vv2 correspond to the inertia effects due to
752
Transport Phenomena in Multiphase Systems
acceleration of fluid. In the present analysis, the vapor flow through a porous medium is described by Darcy’s momentum equation. However, it is assumed that since the velocity profile of the vapor flow along the x-coordinate, uv , is nearly uniform, eq. (9.329) can be used for the description of the inertia effects at the liquid-vapor interface due to acceleration of the vapor. The liquid-vapor interface can be stable provided it has the shape that eliminates the influence of the inertia effects due to acceleration of the vapor flow on the vapor pressure near this interface. While the steady-state situation is analyzed, the liquid pressure along the interface is constant, and the pressure losses in the vapor flow in both directions due to friction and solid obstacles are compensated by the capillary pressure, the vapor pressure gradient along the stable interface due to these inertia effects should be equal to zero. Since the velocity profile of the vapor flow along the x-coordinate is nearly uniform, it follows from eq. (9.329) that ρv uv2 ρ v vv2 + = const (9.330) 2 2 Note that eq. (9.330) is not used for the fluid flow in the porous medium but describes the inertia effects at the adjustable liquid-vapor interface in which the momentum equations for the vapor flow in the porous medium are concerned. Equation (9.330) is necessary in order to find the equilibrium location of the liquid-vapor boundary. Substituting eqs. (9.315) and (9.320) into eq. (9.330) and differentiating it after some rearrangements gives the equation for the vapor blanket thickness, δ v
d δ v2 dx 2 dδ v ° § 2 ®δ v (Tw − Ts ) ζ sin ¨ arctan dx © ° ¯
2 −1 ½ · ª § dδ v · º ° ¸ «1 + ¨ dx ¸ » ¾ ¹« © ¹» ° ¬ ¼¿
ª = (Tw − Ts ) « ¬ dδ −v dx ª§ «¨ «© ¬
³
x
Tw − Ts
0
δv δv
dT · º § dT dx + δ vζ 2 ¨ w − s ¸ » dx ¹ ¼ © dx
(9.331)
³
δv
Tw − Ts
0
2 º · dx ¸ + ζ 2 (Tw − Ts ) 2 » » ¹ ¼
where all of the terms containing (Tw − Ts ) can be calculated in the numerical procedure using the functions Tw(x) and Ts(x) determined at the previous iteration. The second-order differential equation (9.331) should be solved with the two boundary conditions for the variables δ v and d δ v / dx . The first boundary condition is δ v x =0 = C1 (9.332) where C1 should be chosen from the constitutive condition as the value of δ v x =0 , which satisfies the following boundary condition
Rmen
x = Lvb
= Rmen ,o
(9.333)
The second boundary condition is due to the symmetry of the considered element
Chapter 9 Evaporation
753
(Fig. 9.29). Since at the point x = 0, dTw/dx = 0, dTs/dx = 0, and dvvδ / dx = 0 due to the physical reasons, it follows from equation (9.320) dδ v =0 (9.334) dx x =0 Thus we have six main variables (or unknown functions): pvδ , pv , Rmen, Ts, vv and v which should be found from the six equations: (9.248), (9.320), (9.322), (9.325), (9.326), and (9.331). These six equations should be solved along with those presented in the previous sections for variables he , p ( x) and ′′ Tw(x). Note that the value q0 which is needed for eq. (9.311) can now be found as: T ( x) − Ts ( x) 1 Lvb ′′ q0 = keff w dx (9.335) δ v ( x) W0 ′′ For the first configuration, q0 is the heat flux in the solid fin corresponding to the porous structure-vapor channel plane. The heat flux on the outer surface of the evaporator (and the corresponding effective heat transfer coefficient) can be recalculated taking the geometry of the evaporator into consideration. Although the vapor leaving the dry zone of the porous structure is superheated, it is convenient to relate the local effective heat transfer coefficient to the vapor
³
Figure 9.31 Performance characteristics of the modeled evaporator element along the heated fin surface: (a) temperatures of the fin surface and of the porous structure at the liquid-vapor interface and (b) local heat flux across the vapor blanket (Khrustalev and Faghri, 1995a; Reprinted with permission from Elsevier).
754
Transport Phenomena in Multiphase Systems
saturation temperature, because c p ,v (Tw − Ts ) hAv . Thus, the local effective heat transfer coefficient corresponding to the point x = Lvb (outlet of the vapor flow) is defined as: Lvb T ( x) − Ts ( x) 1 heff = keff w dx (9.336) 0 δ v ( x) W (Tw − Tv )o Khrustalev and Faghri (1995a) modeled heat transfer during evaporation from a pore (see Section 9.6.2), heat conduction in the solid fin or wall, vapor flow and heat transfer in the dry region of the porous structure. Numerical solution was performed and the results were obtained with constant thermal physical properties corresponding to the saturation temperature Tsat ( pAδ ) = 100 o C . The numerical results were obtained for the first
³
configuration [Fig. 9.30(b)] for the case of a miniature evaporator: γ = 30D , t pen = 0.2 mm, R p = 20 m, θ men, min = 33D , α = 0.05, keff = 10 W/m-K, k w = 438 W/m-K, ε = 0.5, ε s = 0.5, K = 0.5 × 10−12 m 2 , and pAδ = 1.013 × 105 Pa. Since the longitudinal circulation of the fluid in the heat pipe was not considered, the numerical results were obtained for several fixed Rmen,o. The temperatures of the fin surface and the porous structure at the liquid-vapor interface, as well as local heat flux across the vapor blanket, are shown in Fig. 9.31. It can be seen that the temperature drops at the solid heated surface were significantly smaller than those corresponding to the porous skeleton at the liquid-vapor interface. The real superheat of the liquid, Ts − Tv , which could initiate the boiling, was significantly smaller than the superheat of the heated solid surface: Ts − Tv < Tw − Tv . The local heat fluxes across the dry zone had their maximums at the point x = 0, as can be seen in Fig. 9.31(b). For extremely high heat fluxes, ql′′ , when the temperature drop ( Ts − Tv ) is oc large, boiling of the liquid at the liquid-vapor interface can occur. This can cause instabilities at the liquid-vapor interface. However, boiling of the liquid does not ′′ ′′ necessarily result in dryout for an evaporator of this type. Since qo > qloc due to Lvb > W , the triangular geometry of the solid fin helps to postpone boiling; this may make the fin more desirable than the flat wall in the second configuration shown in Fig. 9.30 (c). In other words, the triangular geometry of the fin provides a higher value for the heat flux on the outer surface of the evaporator, which corresponds to the beginning of liquid boiling at the liquid-vapor interface.
References
Alhusseini, K.A., Tuzla, K., and Chen, J.C., 1998, “Falling Film Evaporation of Single Component Liquids,” ASME Journal of Heat Transfer, Vol. 41, pp. 16231632.
Chapter 9 Evaporation
755
Carey, V.P., 1992, Liquid-Vapor Phase-Change Phenomena: An Introduction to the Thermophysics of Vaporization and Condensation Processes in Heat Transfer Equipment, Hemisphere Publishing Corp., Washington, D. C. Chakraborty, S., and Som, S.K., 2005, “Heat Transfer in an Evaporating Thin Liquid Film Moving Slowly Along the Walls of an Inclined Microchannel,” International Journal of Heat Mass Transfer, Vol. 45, pp. 2801-2805. Chun, K. R. and Seban, R. A., 1971, “Heat transfer to evaporating liquid films,” ASME Journal of Heat Transfer, Vol. 93, pp. 391-396. Chung, M., and Catton, I., 1993, “Steam Injection into a Slow Water Flow through Porous Media,” ASME Journal of Heat Transfer, Vol. 115, pp. 734-743 Coquard, T., Camassel, B., and Prat, M., 2005, “Evaporation in Capillary Tubers of Squared Cross Section,” Proceedings of 2005 ASME Summer Heat Transfer Conference, San Francisco, CA. Debbissi, C., Orfi, J., and Nasrallah, S.B., 2003, “Evaporation of Water by Free or Mixed Convection into Humid Air and Superheated Steam,” International Journal of Heat and Mass Transfer, Vol. 46, pp. 4703-4715. Dombrovsky, L.A., and Sazhin, S.S., 2003, “A Simplified Non-Isothermal Model for Droplet Heating and Evaporation,” International Communication on Heat and Mass Transfer, Vol. 30, pp. 787-796. Elperin, T., Fominykh, A., and Krasovitov, B., 2005, “Modeling of Simultaneous Gas Absorption and Evaporation of Large Droplet,” Proceedings of 2005 International Mechanical Engineering Congress and Exposition, Orlando, FL, Nov. 5-11, 2005. Faghri, A., 1995, Heat Pipe Science and Technology, Taylor and Francis, Washington, D.C. Faghri, A. and Seban, R., 1985, “Heat Transfer in Wavy Liquid Films,” International Journal of Heat Mass Transfer, Vol. 28, pp. 506-508. Hewitt, G.F., Shires, G.L., and Bott, T.R., 1994, Process Heat Transfer, Begell House, New York, NY. Hirschburg, R.F., and Florschuetz, L.W., 1982, “Laminar Wavy Film Flow, Parts I and II,” ASME Journal of Heat Transfer, Vol. 104, pp. 452-464. Holm, F.W., and Goplen, S.P., 1979, “Heat Transfer in the Meniscus Thin-Film Transition Region,” ASME Journal of Heat Transfer, Vol. 101, No. 3, pp. 543547. Kaptiza, P.L., and Kaptiza, K.P., 1975, “Wavy flow of Thin Layers of Viscous Fluid,” Collected Papers of P.L. Kapitza, Vol. 2, pp. 662-709. Pergamon, New York.
756
Transport Phenomena in Multiphase Systems
Kays, W.M., Crawford, M.E., and Weigand, B., 2004, Convective Heat Transfer, 4th ed., McGraw-Hill, New York, NY. Khrustalev, D. and Faghri, A., 1995a, “Heat Transfer in the Inverted Meniscus Type Evaporator at High Heat Fluxes,” International Journal of Heat and Mass Transfer, Vol. 38, pp. 3091-3101. Khrustalev, D.K. and Faghri, A., 1995b, “Heat Transfer during Evaporation and Condensation on Capillary-Grooved Structures of Heat Pipes,” ASME Journal of Heat Transfer, Vol. 117, August, No. 3, pp. 740-747. Khrustalev, D.K. and Faghri, A., 1996, “Fluid Flow Effects in Evaporation from Liquid/Vapor Meniscus,” ASME Journal of Heat Transfer, Vol. 118, pp. 725730. Khrustalev, D., and Faghri, A., 1997, “Thick Film Phenomenon in High Heat Flux Evaporation from Cylindrical Pores,” ASME Journal of Heat Transfer, Vol. 119, No. 2, pp. 272-278. Kutateladze, S.S., 1963, Fundamentals of Heat Transfer, Academic Press Inc., New York. Kutateladze, S.S., 1982, “Semi-Empirical Theory of Film Condensation of Pure Vapors,” International Journal of Heat Mass Transfer, Vol. 25, pp. 653-660. Lock, G.S.H., 1994, Latent Heat Transfer, Oxford Science Publications, Oxford University, Oxford, UK. Patankar, S.V, 1980, Numerical Heat Transfer and Fluid Flow, McGraw-Hill, New York. Rahman, M.M., and Faghri, A., 1992, “Analysis of Heating and Evaporation from a Liquid Film Adjacent to a Horizontal Rotating Disk,” International Journal of Heat and Mass Transfer, Vol. 35, pp. 2644-2655. Rahman, M.M., and Faghri, A., 1993a, “Transport in a Thin Liquid Film on the Outer Surface of a Wedge or Cone Embedded in a Porous Medium, Part I: Mathematical Analysis,” International Communications in Heat and Mass Transfer, Vol. 20, pp. 15-27. Rahman, M.M., and Faghri, A., 1993b, “Transport in a Thin Liquid Film on the Outer Surface of a Wedge or Cone Embedded in a Porous Medium, Part II: Computation and Comparison of Results,” International Communications in Heat and Mass Transfer, Vol. 20, pp. 29-42. Rice, J., Faghri, A., and Cetegen, B.M., 2005, “Analysis of a Free Surface Flow for a Controlled Liquid Impinging Jet over a Rotating Disk Including Conjugate Effects, with and without Evaporation,” International Journal of Heat and Mass Transfer, Vol. 48, pp. 5192-5204.
Chapter 9 Evaporation
757
Rogovan, I.A., Olevski, V.M., and Runova, N.G., 1969, “Measurement of the Parameters of Film Type Wavy Flow on a Vertical Plane,” Thm. Found. Chem. Eng., Vol. 3, p. 164. Seban, R. A., and Faghri, A., 1976, “Evaporation and Heating with Turbulent Falling Liquid Films,” ASME Journal of Heat Transfer, Vol. 98, pp. 315-318. Solov’ev, S.L. and Kovalev, S.A., 1987, “Heat Transfer and Hydrodynamics in the Inverted Meniscus Evaporator of a Heat Pipe,” Proceedings of the 6th International Heat Pipe Conference, Grenoble, France, Vol. I, pp. 116-120. Stephan, K., 1992, Heat Transfer in Condensation and Boiling, Springer-Verlag, New York. Swanson, L.W., and Peterson, G.P., 1994, “Evaporating Extended Meniscus in a V-Shaped Channel,” AIAA Journal of Thermophysics and Heat Transfer, Vol. 8, No. 1, pp. 172-180. Wulz, H., and Embacher, E., 1990, “Capillary Pumped Loops for Space Applications Experimental and Theoretical Studies on the Performance of Capillary Evaporator Designs,” Proceedings of the AIAA/ASME 5th Joint Thermophysics and Heat Transfer Conference, Seattle, WA, AIAA Paper No. 90-1739. Yan, W.M., Lin, T.F., and Tsay, Y.L., 1991, “Evaporative Cooling of a Liquid Film through Interfacial Heat and Mass Transfer in a Vertical Channel – I. Experimental Study,” International Journal of Heat and Mass Transfer, Vol. 34, pp. 1105-1111. Yan, W.M., and Lin, T.F., 1991, “Evaporative Cooling of a Liquid Film through Interfacial Heat and Mass Transfer in a Vertical Channel – II. Numerical Study,” International Journal of Heat and Mass Transfer, Vol. 34, pp. 1113-1124. Yang, W.J., Zhang, N., Chai, A., Guo, K.H., and Sakamoto, M., 1997, “Evaporation induced Benard Convection – A New Type and its Mechanism,” ASME HTD-Vol. 349, pp. 37-50, The 32nd National Heat Transfer Conference, August 8-12, 1997, Baltimore, MD.
Problems
9.1. A liquid drop is on the top of a substrate at temperature Ts. Suppose the liquid drop is disc-like and the initial thickness is δ 0 . Find the time required to completely evaporate the liquid drop. Assume that heat transfer across the disc-like liquid drop is one-dimensional conduction. A small amount of water is contained in a tray, as shown in Fig. P9.1. The bottom of the tray is maintained at a temperature, Tw, which is above the
9.2.
758
Transport Phenomena in Multiphase Systems
ambient temperature, T∞. The initial bulk temperature of the water in the tray, Tb, is between Tw and T∞. The rate of evaporation at the surface of the water is mevp (kg/s). Find the bulk temperature of the water as a function of time. m T∞ Tb Tw
Figure P9.1
9.3.
A thin layer of liquid film sits on top of a heated surface at temperature T0 (see Fig. P9.2). The top surface of the liquid film is at temperature Tδ ( Tδ < T0 ). The mass flux on the top surface of the liquid film due to ′′ evaporation is mδ . Yang et al. (1997) suggested that the liquid layer could be divided into upper and lower regions with linear temperature distribution in each region. Build an appropriate model to determine the location and temperature of the boundary between the two regions.
Figure P9.2
9.4.
Solve the Stefan’s diffusion problem (see Example 9.2) in terms of molar concentration by assuming the molar fractions of species A at liquid/vapor interface, x A,δ , and at the top of the tube, x A,H , are known. The mixture concentration can be assumed constant. Do the Stefan diffusion problem (Problem 9.4) as quasi-steady rather than steady condition and obtain the change in height with time in terms of molar concentration of mixture and initial height H.
9.5.
Chapter 9 Evaporation
759
9.6.
Evaporation from a horizontal thin film is described by three ordinary differential equations (9.14) – (9.16) and the corresponding boundary conditions, eqs. (9.18) – (9.21). Write a computer program to solve the set of ordinary differential equations and discuss the effects of various parameters, such as Pr, Sc, and Ja on Nu x / Re1/ 2 . x Evaporation from a horizontal thin liquid film is described by eqs. (9.1) – (9.6), (9.8) and (9.10). Exact solution of the problem using similarity variables was introduced in subsection 9.3.1. To solve the same problem using integral method, the profiles of velocity, temperature, and concentration need to be determined. If these profiles can be assumed to be second-degree polynomial functions, find velocity, temperature, and concentration distribution in the boundary layers. Solve the evaporation problem using integral approximate method. For the sake of simplicity, assume Pr=Sc=1 so that the velocity, thermal, and concentration boundary layer thicknesses are the same. Use the velocity, temperature, and concentration distributions in the boundary layers obtained from the previous problem. The blowing velocity at the surface due to evaporation can be neglected. Sweat cooling is a technique used to protect an intensely-heated porous wall by evaporating a liquid injected from the pores to the wall surface. In order to cool a 0.5-m long porous stainless surface, water is injected to the surface at a rate that is just sufficient to keep the surface wet. Dry air at 600K and 1 atm flows over the plate at a velocity of 10 m/s. Find the rate of water supply needed to maintain the wall temperature at 400 K.
9.7.
9.8.
9.9.
9.10. Falling water film evaporation occurs on the inner sides of two parallel plates with a length of 1 m, and they are separated by a distance of 2b = 3 cm . The inlet and outlet liquid film temperatures are 50 ÛC and 45 ÛC, respectively. Find the mass flow rate of the liquid film per unit width,
Γ.
9.11. Air enters the bottom of a vertical channel formed by two parallel walls (see Fig. P9.3). A liquid water film is present on the left wall, which is maintained at a constant temperature. The right wall is adiabatic and has no liquid film on it. The liquid film on the left wall is so thin that its thickness can be neglected. The left wall can be considered as a surface where the concentration of the vapor equals the saturation concentration of the vapor at its partial pressure. Assuming that the air flow in the vertical channel is laminar and the air inlet velocity is uniform, find the governing equations and boundary conditions to describe this conjugate heat and mass transfer process.
760
Transport Phenomena in Multiphase Systems
Figure P9.3
9.12. Air flows into a tube with radius R and constant wall temperature Tw. Evaporation occurs on a thin liquid film on the inner wall of the tube. The thickness of the liquid film can be neglected in comparison with the radius of the tube. The inner wall can be considered as a surface where the concentration of the vapor equals the saturation concentration of the vapor at its partial pressure. The air flow at the inlet of the tube is fully developed, but the temperature and concentration distributions at the inlet are uniform. Specify the governing equations and corresponding boundary conditions of the problem. 9.13. Solve the energy and species equations of the above problem to obtain the temperature and concentration distributions in the tube. The Schmidt number, Sc, can be assumed to be equal to one. Obtain the local convective heat transfer coefficient for the system. 9.14. A two-phase thermosyphon, shown in Fig. P9.4, is needed to deliver a steady-state heat rate of Q = 1kW at an operating temperature of Tv = 333 K. An acceptable wall superheat is T = 5 K. The inside diameter of the circular channel is D = 0.1 m, and the working fluid is water. Use classical Nusselt analysis to determine how long the evaporator should be to avoid dryout. 9.15. Derive eq. (9.125) from eqs. (9.122) – transformation of variables ( x, y ) → (ξ ,η ) η = y / δ ( x) . (9.123), using the where ξ = x and
9.16. Show that eq. (9.121) can be transformed into eq. (9.135) by using the variables defined in eqs. (9.124), (9.127), and (9.133) – (9.134).
Chapter 9 Evaporation
761
Figure P9.4
9.17. In a falling liquid film evaporation experiment, the temperature of the heating wall is 106 °C and the vapor pressure is at one atmosphere. The inlet mass flow rate of the liquid film per unit width is Γ 0 = 2.24 × 10−2 kg/s-m . Is the film flow laminar, laminar with waves, or turbulent? If the height of the vertical plate is 0.1 m, determine the average heat transfer coefficient. Please also check the overall energy balance by using eq. (9.109). 9.18. A uniformly-distributed liquid water film at 0.06 kg/s is introduced at the top of a tube with an inner diameter of 4 cm and height of 3 m. The inner surface temperature of the tube is uniformly 105 °C, and the saturation temperature of the steam is 100 °C. Is the liquid film at the inlet and outlet laminar, laminar with waves, or turbulent? Find the evaporation rate and the average heat transfer coefficient. 9.19. Falling film evaporation on a vertical wall is analyzed in great detail in Section 9.4. How would you apply the solutions of falling film evaporation on a vertical heated wall to predict heat transfer of falling film evaporation on an inclined heating wall? The inclination angle between the heating wall and the vertical direction is θ .
762
Transport Phenomena in Multiphase Systems
9.20. You are asked to select the most weight-efficient evaporant for a particular evaporator. The temperature at which the evaporator will operate has not been set. Go through Appendix B.4 and rank each working fluid based on its heat of vaporization ( hAv ). Also, for each working fluid, note what the operating temperature of the evaporator would be if the operating pressure was 1 atm. Which evaporants stand out in the different temperature ranges? 9.21. Consider evaporation for a constant temperature wall rotating disk with controlled liquid jet impingement at its center (see Fig. P9.5). Rotational forces cause the film to spread across the disk and then off the outer edge. The flow field is similar to a laminar, steady falling film, except that rotation replaces gravity as the driving force and the system must be solved in cylindrical coordinates. Obtain the heat transfer coefficient for the thin liquid film evaporation on the rotating disk by neglecting inertia and convective terms.
z
r
Figure P9.5
9.22. For evaporation over a wedge embedded in porous media, as discussed in Section 9.4.5, the film thickness is not a function of x. How would you simplify the transport eq. (9.198)? 9.23. A spray of liquid water droplets with mean diameter d = 50 m at a temperature 373 K is injected into a tank of superheated water vapor at a pressure p = 1.013× 105 Pa held at constant T = 423 K. The convective heat transfer coefficient on the surface of the droplet is 4000 W/m2-K. If the nozzle is 0.2 m from the nearest surface of the tank, what is the
Chapter 9 Evaporation
763
maximum velocity the droplets can have without making contact with the wall? 9.24. The analysis presented in Section 9.5.1 assumed that the heat transfer coefficient on the surface of the liquid droplet was a constant. When a liquid droplet is spread into the hot gas, the following empirical correlation can be used to estimate the heat transfer coefficient: Nu = 2 + 0.6 Re1/ 2 Pr1/ 3 . Since the radius of the droplet is changing in the evaporation process, the heat transfer coefficient is not a constant. How will you modify eq. (9.219) when heat transfer coefficient is not a constant? 9.25. Slurry droplets and hot air are brought into contact in a slurry-drying process. During the early stage, the slurry droplet can be considered as a porous core covered by pure water (see Fig P9.6). Evaporation takes place on the surface of the slurry droplet, the temperature of which can be assumed to be uniform all the time. The first phase of evaporation is referred to as evaporation of the pure water around the porous core. Build an appropriate model to predict the transient radius and temperature of the slurry droplet during the first phase of evaporation.
Figure P9.6
9.26. A stationary liquid fuel droplet with an initial temperature of T0 is introduced into hot gas at a temperature of T∞ . The initial pressure in the liquid droplet and the hot gas are p0 and p∞ respectively. The initial mass fraction of the fuel vapor in the hot gas is ω∞ . Evaporation takes place on the surface of the liquid fuel droplet, which is assumed to be in spherical shape throughout the process. Since the temperature of the fuel droplet is below the temperature of the hot gas, a laminar downward motion of hot gas is formed near the liquid fuel droplet. Write the governing equations and corresponding boundary conditions for liquid and gas flow, and heat and mass transfer in both liquid fuel droplet and hot gas. The spherical coordinate system with origin located at the center of the liquid fuel droplet should be used.
764
Transport Phenomena in Multiphase Systems