Subscriber access provided by UNIV OF LOUISIANA
Thermodynamics, Transport, and Fluid Mechanics
A semi-analytical method for computing the dynamics of bubble growth: the effect of superheat and operating pressure Jyoti Bhati, and Swapan Paruya Ind. Eng. Chem. Res., Just Accepted Manuscript • DOI: 10.1021/acs.iecr.8b02231 • Publication Date (Web): 05 Oct 2018 Downloaded from http://pubs.acs.org on October 5, 2018
Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.
is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.
Page 1 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Industrial & Engineering Chemistry Research
A semi-analytical method for computing the dynamics of bubble growth: the effect of superheat and operating pressure Jyoti Bhati, Swapan Paruya* Department of Chemical Engineering National Institute of Technology Durgapur, West Bengal-713209, India ABSTRACT In this work, we propose a semi-analytical method for computing the dynamics of a sphericalbubble growth in an infinite pool of superheated water. The dynamics shows the temporal variation of bubble radius in superheat of 1.4K-17.89K. The method accurately computes the bubble growth controlled by surface tension, particularly in very low superheats, in which the existing theories are not adequate. The present method consistent with the existing theories computes the bubble growth controlled by heat transfer in high superheats. The semi-analytical method with a modification of heat transfer at the interface is capable of predicting the bubble growth in pool boiling experiments and numerical simulations. The effect of increasing operating pressure on the bubble growth decreases the growth constant following a power law. The variations of energy components due to pressure difference, inertia, surface tension and thermal diffusion required reveal some interesting insights of the bubble-growth mechanism. KEYWORDS: Bubble growth; heat transfer; numerical analysis; phase change; pool boiling 1. INTRODUCTION Insight of growth history of a vapour bubble is essential in estimating the interfacial heat flux and wall heat flux in nucleate boiling. The nucleate boiling is very fast heat transfer process and is very common in many applications including heat exchangers, chemical reactors, thermal power plants, nuclear reactors, rocket engines, electronic cooling, etc. Bubble growth is due to the evaporation of microlayer at bubble surface. The complex dynamical features of the growth
*
[email protected] ACS Paragon Plus Environment
Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
process, detachments, bubble merger, the collapse and the associated interfacial instabilities are also highly interesting and challenging research issues (Westwater1, Dhir2, Tomar et al.3, Uguz and Narayanan4, Ganguli et al.5). Exhaustive experimental and numerical investigations on pool boiling dynamics and collapse dynamics of a single vapour bubble reveals the complexity of bubble dynamics and the transport processes in nucleate boiling. Thus, the research on the bubble dynamics continues to receive serious attentions to clearly understand the mechanism of nucleate boiling. The growth occurs as a result of evaporation around the bubble surface as well as at the bubble base (Dhir2). Tomar et al.3 concluded based on their numerical simulation that high superheat prevented the film rapture due to Rayleigh-Taylor instability at initial stage and eventually caused instability that facilitates bubble formation. Surface tension plays an important role in evaporation. Uruz and Narayanan4 demonstrated that the evaporation was destabilized, when heated from the liquid side, if the surface tension gradient due to the temperature variation leading to the Marangoni convection is not considered. Nevertheless, the pioneering work of Lord Rayleigh6 on vapour bubble collapse also inspires subsequent development of a number of theories on vapour bubble growth in superheated liquid. The Rayleigh theory on bubble collapse due to the pressure increase, which included only the effect of inertia, was used many authors to compute bubble radius and growth rate. Many other authors reported the growth history of a vapour bubble. The problem of the growth of a spherical vapour bubble was addressed by many authors2,7-27. The studies assumed a uniform pressure inside the bubble. Plesset7 modified the Rayleigh equation to include the effect of viscosity and surface tension for describing the growth history of a vapour bubble. The modified equation is known as Rayleigh-Plesset equation. A major assumption of “thin thermal boundary layer” of superheated liquid around the bubble by Plesset-Zwick8-9 simplified the physics of bubble growth to a great extent. The approximations of zero viscosity and linear change of vapour pressure with temperature readily solved RayleighPlesset equation for the bubble growth rate. The small thickness of the layer results from very low thermal diffusivity of the liquid. They also approximated the asymptotic (series) solutions of bubble growth rate to the leading terms. Their theoretical bubble radius agreed well with the experimental result of Dergarabedian0, particularly in high superheats that favor the growth dominated by heat diffusion. Forster and Zuber formulated integro-differential equation for the growth of a vapour bubble in a superheated liquid and identified two distinct time domains - one,
ACS Paragon Plus Environment
Page 2 of 39
Page 3 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Industrial & Engineering Chemistry Research
the time domain of the order of 10-4 s during, in which the effect of the hydrodynamic forces becomes important; two, the time domain in which this effect is unimportant. Until the Scriven theory12, the effect of mass convection was ignored. Scriven first proposed a unified theory based on the thermodynamic principle and first introduced the effect of radial convection of mass in continuity equation. Theofanous et al.4 computed bubble growth in constant and time-dependent pressure fields taking into account the effects of surface tension, inertia of the liquid, heat conduction in the liquid, non-equilibrium at the liquid-vapour interface, and the varying density of the vapour bubble. Mikic et al.6 derived the two limiting solutions - inertia controlled growth (Rayleigh solution) at small time scale and heat diffusion controlled growth (Plesset-Zwick equation8) at large time scale; they predicted the growth rate of spherical bubble in a uniformly superheated liquid by balancing the work done on the liquid and the kinetic energy along with the use of Clausius – Clapeyron equation. Board and Duffey17 extended the Plesset-Zwick theory to determine the condition that exhibits a similar growth-profile with time for water and sodium; They numerically solved four equations - Rayleigh equation, Clausius-Clapeyron equations, energy balance at bubble wall and ideal gas law. Prosperetti and Plesset8 numerically established that the approximation of a thin thermal boundary layer to describe the spherical vapour bubble growth was valid in high superheats which allow the inertia effect to be dominant at the initial stage. Their approximation of linear variation of vapour pressure with temperature leading to a scaling law can adequately characterize the growth phenomena in its latter stage in which the thermal effect dominates. Lee and Merte9 presented a numerical procedure for the solution of vapour bubble growth with radial symmetry from the thermodynamic critical size in an initially uniformly superheated liquid, which includes the influences of surface tension, liquid inertia and heat diffusion. They examined the effect of the disturbance on the solution required to initiate the growth. Robinson and Judd20 eliminated the “thin thermal boundary layer” assumption of Plesset-Zwick by using finite difference scheme with dense gridding to solve the heat diffusion equation for the liquid surrounding the bubble to calculate temperature gradient
( ) 𝑑𝑇
𝑑𝑟 𝑟 = 𝑅.
Kudryashov and Sinelshchikov22 derived a general analytical solution of Rayleigh
equation for both cases - empty hyperspherical vapour bubble and gas–filled hyperspherical vapour bubble. The growth, collapse and collapse motion of the empty hyperspherical bubble were analyzed. The effect of the surface tension on the bubble motion was included in the
ACS Paragon Plus Environment
Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
corresponding general solution of the Rayleigh equation in the three-dimensional space. Sharipov23 made a semi-analytical approximation to the solution of the radial Fourier equation describing liquid temperature dynamics in the vicinity of a spherical bubble. This approximation opens a possibility to construct a computationally efficient bubble model that is flexible enough to simulate different bubble dynamics behavior like bubble growth, collapse, and oscillations. The scenarios like strong bubble parameter oscillations in highly subcooled water and abrupt liquid pressure change were considered. Bhati et al.24 reported some preliminary results of the bubble growth in superheated liquid obtained by the approximation to Plesset-Zwick theory. Their approximation requires further validation of their computational results. The authors25-27 reported the parametric effects on the interfacial instability of the bubble growth and recommended the suitable operating conditions for the growth processes. From the survey of literatures, we note that the studies on bubble growth in superheated water mostly rely on Rayleigh equation and Plesset-Zwick theory8. The theory neglected the effect of vapour-density change with temperature, vapour-pressure variation with time and approximated the asymptotic (series) solutions to the leading terms for bubble growth radius, which essentially suppresses the growth rate. The theory is not adequate for describing the bubble growth in the range of low superheats where the effect of surface tension dominates; the other existing theories significantly deviate from the experiments too. A computationally-efficient and unified method for computing reasonably-accurate bubble growth in low and high superheats is essential. In the present work, we propose a semi-analytical method for computing the history of bubble growth with time in an infinite pool of superheated water and compare the method with the analytical methods of Plesset-Zwick8, Mikic et al.6 and Prosperetti and Plesset8. Prosperetti and Plesset made the assumptions similar to those of Plesset and Zwick. Although our semi-analytical method relies on “thin boundary layer” assumption, it significantly differs from the analytical methods. It makes a leading-order approximation in temperature profile and includes the dependency of phase densities and surface tension on temperature and of variable pressure difference across the bubble wall under varying vapour temperature computed by the energy balance equation. The semi-analytical method determines the temperature profile and the temperature gradient analytically at vapour-water interface and computes heat transfer across the interface, while it numerically solves the equations for the dynamics of the bubble. Our leadingorder approximation in the temperature profile of the surrounding liquid results in a linear
ACS Paragon Plus Environment
Page 4 of 39
Page 5 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Industrial & Engineering Chemistry Research
temperature profile in the “thin boundary layer”. The approximation leads to a significant improvement of predicting bubble growth valid in low superheats and stronger inertia controlled growth in low pressures. We have considered five test cases of bubble growth in superheats from 1.4K to 17.89K for water vapour-bubble growth at a constant system pressure of 1.0 atm and the numerical results of the test cases have been compared with the experiments of Dergarabedian0 and Abdelmessih23. On the other hand, our method has been extended to explain the bubble growth in the pool boiling experiments of Dhir2, Mukherjee and Dhir26 and Akiyama et al.28. One test case of bubble growth at a very low pressure of 0.395 atm, experimented by Board and Duffey17, has been studied. We have also studied the effect of pressure on the bubble growth and examined the role of various energy components in bubble growth. 2. MODEL DESCRIPTIONS We describe the model equations based on the schematic of a growing vapour bubble in an infinite pool of superheated water shown in Figure . R(t) is the bubble radius at time t. The bubble surface, interface and bubble wall are indicated by r=R(t).
𝑑𝑅 𝑑𝑡
is the growth rate of bubble
that produces liquid inertia force. A long radial distance of 0.05 m from the bubble is shown to be R. T is the temperature of superheated water and p or p is the operating or ambient pressure. The surface tension r=R at bubble wall and pressure drop p (=pv-p) across the interface are indicated in Figure . Initially, the effect of surface tension dominates. When the vapour bubble is growing, the effect of surface tension gradually decreases.
ACS Paragon Plus Environment
Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 6 of 39
Figure 1. Schematic diagram of vapour bubble in superheated liquid 2.1 Numerical models Model equation for bubble growth: We make the following simplifications before describing the model equations: (1) Spherical geometry of the bubble regardless its size (2) Thermal equilibrium condition inside bubble and at bubble wall (vapour-liquid interface) (3) Negligible viscous effect of water and liquid sodium (Prosperetti and Plesset8) Modified Rayleigh-Plesset equation used by Plesset-Zwick8 is given by equation (1) for describing the growth of a spherical vapour-bubble in an infinite pool of superheated liquid. The equation computes R(t). The LHS of equation (1) is the pressure difference across the bubble wall. 𝑃𝑣(𝑇𝑣) - 𝑃∞ 𝜌𝑙
𝑑2𝑅
3 𝑑𝑅 2 𝑑𝑡
( )
= 𝑅 𝑑𝑡2 + 2
2𝜎
()
+ 𝜌𝑙𝑅
ACS Paragon Plus Environment
Page 7 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Industrial & Engineering Chemistry Research
Initial conditions are: R(0)=Rc
|
𝑑𝑅 𝑑𝑡 (0)
(2a) (2b)
=0
Young-Laplace equation given by equation (2c) calculates the critical radius Rc for equation (2a). 2𝜎
(2c)
𝑅𝑐 = 𝑃𝑠𝑎𝑡(𝑇∞) ― 𝑃∞ Energy balance at vapour-liquid interface:
The energy required for evaporating the liquid is supplied by thermal diffusion through the liquid. Uniform temperature Tv from the centre of bubble to the bubble wall is assumed. The degree of superheat in the present study is not enough, for which we ignored the effect of the Marangoni convection in interfacial heat transfer. The resistance due to Marangoni force on the bubble is negligible as well. For the spherical bubble, the energy balance at vapour-liquid interface is given by equation (3a):
(
) = ∫𝐴 - 𝑘𝑙(∂𝑇∂𝑟)𝑟 = 𝑅𝑑𝐴
𝑑 4 3 𝑑𝑡 3𝜋𝑅 𝜌𝑣h𝑓𝑔
(3a)
𝑠
As is the surface area of bubble. Considering ρv = ρv(Tv), the initial condition for solving equation (3a) is: (3b)
𝑇𝑣(0) = 𝑇∞
Energy balance in the surrounding liquid: The temperature gradient at the vapour–liquid interface is obtained analytically by solving the one–dimensional energy equation for the moving liquid in spherical coordinate given by following equation. ∂𝑇 ∂𝑡
∂𝑇
( (𝑟 ))
+ 𝑢 ∂𝑟 = 𝛼𝑙
1∂
2∂𝑟
𝑟
2∂𝑇 ∂𝑟
(4)
u is the velocity of bubble wall and is defined by
𝑑𝑅 𝑅 2 𝑑𝑡 𝑟 .
()
r is the radial distance from the centre
of the bubble. The Eulerian coordinate (r, t) is chosen with origin r=0 at the center of the bubble. Initial and boundary conditions for solving equation (4) are: Initial condition: 𝑇(𝑟, 0) = 𝑇∞ (5a)
ACS Paragon Plus Environment
Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 8 of 39
Boundary conditions: T(R, t) = Tv
(5b)
T(R∞, t = T∞
(5c)
2.2 Solution for temperature gradient While solving equation (4) for the temperature gradient at r=R(t), it is advantageous to rewrite equation (4) in Lagrangian coordinate (Prosperetti and Plesset8) to ‘formally’ eliminate the convective term (the second term on the LHS of equation (4)). We also assume that the liquidmotion possesses a spherical symmetry. The Lagrangian coordinate h is defined by: 1
h = 3(𝑟3 - 𝑅3(𝑡))
We assume
𝛿 𝑅
(6)
to be very less than unity because of very small thickness of thermal boundary
layer compared to R. Thus, one can write the condition: h 𝑅3
𝛿
(7)
≈𝑅≪1
Therefore, equation (4) reduces to the following form:
( )
∂𝑇
∂2𝑇
4 ∂𝑡 = 𝛼𝑙𝑅
(8)
∂h2
The derivation of equation (8) is described in detail in Appendix-I. At t=0 and r=R, T=T. It is convenient to introduce a function U (h, t) such that: ∞
(9)
𝑈 = ∫h 𝑇 - 𝑇∞𝑑ℎ
And equation (8) is then written in following form: ∂𝑈 ∂𝑡
( )
= 𝛼𝑙𝑅4
∂2𝑈
(10)
∂h2
It is convenient to introduce a new time variable in equation (10): 𝑡
(11)
𝜃 = ∫0𝑅4(𝑡)𝑑𝑡
equation (11) therefore reduces to the general form of heat diffusion equation: ∂𝑈
( )
∂𝜃 = 𝛼𝑙
∂2𝑈
(12)
∂h2
ACS Paragon Plus Environment
Page 9 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Industrial & Engineering Chemistry Research
The appropriate solution of equation (12) is given by equation (13) taking its Laplace transform subject to the modified initial and boundary conditions (see Appendix-I): 𝑈(h,𝑠) =
𝑈0 ―
𝑠h 𝛼𝑙
𝑆𝑒
(13)
Where, U0=Tv-T Taking Laplace inverse of equation (13), we obtain: 𝑈(h,𝜃) = 𝑈0𝑒𝑟𝑓𝑐
[ ] h 4𝛼𝑙𝜃
(14)
Expanding the erfc function and ignoring the higher order terms in erfc function, we get:
(
𝑈(h,𝜃) = 𝑈0 1 ―
2
h
)
(15)
𝜋 4αlθ
The above leading-order approximation of the higher-order terms is based on the argument that h is very small and r is very close to the bubble wall. This assumption is believed to be appropriate in low superheats. Therefore, our new temperature gradient at bubble wall (r=R) is:
( ) 𝑑𝑇
𝑑𝑟 𝑟 = 𝑅
= ―𝑈0
2 𝑅2
(16)
4𝜋𝛼𝑙𝜃
2.3 Semi-analytical method In our semi-analytical approach, we solved analytically heat diffusion equation applied to the liquid surrounding the bubble, equation (4), to compute temperature gradient
(𝑑𝑇𝑑𝑟)𝑟 = 𝑅at bubble
wall. Then we numerically solved equation (1) and (3) simultaneously for R and Tv at bubble wall. As evident, Tv (boundary condition of equation (4)) and R are also required for computing
(𝑑𝑇𝑑𝑟)𝑟 = 𝑅. This temperature gradient is used in the energy balance to calculate the heat flow (RHS of equation (3)) for evaporating liquid. Thus, the approach avoids the use of dense gridding (Robinson and Judd20) in solving the heat diffusion equation numerically. Specifying the vapour density and pressure as functions of the vapour temperature effectively reduces the problem to determining the instantaneous values of R, ( = 𝑅) and 𝑦3( =
𝑑𝑅 𝑑𝑡
𝑑𝑅 𝑑𝑡
and Tv. In order to do so, the variables 𝑦1( = 𝑇𝑣 ),𝑦2
) are introduced. equation (2) and (4) are rearranged and solved
simultaneously in the form of equation (17) (Robinson and Judd
20)
by fourth order Runge –
Kutta method with a time-step size of 10-9. 𝑑𝑦1 𝑑𝑡
=
𝑑𝑇𝑣 𝑑𝑡
= 𝑓1(𝑦1,𝑦2,𝑦3)
(17a)
ACS Paragon Plus Environment
Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
𝑑𝑦2 𝑑𝑡 𝑑𝑦3 𝑑𝑡
= =
𝑑𝑅 𝑑𝑡
= 𝑓2(𝑦3)
𝑑2𝑅 𝑑𝑡2
Page 10 of 39
(17b)
= 𝑓3(𝑦1,𝑦2,𝑦3)
(17c)
Here, 𝑓1 is determined by the rearranging equation (3).𝑓2 is bubble growth rate.𝑓3is determined by the rearranging equation (1). The initial conditions are y(0)=Tv(0)=T, y2(0)= R(0)=Rc and y3(0)=R(0)=0. 3. NUMERICAL RESULTS AND DISCUSSIONS 3.1. Bubble growth computation In this section, we describe the results of variations of bubble radius with time for five sample test cases of superheats shown in Table .We have presented the results in low superheat and high superheat using our new semi-analytical method, Plesset-Zwick theory8 and Mikic et al. theory16, Prosperetti and Plesset8 and Robinson and Judd20. We also compare our method with the experimental results of Dergarabedian0, Abdelmessih23, Dhir and his group2, 26, and Board and Duffey7. Table . Five test cases of superheats Case #
System/ambient pressure
Superheating
∆𝑻 Experimental
sup
references
1
1.40K (low ∆Tsup)
Dergarabedian0
2
3.10K (low ∆Tsup)
3
1.0 atm
6.83K (low ∆Tsup)
4
10.17K (high ∆Tsup)
5
17.89K (high ∆Tsup)
Growth characteristics Surface
Abdelmessih23
tension,
inertia and thermal diffusion
3.. Bubble growth in low superheats Figures 2-4 show the variations of radius of growing vapour bubble with time in low superheats of 1.40 K to 6.83 K.The figures also show the comparison of our semi-analytical predictions, Plesset-Zwick theory, Mikic et al. theory and the experimental results of Dergarabedian0 and Abdelmessih23.
ACS Paragon Plus Environment
Page 11 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Industrial & Engineering Chemistry Research
In Figure 2, we see that at a low superheat of 1.40 K, our semi-analytical prediction (new solution method) and the experimental result are very close to each other. We also see that the prediction using Plesset-Zwick theory is higher than the experimental results. The reason of overprediction by Plesset-Zwick theory is that they approximated the asymptotic solution of bubble radius given by equation (18) to a leading order term and ignored the higher order of time. This approximation is valid for the latter stage of bubble growth and results in some uncertainties at the initial stage of bubble growth. However, when the bubble is growing over a long time, the growth is heat diffusion-controlled. 𝑅=
12𝑘(𝑇∞ ― 𝑇𝑠𝑎𝑡) 1/2 𝜋 ℎ𝑓𝑔𝑣 𝛼 𝑡
(18)
The bubble radius predicted by Mikic et al. theory (equation (19)) is significantly higher than the experimental results in low superheats. 3
𝑑𝑅 + 𝑑𝑡
3
= (𝑡 + + 1)2 ― (𝑡 + )2 (19a) Their formulation results in two following limiting solutions: R=At
for t+≤1 (inertia-controlled growth)
(19b)
R=Bt/2
for t+>1 (heat diffusion-controlled growth)
(19c)
equation (19) is known as MRG (Mikic, Rohsenow and Griffith) solution. We see in equation (20) that the constant A becomes very small at an extremely low superheat, and the first term dominates over the second term. Consequently, the growth becomes inertia-driven and very fast compared to the heat diffusion-controlled growth. Therefore, Mikic et al. theory shows higher bubble radius at low superheat of 1.4 K. 1 𝑑𝑅 2 𝐴2 𝑑𝑡
( )
+
𝑡𝑑𝑅 𝐵 𝑑𝑡
(20)
-1=0
Mikic et al. assumed that the bubble started to grow from R = 0 at t = 0 rather than R=Rc at t = 0. We also note that the theory predicts higher initial growth rate in all cases of superheat. It is due to the fact that the theory ignores the effect of surface tension.
ACS Paragon Plus Environment
Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
In Figure 3, we see that our semi-analytical results are very close to the experimental results. We also note that Plesset-Zwick theory and Mikic et al. theory (MRG solutions) predict appreciably high bubble radius compared to the experimental results, particularly in longer time. This is due to the fact that they ignored the effect of higher order term of time. Figure 4 shows an excellent agreement between our semi-analytical predications and the experimental results, and a slightly higher prediction by Plesset-Zwick theory and Mikic et al. theory at a low superheat 6.83 K.
Figure 2. Variation of bubble radius with time for Case 1 at Tsup=1.40 K
ACS Paragon Plus Environment
Page 12 of 39
Page 13 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Industrial & Engineering Chemistry Research
Figure 3. Variation of bubble radius with time for Case 2 at Tsup=3.1 K
Figure 4. Variation of bubble radius with time for Case 3 at Tsup=6.83 K
ACS Paragon Plus Environment
Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
3..2 Bubble growth in high superheats We discuss here two test cases of bubble growth in high superheats (Tsup=10.17 K and 17.89 K). The computational results were compared with Plesset-Zwick theory, Mikic et al. and the experimental results of Abdelmessih. In Figure 5, the bubble radius predicted by our method and the analytical models of Plesset-Zwick and Mikic et al. are significantly higher than the experimental results at a high superheat of 10.17 K. We also see the same trend at high superheat of 17.89 K in Figure 6. The reason of higher predictions by Mikic et al. follows. We see in equation (20) that A will be very large at a high superheat and the second term will dominate over the first term, leading to heat-diffusion controlled growth (R t). It clearly indicates that Plesset-Zwick diffusion theory is valid in high superheats. Overall, a very good agreement of our predictions for bubble radius with the experimental results is found in low superheats.
Figure 5. Variation of bubble radius with time for Case 4 at Tsup=10.17 K
ACS Paragon Plus Environment
Page 14 of 39
Page 15 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Industrial & Engineering Chemistry Research
Figure 6. Variation of bubble radius with time for Case 5 at Tsup=17.89 K
3.2 Comparison of our semi-analytical method and other methods and experiments In this section we compare our semi-analytical method with the analytical method of Prosperetti and Plesset18, the computational fluid dynamics (CFD) method (numerical simulation) of Dhir2 and the pool boiling experiment of Dhir2 at Tsup=7 K in Figure 7(a). Prosperetti and Plesset derived equation (2) for the scaled bubble radius R* over scaled time t*, which is very similar to MRG solution. ∗
𝑅 =
2
( )( ) [( 2 3
2
𝜋2
∗
𝜋2 ∗ 2𝑡
2
3 2
3 2
) ( ) ― ]
+ ―
𝜋2 ∗ 2𝑡
∗
2
Where, 𝑅 = 𝜇 𝑅/𝑅𝑐; 𝑡 = 𝛼𝜇 𝑡; 𝜇 =
(2)
2𝜎𝛼𝑙 2 𝜌𝑔ℎ𝑓𝑔 3 𝜋 𝑘𝑙(𝑇∞ ― 𝑇𝑣)
( )
[𝜌𝑓{𝑝𝑣(𝑇∞) ― 𝑝∞}] ― 4
3.2. Comparison with pool boiling experiments In the pool boiling experiment, Dhir2 used wafer of diameter of 0 cm and thickness of mm. A cylindrical cavity of 0 m was created on the wafer surface. The surface was heated by thin-
ACS Paragon Plus Environment
Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 16 of 39
film strain gage heaters (6.5 mm 6.5 mm area was covered by each heater) at bottom. Temperature was measured by thermocouples. A CCD camera captured the boiling phenomenon at a frame rate of 020 fps. A comparison of our semi-analytical method (without modification) with the experiments have been presented in Figure 7(a), which shows the variation of bubble diameter with time at Tsup = 7.0 K and p = .0 atm. The analytical methods of Plesset-Zwick, Mikic et al. and Prosperetti and Plesset are also presented in the figure. All the methods significantly overpredict the bubble diameter compared to the experimental results due to large value of heat flux (very thin boundary layer) from liquid to vapour. To avoid the overprediction, we modified the heat flux term in our method by incorporating the heat transfer coefficient in terms of Nu calculated from the following correlations of Dhir2 and Berenson28:
𝑁𝑢 = 0.264(𝐺𝑟𝑃𝑟) × 4
(
𝐽𝑎𝑣
Nu = 0.264(GrPr/Jav∗ )4
)
+ 0.5 + .3𝐽𝑎𝑣
4
(Dhir2) (Berenson28)
(22) (23a)
𝐽𝑎𝑣
𝐽𝑎𝑣∗ = + 0.5𝐽𝑎
𝑣
(23b) Nu, Gr and Pr are the usual symbols meant in the section of Nomenclature. Jav is Jakob number based on the vapour phase. The modification led to a remarkable improvement of the prediction by our semi-analytical method. The prediction of bubble diameters using the Nu obtained from the Dhir’s correlation (equation (22)) very closely agrees with the experimental bubble diameters compared to that using the Nu obtained from the Berenson’s correlation (equation (23)). Evaluated by Dhir2, the Berenson’s correlation calculates relatively high Nu. In addition, our prediction is also relatively close to the CFD results of Dhir2. We carried out the computation upto 40.0 s to avoid the bubble departure, because we have not included the force balance for the bubble departure in the present model. It is required to mention that there is no provision for the inclusion of the desired heat transfer model in the analytical model of Plesset-Zwick, Mikic et al. and Prosperetti and Plesset. At a relatively low superheat of 5.0 K, the prediction of bubble diameter by our method due to the modification, shown in Figure 7(b), is consistent with and close to the experimental results of Mukherjee and Dhir29, and the CFD results of Dhir30. The under-prediction by our modified model is due to low heat transfer from the superheated liquid calculated by equation (22) and
ACS Paragon Plus Environment
Page 17 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Industrial & Engineering Chemistry Research
more importantly, the physics of bubble merger is not included in our method. Overall, our semimethod closely predicts the experimental bubble diameter in pool boiling.
Figure 7(a). Comparison of our semi-analytical method and other methods for bubble diameter at Tsup=7.0 K and p=1.0 atm 3.2.2 Comparison with other experiments We also compare our method with the experiments of Board and Duffey17 at Tsup=4.3 K and p = 0.395 atm, and the other analytical methods of Plesset-Zwick, Mikic et al., and Prosperetti and Plesset and the numerical simulation of Robinson and Judd. The predicted bubble radii by all the methods are higher than the experimental results as we noted the large gap between prediction and experiment for bubble growth in high superheats presented in Section 3. and 3.2. It is interesting to note that despite the low superheat, the large gap also appears for the bubble growth at very low pressure (p = 0.395 atm). It is because of the effect of inertia in the computations dominates over the heat transfer control. We discussed this issue in detail in Section 3.3. Nevertheless, the predictions using the methods are very close to each other. The bubble radii after 9.07 ms as evident in Figure 7(c) are 2.32 mm (our method), 2.37 mm (Plesset-
ACS Paragon Plus Environment
Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Zwick), 2.35 mm (Mikic et al.), 2.34 mm (Prosperetti and Plesset), 2.377 mm (Robinson and Judd). The discrepancy is within 2%. The experimental radius is .763 mm. As the time progresses, all the above methods yield further close to bubble radii. Although this discrepancy is varying in nature depending on superheat and operating pressure, our proposed computationallyefficient semi-analytical method calculates reasonably accurate bubble radius.
Figure 7(b). Comparison of our semi-analytical method and other methods for bubble diameter at Tsup=5.0 K and p=.0 atm
ACS Paragon Plus Environment
Page 18 of 39
Page 19 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Industrial & Engineering Chemistry Research
Figure 7 (c). Comparison of our semi-analytical method and other methods for bubble radius at Tsup=4.3 K and p = 0.395 atm In Figure 7b and 7c, the experimental observation could not be adequately predicted by our semianalytical model. It is due to some uncertainties of heat transfer calculation from the superheated liquid by eq. (22) under the operating conditions. The bubble is assumed to grow to a perfect sphere. The contribution of microlayer evaporation is not considered in our model. Micro-flows including the Marangoni convection would improve the computation of bubble growth, although the temperature dependency of the surface tension is calculated from our water-property codes. In the saturated pool boiling, the effect of surface tension gradient due to the temperature gradient is not much as in the subcooled boiling, except at interface. However, the degree of superheat in the present study is not so high to have an inverse effect by the convection on the interfacial heat transfer. The resistance due to the Marangoni force on the bubble for departing the surface is negligible. The full-scale CFD-model along with micro-flow computations would be appropriate to bridge this gap. 3. 2.3 Evaluation of our linear approximation Now, we examine the performance of our linear approximation of the temperature profile of surrounding liquid. We computed the deviations of the theoretical bubble radius from
ACS Paragon Plus Environment
Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 20 of 39
experimental bubble radius in low superheats and in high superheats. Table 2 presents the percentage deviations from the experimental results. In low superheats, our results are very close to the experimental results while the results of Plesset-Zwick and Mikic et al. show a large deviation. The reason for relatively higher deviations of our results from the experimental result and the Plesset-Zwick theory at high superheat and is that we have approximated the asymptotic (series) solution for temperature profile of liquid surrounding the bubble to be a linear one (see equation (15)). From Figure 1 and equation (6), we see that the approximation is valid when h is very small. It occurs when r is very close to bubble wall. And the surrounding temperature 𝑇∞ is also very close to Tv at bubble wall, allowing the approximation to be applicable in low superheats. In high superheat, the approximation permits slightly higher heat diffusion from the surrounding liquid. In high superheats, all the theoretical results are consistent and are close to each other. We also note that at higher superheat, the theories predict higher radius compared to experiments. The theories follow inertia- and heat diffusion-controlled growth at high superheat. An illustration has been made. In addition, the theories and our semi-analytical method compute high inertia force to give rise to an increased bubble radius. Prosperetti and Plesset18 concluded that the rate of bubble growth for large superheats was somewhat overestimated in the intermediate stage in which both inertial and thermal effects concurrently played a vital role. However, the overestimate does not practically yield a serious error in the radius-time behavior. Table 2 Accuracy of the predicted results Case
Experiments
Superheats
Theories
∆𝑻𝒔𝒖𝒑
%
Relative
error
of
bubble radius 1
1.40K Dergarabedian0
2
3.1K
Our semi-analytical method
2.21
Plesset-Zwick
75.77
Mikic et al. equation (19b-c)
71.196
Our semi-analytical method
7.17
Plesset-Zwick
12.59
Mikic et al. equation (19b-c)
13.60
ACS Paragon Plus Environment
Page 21 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Industrial & Engineering Chemistry Research
3 Abdelmessih23
6.83K
4
10.17K
5
17.89K
Our semi-analytical method
1.23
Plesset-Zwick
9.86
Mikic et.al equation (19b-c)
10.21
Our semi-analytical method
51.15
Plesset-Zwick
42.98
Mikic et al. equation (19b-c)
59.97
Our semi-analytical method
57.77
Plesset-Zwick
62.17
Mikic et al. equation (19b-c)
31.273
In order to test the computational accuracy of our semi-analytical method, we performed the computation of bubble radius using three different time steps of 10-7 s, 10-8 s and 10-9 s at Tsup = 1.40 K. All the results were found to be almost equal; the change in sampling time from 10-7 s to 10-9 s did not affect the vapour bubble growth. We used a time step of 10-9 s in the present study. We also demonstrate here the convergence and the stability of our method by deriving an asymptotic solution. From equation (6), we have successfully derived an asymptotic solution for
( ) 𝑑𝑇
𝑑𝑅
𝑑𝑟 𝑟 = 𝑅 ( 𝑑𝑡 0
(𝑑𝑇𝑑𝑟)𝑟 = 𝑅 = ―2𝑈
for 𝑡) in heat-diffusion controlled growth:
0 4𝜋𝛼𝑙
𝑡
―2
(25)
The above asymptotic solution is comparable with Plesset-Zwick theory8 of heat-diffusion controlled growth having the growth rate: 𝑑𝑅 𝑑𝑡
=𝐶
𝛼𝑙
𝜋 𝐽𝑎𝑡
―2
(26)
with Ja = Jakob number (=cplρl∆Tsup/ρvhfg), C=3 for Plesset-Zwick theory and C=1 for our asymptotic solution. The asymptotic solution shows a good validation of our approximation. 3.3 Effect of operating pressure on bubble growth It has been noted in the experiments of Akiyama et al.3 that the operating pressure has a significant influence on the growth of vapour bubble. The bubble diameter decreases remarkably
ACS Paragon Plus Environment
Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
as we increase the operating pressure. It is due to very fast decrease of bubble growth rate by increasing the operating pressure. Figure 8 shows the variation of bubble diameter with time at different pressures of .0 atm, 2.0 atm and 5.0 atm. The required heat flux is 40.706 kW/m2. The predictions by our method, and the analytical methods of Plesset-Zwick, Mikic et al. and Prosperetti and Plesset agree well with the experiments of Akiyama et al. (variation of bubble diameter with time). The agreement at higher pressure is better. It has been noted that the inertia force is greatly decreased by the increase of pressure. Consequently the bubble growth rate significantly decreases. The growth regime quickly switches to heat-diffusion controlled at high pressure. To understand how the bubble-growth kinetics (R=A0tn) is affected by the operating pressures (0.395-3.0 atm), a numerical simulation has been performed at Tsup= 4.3 K using our method and the analytical methods of Plesset-Zwick, Mikic et al. and Plesset and Prosperetti. The numerical results of bubble radius varying with time have been presented in Figure 9(a). In low pressures, the prediction by our method differs from the other methods. As the pressure is increased, the growth exponent n computed by semi-analytical method decreases from 0.72 (at 0.395 atm) to a limiting value of about 0.58 while n computed by the other analytical methods remains at about 0.5 regardless of pressure variation. In addition, the growth constant (A0) is appreciably low in our method. It is because of the fact that the boundary layer thickness around the bubble in our method is much affected by low pressures. However, the decrease of bubble radius with pressure follows a power law with an exponent of 0.90, as shown Figure 9(b).
ACS Paragon Plus Environment
Page 22 of 39
Page 23 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Industrial & Engineering Chemistry Research
Figure 8. Effect of operating pressure on bubble growth and comparison with experiments To compare the predicted growth kinetics at Tsup= 4.3 K with the experimental growth kinetics of Akiyama et al. at heat flux of 40.706 kW/m2 in the range of pressures of 0.395-3.0 atm, A0 in of the predicted growth kinetics and of the experimental growth kinetics have been plotted against operating pressure in Figure 9(c). It is evident in Figure 9(c) that the predicted A0 and the experimental A0 are consistent. The decrease of A0 for prediction and experiment follows a power-law variation with pressure with exponent of 0.85 and 0.86 respectively. We also noted that the experimental n varies from about 0.8 to about 0.46 as the pressure is increased. The variation of n with pressure in the experiments closely agrees with that in our semi-analytical method. n obtained by our method varies from 0.72 to about 0.58 with increase in pressure (see Table 3). The high value of n in low pressures indicates that unlike the analytical methods, the semi-analytical method yields inertia-controlled bubble growth in low pressures. This is due to our leading-order approximation in temperature profile and the inclusion of the dependency of phase densities and surface tension on temperature and of variable pressure difference (pv-p) across the bubble wall under varying vapour temperature (Tv).
ACS Paragon Plus Environment
Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 24 of 39
Figure 9(a). Effect of operating pressure variation on the bubble growth at Tsup = 4.3 K 3.5 Interaction of energy components The interaction of various energy components during bubble growth plays an important role in generating feedbacks. Time-variation of these components reveals some interesting dynamical features of bubble growth. For convenience, we define the following energy components (J/kg): (1) Pressure difference component across the interface (PDC) = 3 𝑑𝑅 2 (2) Inertia component (IC) = 2 𝑑𝑡
( )
2𝜎
(3) Surface tension component (STC) = 𝜌 𝑅 𝑙 (4) Thermal diffusion component (TDC) =
𝑑𝑇
―𝑘𝑙𝑑𝑟
|𝑟 = 𝑅4𝜋𝑅2
4𝜋𝑅3𝜌𝑙/3
ACS Paragon Plus Environment
𝑃𝑣(𝑇𝑣) ― 𝑃∞ 𝜌𝑙
Page 25 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Industrial & Engineering Chemistry Research
Figure 9(b). Power law variation of bubble radius with pressure at Tsup = 4.3 K
Figure 9(c). Effect of pressure on the growth constant 3.5. Energy components in low superheats
ACS Paragon Plus Environment
Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 26 of 39
Figure 0(a-b) shows the time-variations of PDC, IC, STC and TDC during bubble growth at a low superheat of 1.4 K obtained using our semi-analytical method and Plesset-Zwick theory. PDC and STC vary almost identically and gradually decrease with time in the semi-analytical method (Figure 0(a)). They do not vary identically in Plesset-Zwick theory (Figure 0(b)); PDC remains unchanged with time while STC decreases with time. The constant PDC is one of the approximations in the Rayleigh theory and Plesset-Zwick theory, when R is very large compared to Rc. The rate of growth or collapse following the theories is essentially heat transfer controlled. STC decreases with time because of increasing bubble radius. Table 3 Effect of operating pressure on growth exponent at ∆Tsup=4.3 K Operating
n
n
n
n
pressure
(Semi-analytical
(Plesset-Zwick)
(Mikic et al.)
(Prosperetti
(atm)
method)
0.395 .00 .50 2.00 2.50 3.00
and
Plesset) 0.720
0.486
0.532
0.584
0.613
0.487
0.502
0.525
0.597
0.488
0.500
0.515
0.590
0.488
0.500
0.510
0.587
0.488
0.500
0.507
0.584
0.488
0.500
0.506
The reason of decreasing PDC with time in the semi-analytical method can be explained based on the energy equation at vapour-liquid interface (see equation (3)), which is rewritten as: 𝑑𝜌𝑣 𝑑𝑡
=
( ― 𝑘𝑙∂𝑇∂𝑟)𝑟 = 𝑅𝑑𝐴 ―
3 0 ∫ 4𝜋𝑅2ℎ𝑓𝑔𝑅 𝐴𝑠
3𝜌 𝑣ℎ𝑓𝑔𝑑𝑅 ℎ𝑓𝑔 𝑅 𝑑𝑡
(25)
In the above equation, the first term on RHS is the heat flow due to thermal diffusion to vapourliquid interface and the second term is the latent heat required for the phase change. The LHS is due to the energy required to maintain the saturation temperature. The latent heat is always significantly higher than the heat flow and the deficient amount for the latent heat is compensated by the heat released from due to the drop in vapour temperature. Thus, the vapour
ACS Paragon Plus Environment
Page 27 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Industrial & Engineering Chemistry Research
temperature decreases with time (with increasing vapour bubble radius) and so are the vapour pressure and the density of saturated vapour. The variation of IC and TDC are consistent and advances through two maxima. After the second maximum (t=7.76 ms for IC and t=5.606 ms for TDC), they start to asymptotically decrease (noted in a growth over 40 ms) in the semi-analytically method (Figure 0(a)). It indicates that both inertia and thermal diffusion play significant role in the latter stage of the bubble growth in low superheats (the effect of surface tension starts to decrease approximately t=5.606 ms), as evident in Figure 0(a). The inertia then decreases due to the decrease of the thermal diffusion. It is also noted the semi-analytical method calculates large STC compared to IC. In contrary, STC, IC and TDC decrease asymptotically in Plesset-Zwick theory from t=0 (Figure 0(b)). The leading-order time approximation made by Plesset-Zwick causes a very fast asymptotic decrease of thermal energy and energy of inertia from t=0. The effect of thermal diffusion and inertia is significantly stronger in Plesset-Zwick theory (Figure 0(b)) from t=0 compared to those in our semi-analytical method (Figure 0(a)). It clearly explains why the approximation of PlessetZwick calculates higher bubble radius in very low superheat. STC and IC are almost equal at t=0 in Plesset-Zwick theory. As the time progresses, the effect of STC in both the methods continues to be significant for a sufficient time and is greater than the inertia effect at a low superheat.
Figure 0(a). Variations of energy components at Tsup=1.40 K (semi-analytical method)
ACS Paragon Plus Environment
Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
3.5.2 Energy components in high superheat Increasing the superheat to a high value causes the effect of STC, IC and TDC to decrease asymptotically. Figure (a) shows the time-variations of PDC, STC, IC and TDC during bubble growth at a high superheat of 17.89 K in our semi-analytical method. The difference between STC and IC decreases with increase in superheats. Therefore, the combined effect of inertia and diffusion quickly dominates over the effect of surface tension. Unlike in low superheats, the variations of PDC and STC with time differ significantly. IC is less than PDC and STC for a short duration from t=0. Interestingly, IC exceeds STC and then exceeds PDC after some time and continues to be higher than STC and PDC for a long time. Then, IC becomes less than PDC. Our semi-analytical method indicates the growth is inertia-controlled and thermal diffusioncontrolled after a short while from t=0, in which the effect of surface tension dominates.
Figure 0(b). Variations of energy components at Tsup=1.40 K (Plesset-Zwick theory) In the Plesset-Zwick theory for high superheats, STC, IC and TDC almost decrease asymptotically from t=0 (the initiation of bubble growth). PDC remains constant in PlessetZwick theory. IC is always greater than STC. IC becomes greater than PDC for a short while from t=0 and then continues to be less than PDC. Plesset-Zwick theory indicates that the growth
ACS Paragon Plus Environment
Page 28 of 39
Page 29 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Industrial & Engineering Chemistry Research
is inertia-controlled and thermal diffusion-controlled from t=0, and the effect of surface tension is negligible in high superheats. To confirm the above illustration in high superheat, the growth rate computed by the semianalytical method has been compared with that by the analytical methods of Plesset-Zwick, Mikic et al., Prosperetti and Plesset, and the numerical simulation of Robinson and Judd at a high superheat of Tsup= 15 K and p = 1.0 atm in Figure 12. Robinson and Judd solved the heat diffusion problem for the liquid around the bubble using implicit finite difference method.
Figure (a).Variations of energy components at Tsup=17.89 K (semi-analytical method)
ACS Paragon Plus Environment
Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Figure (b). Variations of energy components at Tsup=17.89 K (Plesset-Zwick theory) The thin right-arrows mean that the effect of inertia becomes stronger in the semi-analytical method and the numerical simulation of Robinson and Judd. The thin left-arrows indicate that the effect of inertia become stronger upto t=0 in the analytical methods of Mikic et al. and Prosperetti and Plesset. The bold left-arrows mean that the effect of surface tension becomes important upto t=0 in our semi-analytical method and the numerical simulation of Robinson and Judd. Figure 12 shows that all the methods compute almost equal growth rate in the heatdiffusion controlled growth (the decrease of growth rate follows a power law with an exponent close to about 0.50), although the prediction by Prosperetti and Plesset is slightly low. The constant growth rate computed by Mikic et al. and Prosperetti and Plesset indicates the inertiadominant growth, which continues for a long time compared to the short duration of inertiadominant growth predicted by our semi-analytical method and Robinson and Judd. The simplifications made by Mikic et al. and Plesset and Prosperetti cannot capture the surface tension dominated growth-zone. Our semi-analytical method and the numerical simulation of Robinson and Judd clearly and consistently identify the surface tension dominated growth-zone (increasing section of the growth rate); the semi-analytical method does not compute the effect of
ACS Paragon Plus Environment
Page 30 of 39
Page 31 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Industrial & Engineering Chemistry Research
surface tension as much as that calculated by Robinson and Judd. The zone persists for a short duration. This numerical uncertainty due to low-dimensional model does not practically affect the computation of eventual bubble radius; all the methods compute almost same bubble radius (t 0.084 ms) despite the different growth mechanisms proposed by the methods.
Figure 12. Comparison of our semi-analytical method with other methods from literature for bubble growth rate (inertia) at Tsup=15 K and p = 1.0 atm Figure 3 presents the temporal variations of various energy components depending Tsup, as computed by the semi-analytical method. We note that in low superheats, the inertia effect is very negligible compared to the effect of surface tension and pressure difference across the bubble wall. After a long time, the growth becomes heat-diffusion controlled. As we increase the superheat gradually, the inertia effect overcomes the effect of surface tension after a time t=t1 and then exceeds the effect of pressure difference after a time t=t2 > t (see Figure 3(c-d)). 3 𝑑𝑅 2 2 𝑑𝑡
( )
=
3 𝑑𝑅 2 2 𝑑𝑡
=
( )
2𝜎 𝜌𝑙𝑅 𝑝𝑣 ― 𝑝∞ 𝜌𝑙
(at t=t)
(26)
(at t=t2)
(27)
After t=t2, the growth is controlled by both inertia and diffusion, leading to a large bubble radius.
ACS Paragon Plus Environment
Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Figure 3. Temporal variation of energy components depending Tsup at p=.0 atm (computed by semi-analytical method) We also see in Figure 3 that the peak values of IC, TDC and the peak time for the growth rate computed by our semi-analytical method increase with increase in superheat. As we increase the superheat, the bubble growth rate increases and so are the effect of inertia and thermal diffusion. It is evident from the figure that the peak time corresponding to the time, at which the growth rate is maximum, and it varies from 5.606 ms to 0.003 ms as we increase the superheat from 1.4 K to 17.89 K. This peak time indicates the time, at which the effect of inertia and thermal diffusion become important. Interestingly, increasing the superheat makes the effect of inertia and thermal diffusion to occur much earlier. At Tsup ≥ 10.17K, the effect of inertia and thermal diffusion becomes important and the effect of surface tension becomes negligible.
4. CONCLUSIONS We have established that our proposed computationally-efficient semi-analytical method adequately describes the bubble growth in varying superheat at atmospheric pressure, low pressure and high pressure. The numerical results have been validated with the experiments, the other analytical methods and higher-order numerical simulations from literature. Unlike the analytical methods, the semi-analytical method makes a leading-order approximation in temperature profile and includes the dependency of the phase densities and surface tension on
ACS Paragon Plus Environment
Page 32 of 39
Page 33 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Industrial & Engineering Chemistry Research
temperature and of variable pressure difference across the bubble wall under varying computed vapour temperature. This solution could be useful for CFD-result verifications, when we do not have experimental results. The present semi-analytical method can be successfully extended by modifying heat flux computed from Nu using Dhir’s correlation to predict the experimental bubble diameter in pool surface-boiling. The method has been found to be reasonably accurate before the bubble departs. The modification led to a remarkable improvement of the prediction by our semi-analytical method. Analysis of the interactions of various energy components such as pressure difference, inertia (hydrodynamics), surface tension and thermal diffusion has also been carried out. We summarize the major findings in the following: 1. The semi-analytical method is capable of computing the bubble growth with a high accuracy in very low superheats (the effect of surface tension dominates) compared to the analytical methods. In high superheats, all the theoretical results are consistent and are very close to each other. 2. The effect of pressure has a significant impact on the bubble-growth kinetics. The decrease of growth constant A0 follows a power-law variation with pressure with an exponent of 0.85. Unlike the existing analytical methods, the semi-analytical method yields inertia-controlled bubble growth in low pressures; the analytical methods essentially yield heat diffusion-controlled growth regardless of pressure variation. In high pressures, all the methods predict thermal-diffusion controlled growth. 3. The rigorous analysis of the interactions of various energy components in low and high superheats reveals that the effect of surface tension dominates for a sufficient time, and both inertia and thermal diffusion play a significant role in the latter stage of the bubble growth in low superheats, as predicted by the semi-analytical method. In future, we shall incorporate the model for bubble departure based on force balance in our semi-analytical models and compare the results in low pressures with our experiments and the experiments from literature. We also believe that the higher-order models along with the microflow computations are equally essential for improved computation of more complex bubble dynamics and the related interfacial instability.
ACS Paragon Plus Environment
Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
NOMENCLATURE A
Constant in equation (19)
As
Bubble surface area (m2)
A0
Growth constant
B
Constant in equation (19)
𝐶𝑝
Specific heat (J/kg-K)
Gr
Grashof number
h
Lagrangian coordinate, equation (6)
hfg
Latent heat of evaporation (J/Kg)
Ja
Jakob number, cplρl∆Tsup/ρvhfg
Jav
Jakob number, cpv(T-Tv)/hfg
k
Thermal conductivity, W/m-K
n
Growth exponent
Nu
Nusselt number
p
Pressure (Pa)
𝑝𝑣
Vapour pressure (Pa)
𝑃∞
Bulk liquid pressure (Pa)
Pr
Prandtl number
r
Radial coordinate
R
Bubble radius (m)
𝑅+
Dimensionless bubble radius, 𝑅 + = 𝐵2 ;
𝐴𝑅
𝐴= 𝑏
(
𝜌𝑙𝑇𝑠𝑎𝑡
1 2
1 2
) ; 𝐵 = ( ) 𝐽𝑎; b=(2/3)
ℎ𝑓𝑔𝜌𝑣∆𝑇
12𝑎𝑙 𝜋
0.5
R*
Dimensionless bubble radius, equation (2)
𝑇𝑠𝑎𝑡
Saturation temperature (K)
𝑇𝑣
Vapour temperature (K)
𝑇∞
Bulk liquid Temperature (K)
t
Time (s)
𝑡+
Dimensionless bubble growth time,𝑡 + =
𝐴2𝑡 𝐵2
ACS Paragon Plus Environment
Page 34 of 39
Page 35 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Industrial & Engineering Chemistry Research
t*
Dimensionless bubble growth time, equation (2)
u
Local velocity (m/s)
U
∫ℎ 𝑇 ― 𝑇∞𝑑ℎ , a new variable for temperature, equation (9)
U0
Tv-T (K)
∞
Greek Symbols α
Thermal diffusivity (m2/s)
𝛿
Thickness of Thermal boundary layer (m)
∆
Differences
Far field from bubble
μ
Dimensionless parameter, equation (2)
ρ
Density(kg/m3 )
𝜎
Surface tension (N/m)
θ
New variable for time, equation (1)
Subscripts c
Critical
l
Liquid phase
sat
Saturation
sup
Superheated
v
Vapour phase
Bulk liquid
Abbreviations CFD IC PDC STC
Computational fluid dynamics 3 𝑑𝑅 2
( )
Inertia component = 2 𝑑𝑡
Pressure difference component across the interface = 2𝜎 Surface tension component = 𝜌 𝑅 𝑙 𝑑𝑇
TDC Thermal diffusion component =
―𝑘𝑙𝑑𝑟
|𝑟 = 𝑅4𝜋𝑅2
4𝜋𝑅3𝜌𝑙/3
ACS Paragon Plus Environment
𝑃𝑣(𝑇𝑣) ― 𝑃∞ 𝜌𝑙
Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
ACKNOWLEDGMENT The authors gratefully acknowledge the financial support of Department of Science and Technology (DST), New Delhi, under SERB scheme for this work.
SUPPORTING INFORMATION Appendix-I contains the detailed derivation for temperature gradient at bubble wall. REFERENCES 1. Westwater, J. W. Advances in Chemical Engineering. Adv. Chem. Eng. 1956, 1, 1. 2. Dhir, V. K. Numerical simulations of pool-boiling heat transfer. AIChE Journal 2001, 47(4), 813. 3. Tomar, G.; Biswas, G.; Sharma, A.; Agrawal, A. Numerical simulation of bubble growth in film boiling using a coupled level-set and volume-of-fluid method. Phys. Fluids. 2005, 17, 112103-1. 4. Uguz, E. K.; Narayanan, R. Instability in evaporative binary mixtures, part I- the effect of solutal marangoni convection. Phys. Fluids 202, 24, 094101-1. 5. Ganguli, A. A.; Pandit, A. B.; Joshi, J. B. Bubble dynamics of a single condensing vapor bubble from vertically heated wall in subcooled pool boiling system: Experimental measurements and CFD simulations. Int. J. Chem. Eng. 202, ID 712986-1. 6. Rayleigh, L. On the pressure developed in a liquid during the collapse of a spherical cavity. Philos. Mag. Series 6 1917, 34, 94. 7. Plesset, M. S. The Dynamics of cavitation bubbles. J. Appl. Mech. 1949, 16, 277. 8. Plesset, M. S.; Zwick, S. A. The growth of vapor bubbles in superheated liquids. J. Appl. Phys. 1954, 25, 493. 9. Plesset, M. S.; Zwick, S. A. A nonsteady heat diffusion problem with spherical symmetry. J. Appl. Phys. 1952, 25, 493. 10. Dergarabedian, P. The rate of growth of vapor bubbles in superheated water. J. Appl. Mech. 1953, 20, 537.
ACS Paragon Plus Environment
Page 36 of 39
Page 37 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Industrial & Engineering Chemistry Research
11. Forster, H. K.; Zuber, N. Growth of a vapour bubble in a superheated liquid. J. Appl. Phys. 1954, 25, 474. 12. Scriven, L. E. On the dynamics of phase growth. Chem. Eng. Sci. 1959, 10, 1. 13. Hsieh, D-Y. Some analytical aspects of bubble dynamics. J. Basic Eng. 1965, 87, 99. 14. Theofanous, T. G.; Baisi, L.; Isben, H. S.; Fauske, H. A theoretical study on bubble growth in constant and time-dependent pressure fields. Chem. Eng. Sci. 1969, 24, 885. 15. Abdelmessih, A. H. Spherical bubble growth in a highly superheated liquid pool. Cocurrent gas-liquid flow, Plenum Press. 1969, 485–495. 16. Mikic, B.; Rohsenow, W.; Griffith, P. On bubble growth rates. Int. J. Heat Mass Trans. 1970, 13, 657. 17. Board, S. J.; Duffey, R. B. Spherical vapour bubble growth in superheated liquids. Chem. Eng. Sci. 1971, 26, 263. 18. Prosperetti, A.; Plesset, M. S. Vapour-bubble growth in a superheated liquid. J. Fluid Mech.1978, 85, 349. 19. Lee, H. S.; Merte, Jr. H. Spherical vapor bubble growth in uniformly superheated liquids. Int. J. Heat Mass Trans. 1996, 39, 2427. 20. Robinson, A. J.; Judd, R. L. The dynamics of spherical bubble growth. Int. J. Heat Mass Trans. 2004, 47, 5101. 21. Mohammadein, S. A.; Gouda, S. A. Temperature distribution in a mixture surrounding a growing vapour bubble. Heat Mass Trans. 2006, 42(5), 359. 22. Kudryashov, N. A.; Sinelshchikov, D. I. Analytical solutions for problems of bubble dynamics. Phys. Lett. A 205, 379, 798. 23. Sharipov, V. High Order Bubble Dynamics in incompressible Liquid. J. Heat Trans. 2015, 37, . 24. Bhati, J.; Naik L, J.; Paruya, S. A semi-analytical method for computing bubble growth in superheated water. Presented in NURETH-17 (17th International Topical Meeting on Nuclear Reactor Thermal Hydraulics), September 3-8, 2017, Xian, China. 25. Pandey, V.; Biswas, G.; Dalal, A. Effect of superheat and electric field on saturated film boiling. Phys. Fluids. 2016, 28, 052102-1.
ACS Paragon Plus Environment
Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
26. Hens, A; Biswas, G; De, S. Analysis of interfacial instability and multimode bubble formation in saturated pool boiling using Coupled Level Set and Volume- of- Fluid approach. Phys. Fluids. 2014, 26, 012105-1. 27. Karapetsas, G.; Sahu, K. C.; Matar, O. K. Evaporation of sessile droplets laden with particles and insoluble surfactants. Langmuir 2016, 32, 6871. 28. Berenson, P. J. Film Boiling Heat Transfer from a Horizontal Surface. J. Heat Trans. 1961, 83, 351. 29. Mukherjee, A.; Dhir, V. K. Study of lateral merger of vapor bubbles during nucleate Pool Boiling. J. Heat Trans. 2004, 126, 1023. 30. Dhir, V. K. Mechanistic prediction of nucleate boiling heat transfer–achievable or a hopeless task? J. Heat Trans. 2006, 28, . 31. Akiyama, A.; Tachibana, F.; Ogawa, N. Effect of pressure on bubble growth in pool boiling. Bulletin of JSME 969, 2, 2.
ACS Paragon Plus Environment
Page 38 of 39
Page 39 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Industrial & Engineering Chemistry Research
ACS Paragon Plus Environment