Boiling Crisis as a Critical Phenomenon

We present the first experimental study of intermittency and avalanche distribution during a boiling crisis. To understand the emergence of power law statistics we propose a simple spin model capturing the measured critical exponent. The model suggests that behind the critical heat flux is a percolation phenomenon involving drying-rewetting competition close to the hot surface. Among known hydrodynamic instabilities, a boiling crisis is one of the most mysterious and poorly understood [1,2]. It may be described as a sudden loss of liquid contact with a heating surface (dryout) when a small increase of heat flux transforms nucleate boiling into film boiling. Reaching the critical heat flux (CHF) results in thermal blocking and a temperature increase with often devastating consequences for various nuclear and electronic systems (burnout) [3]. Boiling crisis is of fundamental theoretical interest because it marks a transition between two different nonequilibrium dynamic steady states, each incorporating fluid flow, phase transition, contact line dynamics, and heat transfer. It is associated with thermal bi-stability which is observed in many metals and superconductors [4]. Recently experimental indications have been obtained that the CHF is accompanied by 1=f-noise [5,6]. The implied scale free behavior places boiling crisis in the class of nonequilibrium phase transitions in driven systems and links it to such diverse phenomena as earthquakes, plasticity and magnetism [7,8]. Boiling crisis takes place on a contact between a liquid and an adjacent solid surface. The surface generates a heat flux È which produces an overheat ÁT ¼ T À T t , where T t is the saturation temperature. At subcritical values of the heat flux separate bubbles are released from the hot surface (nucleate boiling), while at the CHF È C a vapor layer suddenly covers the hot surface. The transition is akin to the formation of a macro-crack after a prolonged micro-cracking in the continuously loaded brittle system [9]. At supercritical heat flux values, the system stabilizes again in a film boiling regime where bubbles are released from a detached liquid-vapor surface. The value of È C is known to depend on many physical parameters (roughness, wettability, tangential flow, etc.) and the corresponding correlations are thoroughly documented in the engineering literature [1,10,11]. An assumption that the boiling crisis is a manifestation of a macroscale hydrodynamical instability [12] gives the Kutateladze dimensional analysis for È C [13], which is an agreement with most experiments at normal conditions [1]. For a …

Among known hydrodynamic instabilities, a boiling crisis is one of the most mysterious and poorly understood [1,2].It may be described as a sudden loss of liquid contact with a heating surface (dryout) when a small increase of heat flux transforms nucleate boiling into film boiling.Reaching the critical heat flux (CHF) results in thermal blocking and a temperature increase with often devastating consequences for various nuclear and electronic systems (burnout) [3].Boiling crisis is of fundamental theoretical interest because it marks a transition between two different nonequilibrium dynamic steady states, each incorporating fluid flow, phase transition, contact line dynamics, and heat transfer.It is associated with thermal bi-stability which is observed in many metals and superconductors [4].Recently experimental indications have been obtained that the CHF is accompanied by 1=f-noise [5,6].The implied scale free behavior places boiling crisis in the class of nonequilibrium phase transitions in driven systems and links it to such diverse phenomena as earthquakes, plasticity and magnetism [7,8].
Boiling crisis takes place on a contact between a liquid and an adjacent solid surface.The surface generates a heat flux È which produces an overheat ÁT ¼ T À T t , where T t is the saturation temperature.At subcritical values of the heat flux separate bubbles are released from the hot surface (nucleate boiling), while at the CHF È C a vapor layer suddenly covers the hot surface.The transition is akin to the formation of a macro-crack after a prolonged microcracking in the continuously loaded brittle system [9].At supercritical heat flux values, the system stabilizes again in a film boiling regime where bubbles are released from a detached liquid-vapor surface.
The value of È C is known to depend on many physical parameters (roughness, wettability, tangential flow, etc.) and the corresponding correlations are thoroughly documented in the engineering literature [1,10,11].An assumption that the boiling crisis is a manifestation of a macroscale hydrodynamical instability [12] gives the Kutateladze dimensional analysis for È C [13], which is an agreement with most experiments at normal conditions [1].For a detailed explanation of the Kutateladze argument see Sec.I in Supplemental Material [14].Measurements in low gravity conditions [2,15] and a very detailed observations around the heater surface [3] supported an alternative theory that the instability is related to much smaller scales [10] and is controlled by the processes in a thin layer adjacent to the surface, in particular, by the vapor recoil mechanism [15][16][17].Away from the thermodynamic critical point, the microscale theories were shown to be compatible with the Kutateladze formula with the same order of magnitude for the nondimensional constant [17].The concurrence of very different micro and macro descriptions is a sign that the CHF phenomenon may be lacking a characteristic scale.
To test the hypothesis of criticality for the boiling crisis we study in this Letter intermittent acoustic emission (AE) generated at the CHF conditions.The AE sensors have been widely used to monitor boiling activity in containers and pipes [18,19], however the statistics of avalanches has not been studied before.Our experiments show that during boiling crisis the distribution of avalanches have a powerlaw structure which disappears in the bordering nucleate and film boiling regimes.To understand the origin of criticality we developed a simple theoretical description of the boiling crisis in terms of a lattice spin model with nonlocal interactions.Despite drastic simplifications of the hydrodynamical side of the boiling problem our model captures the most important experimental observations including the value of the main critical exponent.
Experiment.-We used the conventional quenching set up [20] and studied boiling during the immersion of an aluminum cylinder into liquid nitrogen [21].The initial temperature of the sample was much higher than the equilibrium temperature (T t ¼ 77 K) causing film boiling.In the process of cooling, the thickness of the vapor film around the sample decreases and near the CHF conditions, the liquidvapor interface starts breaking and hitting the metallic cylinder.At the CHF the intensity of fluctuations increases and subsequently the system enters the nucleate boiling regime with bubbles detaching directly from the metal surface.
For the detection of AE signal we used the device shown in Fig. 1.During the experiment the liquid nitrogen level was placed in the middle of the neck.We located several thermocouples at different heights in the upper bar to monitor the heat flux through the neck.The signal from the transducer was first preamplified (20 dB) and then analyzed by a Mistras PCI-2 system.Figure 1 shows the behavior of the temperature of the head TðtÞ, its time derivative dT=dt, characterizing the heat flux, and a segment of the bare AE signal (voltage) where the spikes represent the intermittent component of the signal while the background represents the noise of hydrodynamic origin [18].The first decreasing part of the TðtÞ curve corresponds to film boiling with low flux values.The abrupt temperature drop indicates a rapid increase in flux which reaches the CHF (boiling crisis) approximately at the maximum slope of the TðtÞ curve.The subsequent slow decay of the temperature is associated with nucleate boiling regime.
The time evolution of the AE activity is presented in Fig. 2 showing a different realization of the experiment.Here we can identify three standard boiling regimes: (a) t < t 1 ¼ 415 s (film boiling), (b) t 1 < t < t 2 ¼ 430 s (boiling crisis), and (c) t > t 3 ¼ 440 s (nucleate boiling).The maximum of acoustic activity (at t ¼ 427:5 s) coincides with the second inflexion point in the TðtÞ curve (indicated again by an arrow).Knowing the specific heat of the sample CðTÞ [22] we computed the heat flux È ¼ CðTÞdT=dt which is shown as inset in Fig. 2 vs the overheat ÁT (Nukiyama curve [21,23]).The CHF È c ¼ 2:5 W=cm 2 agrees with Kutateladze formula [13] for liquid nitrogen.See Sec.I in Supplemental Material [14] for a more detailed comparison.
The statistics of avalanches is presented in Fig. 3(a) (More detailed graphs can be found in Sec.II of Supplemental Material [14]).In the film boiling regime (red curve), where the distribution is almost flat, the con-siderable presence of high energy events reflects the hits of the sample by the detached liquid-gas interface.In the second regime (black curve), which we associate with the boiling crisis, one can see the emergence of the power-law statistics with the exponent $ 2:1 AE 0:1 obtained by maximum likelihood (ML) fit [24].This behavior disappears in the third regime (blue curve), where the number of hits with high energy is small reflecting conventional gas bubble formation on a hot surface.The universality of the power-law distribution in the crisis regime was tested under different vessel sizes, sample sizes, and materials (Al, Cu); in all cases the exponent remained close to the reported value (see Sec. II in Supplemental Material [14]).
The avalanche distribution in the interval t 2 < t < t 3 resembles nucleate boiling while also exhibiting some To test the scale free behavior of the system further we performed a spectral analysis of the acoustic signal shown in Fig. 1.The averaged power spectrum (PS) in the frequency range 0.001-10 kHz is shown in Fig. 3(b).Higher frequencies are omitted since they are crucially affected by the transducer response.The PS in the film regime (red curve) shows a white noise at low-frequencies revealing uncorrelated events involving the whole system [25].At high frequencies we see the behavior which is qualitatively similar to what has been observed in other hydrodynamic systems [26].At the boiling crisis regime (black curve) the PS exhibits a characteristic 1=f noise in the range 0.001-1 kHz with $ 0:5 (cf.[6]).In nucleate boiling regime (blue curve) we again observe an almost uncorrelated signal at low frequencies indicating statistical independence of the individual bubble formation events.
Modeling.-Thephenomenon of a boiling crisis involves fluid motion, heat transfer, and phase transition.Our model is based on a detailed description of the heat transfer with a schematic representations of the fluid motion (prescribed kinematic flow) and of the phase transition (spin model).In view of the expected scale free character of the boiling crisis regime, we aim at capturing the universality class and do not attempt to describe the structure of the turbulent fluid motion in full detail.
In terms of the velocity field vðr; tÞ, the temperature field Tðr; tÞ, and the phase field ðr; tÞ 2 ½0; 1 we can write the heat transfer equation accounting for a phase transition in the form cdT=dt þ Ld=dt ¼ kÁT, where L is the latent heat, is the fluid density, c is the specific heat, k is thermal conductivity and d=dt ¼ @=@t þ vr.In dimensionless variables " r ¼ r=a, " t ¼ tv=a, " T ¼ Tc=L, " ¼ , where v is the velocity scale and a is the spatial scale, the problem depends only on the Peclet number Pe ¼ vac=k.Considering typical bubble sizes of a $ 10 À3 m, velocities of v $ 1 m=s and physical parameters for N 2 , we obtain Pe $ 10 2 -10 4 .This means that advection occurs much faster than conduction, which enables us to separate the two processes and solve the heat balance equation in two steps.During the fast advection step we need to solve d " T=d " t þ d " =d " t ¼ 0, while during the slow conductive step the process is governed by @ " T=@ " þ @ " =@ " ¼ Á " T where we introduced slow time " ¼ " t=Pe.To complete the formulation we need additional assumptions regarding the velocity field " vðr; tÞ and the phase field " ðr; tÞ.We assume that " is a spin variable taking values " ¼ 0 in the liquid phase and " ¼ 1 in the gas phase.We also assume that the kinetics of our phase transition is instantaneous and that the spin field " is enslaved to the field T through a hysteretic, thresholdbased relation.The origin of the hysteresis is an assumption that the transformation always takes place at the fast time scale and is therefore adiabatic.To compute the size of the hysteresis we assume that after the latent heat is released or absorbed the new phase is at transition temperature T t .Then, the liquid must transform if heated above T L!G ¼ T t þ L=ðc L Þ while the gas transforms if cooled below T G!L ¼ T t À L=c G (see insert in Fig. 4).
We discretize our equations on a 3D cubic lattice with the cell scale a representing the linear size of the smallest bubble.We assume that the gas cells advance along y direction with constant velocity v.During the advection step the spin value is preserved d " =dt ¼ 0 and therefore the temperature is also preserved d " T=dt ¼ 0. In the course of the advance some gas particles will overlap with liquid particles then vacating lattice spots forming ''voids.''The overlapped liquid fills the voids while acquiring the average temperature of the neighboring liquid cells.For a more detailed description of the model, see Sec.III in Supplemental Material [14].
In our numerical experiments the temperature of the bottom surface of the domain in contact with a heater was taken equal to T m with superimposed quenched uncorrelated Gaussian fluctuations while the upper plane was assumed to be in contact with a liquid reservoir at the saturation temperature T t .In the horizontal direction the boundary conditions were taken to be periodic.To reproduce the experimental data we chose the modeling parameters matching liquid nitrogen at $77 K and normal pressure.The disorder was characterized by standard deviation 0.1 in accordance with experimental observations [3,10] and the results were averaged over different realizations of the disorder.
Our numerical simulations produced a realistic shape of the Nukiyama curve (Fig. 2) and generated statistics of avalanches supporting the idea of criticality (Fig. 4).We associate avalanches with the intermittent formation of dry spots represented by the 2D contact area of the gas bubbles adjacent to the hot surface (Fig. 5).The statistics of dry spots, averaged over different time steps, is shown in Fig. 4 for several values of ÁT around the CHF conditions; to avoid size effects we discarded the spanning avalanches [27].One can see that at low temperatures, the distribution is exponentially damped whereas at high temperatures there is a significant contribution of large size events.The analysis of the Kolmogorov-Smirnov distance [24,28] renders an optimal power-law distribution for ÁT=T t ¼ 0:27 with the ML exponent ' 2:1 AE 0:1.A more extensive analysis of the distributions of avalanches, size effects and exponent estimation can be found in Sec.IV of Supplemental Material [14].
The overall parametric dependence of the pdf and the value of the exponent are therefore in excellent agreement with our experimental observations.Similar value of the exponent (Fisher exponent) has been also found in a number of percolation models [29], which suggests that our 3D model can be simplified even further and reformulated as a problem of finding a percolating cluster in an effective 2D model.
Conclusions.-Theexperimental study of AE during boiling crisis established for the first time a power-law statistics of avalanches near the point where the heat flux reaches its maximum.A simple spin automaton model reproduced the basic features of the experiment and captured the value of the main critical exponent.The model suggests that the origin of criticality can be linked to the phenomena close to the heating surface involving a competition between the percolative surface advance of dry spots and the bubble buoyancy.In particular, this means that the hydrodynamic flow far from the hot surface is of secondary importance for the critical dryout.

FIG. 1 (FIG. 2 (FIG. 3 (
FIG. 1 (color online).(color online) AE signal and temperature trace TðtÞ during the cooling of the system.The arrow indicates the maximum slope of the TðtÞ curve.The inset shows a magnification of the AE signal and the 45 dB ( ¼ 1:77 mV) threshold used for the statistical analysis.The scheme on the right shows the sample (the cylindrical head (radius 9 mm, height 15 mm) preceded by a neck (radius 3 mm height, 34 mm)), the position of the broadband PZT piezoelectric transducer (Tr) and the Pt-100 resistor.Lower panel: time derivative of TðtÞ.

FIG. 4 (
FIG. 4 (color online).(color online)The computed dry spot size distribution at different degrees of overheating ÁT=T t .The system size is 100 Â 100 Â 10.The inset shows the adiabatic threshold rule for phase transition.

25 FIG. 5 .
FIG.5.Phase configurations in a vertical plane (upper panels) and in the horizontal layer closest to the hot surface (lower panels) for 3 different temperatures.White areas correspond to gas (e.g.dry spots), whereas areas black correspond to liquid.The system size is 100 Â 100 Â 20.