3.1. Bubble appearance and growth

Figure 3 shows a series of phase-detection (labelled PD) images during flow boiling at 0.2 MPa, a relatively low saturation pressure, obtained from the study by Kossolapov (Reference Kossolapov2021). As revealed by these images, and also shown by Sinha, Narayan & Srivastava (Reference Sinha, Narayan and Srivastava2022), when a vapour bubble grows in low-pressure conditions, it is dragged and rolled in the direction of the flow, creating an asymmetric microlayer (shown as white fringes) and dry spot (shown in solid white). This asymmetric shape leads to contact line surface tension forces that inhibit bubble sliding. This bubble asymmetry can also be seen in the side view images, where there is a clear difference in the advancing and receding contact angles. These observations are in striking contrast to the images of vapour bubbles at high pressures collected in this study and shown in figure 4. At high pressures, we did not observe a microlayer, the optical footprint is rather circular and the physical footprint appears to be close to the centre of the bubble, indicating spherical symmetry. Even at the lowest pressure (i.e. 1.0 MPa), the physical footprint is located near the centroid of the optical footprint up until the bubble leaves the nucleation site. The physical size of the bubbles at the moment of departure from the nucleation site also appears to be very different. For instance, the equivalent departure radius of vapour bubbles at 0.2 MPa is approximately 0.25 mm, while the departure radii of bubbles at high pressures are well below 0.05 mm. The assumption of spherical symmetry can be validated by comparing the magnitude of surface tension and other forces (i.e. buoyancy, viscous and inertial forces) acting on the bubble as pressure increases. Consider the Bond, capillary and Weber numbers given by

(3.1)\begin{gather} Bo = \frac{(\rho_{{l}}-\rho_{{v}})gD_{{d}}^{2}}{\sigma}, \end{gather}

(3.2)\begin{gather} Ca = \frac{\mu_{{l}} U_{{l}}}{\sigma}, \end{gather}

and

(3.3)\begin{equation} We = \frac{\rho_{{l}} U_{{l}}^{2}D_{{d}}}{\sigma}, \end{equation}

respectively, where $\rho$ is the density, $\sigma$ is the liquid-vapor surface tension, and $g$ is the acceleration of gravity. The characteristic length and velocity scales are taken to be the departure diameter $D_{d}$ and local liquid velocity $U_{l}$, respectively. From 0.2 to 1.0 MPa, for example, surface tension decreases by 30 %, liquid density by 6 %, liquid viscosity by 60 %, while vapour density increases by 78 %. Ultimately, the most important effect is the decrease in the bubble departure diameter, which drops by a factor of five, resulting in the Bond number decreasing by a factor of 20. With pressure, the characteristic velocity $U_{{l}}$ will drop according to the bubble diameter; therefore, the capillary number and Weber number will also decrease. Other factors that may cause spherical asymmetry, such as shear (Taylor Reference Taylor1934) or turbulence deformations (Ni Reference Ni2024), should also be negligible in our conditions. Overall, dimensional arguments and experimental observations both indicate that bubbles should become more spherical as pressure increases. Further discussion and more quantitative justification of these assumptions can be found in Appendix A.

Figure 3. Phase-detection (PD) and shadowgraphy side view images of flow boiling at 0.2 MPa. The dashed orange line represents the location of a nucleation site. The solid blue line shows the vertical displacement of a bubble over time.

Figure 4. Phase-detection images of flow boiling at different pressures. The dashed orange line represents the location of a nucleation site. The solid blue line shows the vertical displacement of a bubble over time.

To model bubble departure and sliding, we must account for the change in bubble size with time. Therefore, we measure the bubble radius over time as it grows and slides on the surface. Figures 5–7 summarize the statistical distribution of bubble sizes with time for each operating pressure, coloured as blue, yellow and red for 1, 2 and 4 MPa, respectively, where the last bubble size measurement occurs when it leaves the high-speed camera field of view. Data are grouped by mass flux: 500 kg m$^{-2}$ s$^{-1}$ in figure 5, 1000 kg m$^{-2}$ s$^{-1}$ in figure 6 and 1500 or 2000 kg m$^{-2}$ s$^{-1}$ in figure 7. At each frame, or time step, the equivalent bubble radius is visualized using a box and whisker plot, alongside dots to show any statistical outliers. To the right of the statistical distributions are samples of individual bubble growth histories, depicted using the same colour scheme noted above.

Figure 5. Bubble growth statistics (a,c,e) and selected histories (b,d,f) from 1 to 4 MPa at a mass flux of 500 kg m$^{-2}$ s$^{-1}$.

Figure 6. Bubble growth statistics (a,c,e) and selected histories (b,d,f) from 1 to 4 MPa at a mass flux of 1000 kg m$^{-2}$ s$^{-1}$.

Figure 7. Bubble growth statistics (a,c,e) and selected histories (b,d,f) from 1 to 4 MPa at the highest mass fluxes tested (2000 kg m$^{-2}$ s$^{-1}$ at 1 MPa and 1500 kg m$^{-2}$ s$^{-1}$ at 2–4 MPa).

Unfortunately, due to the of lack wall and local fluid temperature measurements, we cannot formulate a model of the bubble growth. Instead, we provide a simple but physically motivated and tractable semi-empirical correlation that will be used as an input in our mechanistic departure and sliding models. The measurements show that the equivalent bubble radius is mostly proportional to the square root of time, and can be fitted to the equation

(3.4)\begin{equation} R_{{b}} = C_{{RB}}\sqrt{t}, \end{equation}

where $C_{{RB}}$ is an empirical constant. The value of $C_{{RB}}$ is chosen so that the mean absolute error between (3.4) and the instantaneous bubble radius distributions, shown as box and whisker plots in figures 5–7, is minimized. Equation (3.4) is shown as a solid black line over the experimental data in figure 5(a,c,e). Dashed black lines corresponding to an uncertainty in $C_{{RB}}$ equal to ${\pm }35\,\%$ are also shown for illustrative purposes. We note that it could be possible to correlate bubble growth with other power laws (e.g. $R_{b} = a t^{b}$), but there would be no major improvement in the prediction accuracy compared with the square root of time fit. Instead, a square root of time dependency is physically motivated and hints toward a heat-diffusion-controlled growth. As discussed by Mikic, Rohsenow & Griffith (Reference Mikic, Rohsenow and Griffith1970), bubble growth can be limited by two mechanisms, namely inertia and heat diffusion from the superheated liquid. Since the time scale associated with inertial-controlled bubble growth is of the order of nanoseconds for our condition, we expect bubble growth to be heat-diffusion controlled, which should result in the bubble radius increasing with the square root of time. Other effects, such as the translational motion of the bubble, can increase the growth rate (Legendre, Borée & Magnaudet Reference Legendre, Borée and Magnaudet1998), but such effects are only pronounced at high bubble Reynolds numbers $(Re_{b} = \frac{\rho_{l} \left|U_{l} - U_{b}\right|}{\mu_{l}})$, which is not the case in our conditions. With these considerations in mind, we refer to the work of Plesset & Zwick (Reference Plesset and Zwick1954) to provide an interpretation of the fitting coefficient $C_{RB}$. According to Plesset & Zwick (Reference Plesset and Zwick1954), the heat-diffusion-controlled growth of the bubble radius in a uniformly superheated liquid can be predicted by

(3.5)\begin{equation} R_{{b}} \approx B \sqrt{t} = \sqrt{\frac{12\alpha_{{l}}}{\rm \pi}} \frac{\rho_{{l}}}{\rho_{{v}}}\frac{c_{{p,l}}}{h_{{lv}}} (T_{{w}} - T_{{sat}})\sqrt{t}. \end{equation}

Note that our $C_{{RB}}$ is equivalent to the theoretical parameter $B$ in (3.5). Since we could not measure the nucleation temperature $T_{{w}}$, it is impossible to evaluate the parameter $B$ experimentally. However, we can estimate $(T_{{w}} - T_{{sat}})$ based on the correlation by Jens & Lottes (Reference Jens and Lottes1951) and the measured wall heat flux. The empirical $C_{{RB}}$, bubble radius mean absolute error and theoretical $B$ are reported in table 2, where, on average, the ratio $C_{{RB}}/B$ is $0.55 \pm 0.21$. Our bubbles grow slower than predicted by (3.5). This is expected, as (3.5) was obtained for a spherical bubble growing in a uniformly superheated liquid, while our bubbles grow in the presence of subcooled liquid, which slows the growth down.

Table 2. Summary of empirical $C_{{RB}}$ measurements and their comparison against theoretical predictions.

It can be seen both in figures 5–7 that the heat-diffusion-controlled growth fitting (3.4) is less accurate at 1 MPa. This error is typically largest long after the bubble has nucleated, i.e. near the end of when bubbles are tracked. This is because the bubble growth rate at 1 MPa starts to decrease and it eventually appears to stop growing. Considering the fact that bubbles at 1 MPa are significantly larger than bubbles at 2 and 4 MPa, it is plausible that their size makes effects of subcooling more substantial, preventing the bubbles from growing further. The growth statistics for 2 and 4 MPa, on the other hand, generally demonstrate better agreement, but contain more outliers. This is primarily because many more bubbles can be tracked at higher pressures due to their smaller size and shorter time near the nucleation site; therefore, more outliers will inevitably be measured. However, the relative number of outliers is small, as illustrated in the sample of the red growth histories in figures 5–7, where many of the individual bubbles follow the expected growth trend quite well. On average, (3.4) fits the bubble size measurements within 27 % across all operating conditions and, overall, measurements indicate that bubble growth is consistently proportional to the square root of time (i.e. heat-diffusion controlled), particularly during the times when the bubble is still near the nucleation site. However, while we acknowledge that our treatment of bubble growth is not complete, and a more precise and mechanistic growth model will be an important goal of future studies, (3.4) provides a reasonable boundary condition for the sliding and departure models, which is the focus of this work.

In passing, we also observe that, for a given mass flux, an increase of pressure decreases the bubble growth rate. For instance, when increasing the pressure from 1 to 2 MPa, the liquid–vapour density ratio decreases by a factor of two, while the thermal diffusivity and Jakob number $Ja=c_{{p,l}} (T_{{w}} - T_{{sat}})/h_{{lv}}$ only decrease by 3 % and 6 %, respectively. The Jakob number does not change as much because the reduction of latent heat is nearly balanced out by the reduction in nucleation temperature due to decreased surface tension. Therefore, the primary reason why $C_{{RB}}$ decreases with increasing pressure is because of the increase in vapour density, which explains the slower bubble expansion as liquid evaporates.

Mass flux, on the other hand, does not have an obvious effect on the bubble size, which supports the numerical findings of Legendre et al. (Reference Legendre, Borée and Magnaudet1998), but it of course affects the local liquid velocity near the bubble. The secondary axes of the plots in figures 5–7 show the equivalent bubble radius in wall units, where the friction velocity is estimated using the friction factor correlation given by McAdams (Todreas & Kazimi Reference Todreas and Kazimi2011) (see Appendix A for details). Note that, as mass flux increases, the thickness of the viscous sublayer decreases, the local fluid velocity near the bubble increases and creates stronger detaching forces that promote bubble sliding.

3.2. Bubble sliding

We now revisit figures 3 and 4 to discuss the observed sliding phenomenon in detail. It can be seen that the bubble continues to adhere to the boiling surface as it slides vertically upward in the direction of the bulk flow. This is evidenced by presence of a bright patch, indicating a dry spot, in the centre of the bubble optical footprint. At 0.2 MPa, bubbles do not physically lift off from the boiling surface until they are far away (0.75 mm, ${\sim }1.5 D_{{d}}$) from their original nucleation site (Kossolapov Reference Kossolapov2021). At higher operating pressures, the distance at which bubbles lift off from the boiling surface is even larger (relative to the bubble size). This is expected because the bubble asymmetry at 0.2 MPa should yield contact line surface tension forces that reduce bubble acceleration in the vertical direction. One can also see an increase in the number of visible bubbles at a given mass flux as pressure increases from 1 to 2–4 MPa. This behaviour can be explained qualitatively, as increasing pressure decreases surface tension, which in turn decreases the nucleation temperature and wait time. The vertical displacement of individual bubbles is illustrated with solid blue lines in figure 4. It can be seen that the bubble trajectory is nonlinear with time, as the bubble accelerates as it grows. One explanation for this observation is that a larger bubble is exposed to larger local fluid velocities, which aids in accelerating the bubble away from its original nucleation site position. To better understand sliding from a mechanistic point of view, a momentum balance is written for a growing vapour bubble that is surrounded by a forced flow of liquid, shown schematically in figure 6. By assuming axisymmetric growth of a spherical bubble with negligible momentum transfer due to evaporation, the momentum balance can be written as (Bucci Reference Bucci2020)

(3.6)\begin{align} \rho_{{v}}V_{{b}}\frac{{\rm d} \boldsymbol{U}_{b}}{{\rm d}t} = \unicode{x2230}_{V_{{b}}}{\rho_{{v}} \boldsymbol{g}} \,{\rm d}V + \oint_{CL} \sigma \boldsymbol{t}\,{\rm d}\ell - \iint_{S_{{b}}} p_{{v}} \boldsymbol{n}\,{\rm d}S -\iint_{S_{{lv}}} p_{{l}} \boldsymbol{n}\,{\rm d}S + \iint_{S_{{lv}}} \boldsymbol{\tau} \boldsymbol{{\cdot}} \boldsymbol{n}\,{\rm d}S, \end{align}

where $V_{{b}}$ is the volume bubble, p is pressure, $\boldsymbol{\tau}$ is the viscous stress, $\boldsymbol {U}_{b}$ is the bubble velocity, $CL$ is the contact line and the meaning of the other symbols is obvious from figure 8. We have assumed that the only body force is the weight of the bubble, and the only major surfaces forces are the contact line surface tension forces, viscous stresses along the liquid–vapour interface and pressure forces from the liquid and vapour phase. Because the average flow field and gravity are in the vertical direction, only the vertical direction of the momentum balance is needed to predict sliding. This simplification coupled with some common algebraic manipulations yields

(3.7)\begin{align} \rho_{{v}}V_{{b}}\frac{{\rm d}U_{b}}{{\rm d}t}=\oint_{CL}\sigma\boldsymbol{t}\,{\rm d}\ell\boldsymbol{{\cdot}} \boldsymbol{i} + (\rho_{{l}}-\rho_{{v}})gV_{{b}} -\iint_{S_{{lv}}}(p_{{h}}-p_{{c}})\boldsymbol{n}\,{\rm d}S\boldsymbol{{\cdot}} \boldsymbol{i} + \iint_{S_{{lv}}}\boldsymbol{\tau}\boldsymbol{{\cdot}} \boldsymbol{n}\,{\rm d}S\boldsymbol{{\cdot}} \boldsymbol{i}, \end{align}

where the forces (on the right-hand side) from left to right are the contact line surface tension force, buoyancy force and total hydrodynamic force (given by the third and fourth terms). The quantity $p_{{h}}$ is the hydrodynamic pressure of the liquid, and $p_{{c}}$ is a reference pressure, taken to be the liquid pressure at the bubble base. This is done to algebraically shift the effect of hydrostatics from the liquid pressure surface integral to the buoyancy force term (given by the second term in (3.7)) (Bucci Reference Bucci2020). Because most of the bubbles that are observed during sliding appear to be spherical and symmetric, the net contact line surface tension force can be neglected. The total hydrodynamic force (i.e. the last two terms on the right-hand side of (3.7)) is modelled as the combination of a quasi-steady drag force and a virtual (or added) mass force (Favre et al. Reference Favre, Colin, Pujet and Mimouni2023). Quasi-steady drag of a bubble inside a turbulent shear flow near a wall near our conditions has yet to be measured or simulated. However, drag models for rigid spheres and clean bubbles within linear shear flows bounded by a wall (i.e. without turbulence) have been extensively studied and match our operating conditions well because most bubbles remain in the viscous sublayer (see figure 22a in Appendix A). We also note that turbulence typically does not have a systematic effect on mean drag forces (Wu & Faeth Reference Wu and Faeth1994; Sridhar & Katz Reference Sridhar and Katz1995; Bagchi & Balachandar Reference Bagchi and Balachandar2003; Salibindla et al. Reference Salibindla, Masuk, Tan and Ni2020), while the presence of a wall does (Zeng et al. Reference Zeng, Najjar, Balachandar and Fischer2009; Shi et al. Reference Shi, Rzehak, Lucas and Magnaudet2020, Reference Shi, Rzehak, Lucas and Magnaudet2021); therefore, a model based on measurements or simulations of bubble drag within a linear shear flow bounded by a wall should sufficiently capture drag on the bubbles we observe. Additional effects that may impact drag, such as bubble interactions and wake effects, should also be negligible due to our low bubble Reynolds numbers (often below 10, see figure 22b). Typically, wakes from flow over rigid spheres do not appear until ${{Re}}_{b}>20$, while vortex shedding does not occur until approximately ${{Re}}_{b} > 200$ (Clift, Grace & Weber Reference Clift, Grace and Weber2005). More recent studies have also confirmed this fact for flows over spheres in turbulent environments. Wu & Faeth (Reference Wu and Faeth1995), for example, did not observe vortex shedding until ${{Re}}_{b}<300$ and note that velocity fluctuations are the same as the surrounding turbulent fluctuations (i.e. without bubbles) for ${{Re}}_{b}<300$. This indicates that there should not be severe distortions in the flow field for upstream bubbles, provided their ${{Re}}_{b}$ values are low. All this allows us to assume that the drag force on our bubbles and the velocity field around them are approximately the same as the case of a single isolated bubble. Substituting these assumptions into (3.7) yields

(3.8)\begin{equation} \rho_{{v}}V_{{b}}\frac{{\rm d}U_{b}}{{\rm d}t}= (\rho_{{l}}-\rho_{{v}})gV_{{b}} + \frac{1}{2}\rho_{{l}}S_{{p}}C_{{D}}(U_{l}-U_{b})|U_{l}-U_{b}| + F_{{AM,x}}, \end{equation}

where $S_{p}$ is the projected area of the bubble (i.e. ${\rm \pi} R_{{b}}^{2}$), $C_{{D}}$ is the drag coefficient and $F_{{AM,x}}$ is the added mass force. The drag coefficient can be calculated using the model by Mazzocco et al. (Reference Mazzocco, Ambrosini, Kommajosyula and Baglietto2018), which is the drag coefficient for a solid sphere (Zeng et al. Reference Zeng, Najjar, Balachandar and Fischer2009) corrected to account for the different boundary conditions between a solid sphere and vapour bubble (Legendre & Magnaudet Reference Legendre and Magnaudet1997)

(3.9)\begin{equation} C_{{D}} = 1.13\frac{24}{{{Re}}_{{b}}}\left(1+0.104{{Re}}_{{b}}^{0.753}\right), \end{equation}

with the bubble Reynolds number defined as

(3.10)\begin{equation} {{Re}}_{{b}} =\frac{\rho_{{l}}|U_{l}-U_{b}|D_{{b}}}{\mu_{{l}}}. \end{equation}

Finally, we note that thermally induced surface tension gradients (i.e. Marangoni forces) are of little concern in our conditions, especially compared with Marangoni forces from surface contaminants (Clift et al. Reference Clift, Grace and Weber2005). We note that, if the bubbles were contaminated, then the bubble interface would be a rigid boundary with a no-slip boundary condition, making bubbles behave more like rigid spheres (Harper, Moore & Pearson Reference Harper, Moore and Pearson1967), in which case the drag coefficient by Zeng et al. (Reference Zeng, Najjar, Balachandar and Fischer2009) may be used. To elucidate these effects, we compare the Mazzocco et al. (Reference Mazzocco, Ambrosini, Kommajosyula and Baglietto2018) and Zeng et al. (Reference Zeng, Najjar, Balachandar and Fischer2009) correlation in Appendix A and observe that the change in bubble sliding predictions are nearly negligible. The liquid velocity at the bubble centroid $U_{l}$ is calculated using the Van Driest (Reference Van Driest1956) eddy-diffusivity model and McAdams’ friction factor correlation (Todreas & Kazimi Reference Todreas and Kazimi2011) (more details of the liquid flow analysis are provided in Appendix A). The distance from the wall at any given point in time is taken as the bubble radius $R_{b}$, because the bubble will remain adhered to the wall during sliding. To model the added mass force, we use the approach and results by Favre et al. (Reference Favre, Colin, Pujet and Mimouni2023) and van der Geld (Reference van der Geld2009), which consider the effect of a nearby wall on the flow field surrounding the bubble. The vertical component of the added mass force then becomes

(3.11)\begin{equation} F_{{AM,x}} = C_{{AM,x}}\rho_{{l}}V_{{b}}\left(3\frac{\dot{R}_{{b}}}{R_{b}}-\frac{{\rm d}U_{b}}{{\rm d}t}\right), \end{equation}

where $C_{{AM,x}}=0.636$. The second term in (3.11) represents added inertia due to displacing the liquid as the bubble moves, while the first term actually promotes bubble sliding if the liquid velocity is faster than the vapour bubble, which is expected to be the case. With these closure relations in mind, (3.8) can be written as the following nonlinear ordinary differential equation:

(3.12) $$\begin{align} (\rho_{{v}} + C_{{AM,x}}\rho_{{l}})V_{{b}}\frac{{\rm d}U_{b}}{{\rm d}t} &=\underbrace{(\rho_{{l}}-\rho_{{v}})gV_{{b}}}_{F_{{B}}} + \underbrace{\frac{1}{2}\rho_{{l}}S_{{p}}C_{{D}}(U_{l}-U_{b})|U_{l}-U_{b}|}_{F_{{D}}} \nonumber\\ &\quad+ \underbrace{3C_{{AM,x}}\rho_{{l}}V_{{b}}(U_{l}-U_{b})\frac{\dot{R}_{{b}}}{R_{b}}}_{F_{{A}}}. \end{align}$$

This equation implies that bubble acceleration will remain positive during the growth phase, and that there are three terms that contribute to bubble sliding. In an effort to simplify this expression, (3.12) can be adjusted to solve for the bubble acceleration

(3.13)$$\begin{gather} \frac{\rho_{{v}} + C_{{AM,x}}\rho_{l}}{\rho_{l}}\frac{{\rm d}U_{b}}{{\rm d}t}\nonumber\\ =\left(1-\frac{\rho_{v}}{\rho_{l}}\right)g + \frac{3}{8}C_{D}\frac{|U_{l}-U_{b}|}{R_{b}}(U_{l}-U_{b})+3C_{{AM,x}}(U_{l}-U_{b})\frac{\dot{R}_{{b}}}{R_{b}}. \end{gather}$$

Next, it will be assumed, and later justified, that ${{Re}}_{b} < 10$ for most cases, implying

(3.14)\begin{equation} C_{D} \approx \frac{27.12}{{{Re}}_{b}}, \end{equation}

which yields the following simplification:

(3.15) $$\begin{align} \frac{\rho_{{v}} + C_{{AM,x}}\rho_{l}}{\rho_{l}}\frac{{\rm d}U_{b}}{{\rm d}t}&=\left(1-\frac{\rho_{v}}{\rho_{l}}\right)g + \frac{3}{16}\frac{27.12\mu_{l}}{\rho_{{l}}R_{b}^{2}}(U_{l}-U_{b})\nonumber\\ &\quad +3C_{{AM,x}}(U_{l}-U_{b})\frac{\dot{R}_{{b}}}{R_{b}}. \end{align}$$

Then, by relating the fact that $R_{b} = C_{{RB}}\sqrt {t}$, we obtain

(3.16)\begin{equation} \frac{\rho_{{v}} + C_{{AM,x}}\rho_{l}}{\rho_{l}}\frac{{\rm d}U_{b}}{{\rm d}t}= \left(1-\frac{\rho_{v}}{\rho_{l}}\right)g + \left(\frac{81.36\nu_{l}}{16C_{RB}^{2}}+3C_{{AM,x}}\right)\frac{U_{l}-U_{b}}{t}. \end{equation}

Finally, the bubble equation of motion is scaled to wall units, where the velocity is scaled by the friction velocity $U_\tau =\sqrt {\tau _{w}/\rho _{l}}$, and the length is scaled by $\nu _{l}/U_\tau$. Because the bubble remains adhered to the wall and is relatively small in size during sliding, this scaling procedure is both convenient and physically meaningful. Using $y$ to denote the distance from the wall, where $y=R_{b}=C_{{RB}}\sqrt {t}$, one obtains

(3.17)\begin{align} \frac{{\rm d}U_{b}^{+}}{{{\rm d} y}^{+}}=\underbrace{\frac{2\rho_{l}}{\rho_{{v}} + C_{{AM,x}}\rho_{l}}\left(\frac{81.36\nu_{l}}{16C_{RB}^{2}} +3C_{{AM,x}}\right)}_{\varPi_{1}}\frac{U_{l}^{+}-U_{b}^{+}}{y^{+}} +\underbrace{\frac{2g\nu_{l}}{U_\tau^{3}}\frac{1-\rho_{v}/\rho_{l}}{C_{{AM,x}} +\rho_{v}/\rho_{l}} \sqrt{\frac{\nu_{l}}{C_{RB}^{2}}}}_{\varPi_{2}}. \end{align}

Because the liquid velocity profile is strictly a function of $y^+$, the differential equation can be solved analytically with an integration factor, yielding the exact solution

(3.18)\begin{equation} U_{b}^{+}=\frac{\displaystyle\int_{0}^{y^{+}}\varPi_{1}y^{+^{\varPi_{1}-1}}U_{l}^{+}+\varPi_{2}y^{+^{\varPi_1-1}}\, {{\rm d} y}^{+}}{y^{+^{\varPi_{1}}}} \end{equation}

which can be readily solved using a conventional eddy-diffusivity model. It is noted that the gravitational component in this equation is very small for the operating conditions here, i.e. $\varPi _2 \ll 1$, and so it is neglected for the sake of simplicity

(3.19)\begin{equation} U_{b}^{+}=\frac{\displaystyle\int_{0}^{y^{+}}\varPi_{1}y^{+^{\varPi_{1}-1}}U_{l}^{+}\,{{\rm d} y}^{+}}{y^{+^{\varPi_{1}}}} ,\end{equation}

which results in a function that is dependent on just one non-dimensional number, $\varPi _1$. A phenomenological interpretation of $\varPi _1$ can be understood by considering how drag, added mass and inertial forces change as bubbles get smaller. Precisely, smaller bubbles will have less inertia but a larger area-to-volume ratio and drag coefficient, as shown in (3.13). This results in larger bubble acceleration, which makes them more rapidly approach the velocity of the surrounding liquid. Therefore, $\varPi _1$ as defined in (3.17) can be viewed as a ratio of drag (first term in parentheses) and added mass forces (second term in parentheses) to bubble and liquid inertial forces (first and second terms in the denominator outside the parentheses), where bubbles with larger values of $\varPi _1$ will more closely match the liquid velocity, all else being equal. The simplified analytical model for bubble velocity (i.e. (3.19)) is compared against all measured bubble velocity data in figure 9 (black solid line), where the experimentally measured bubble velocities and radii (represented as distance from the wall) for all operating conditions are plotted in wall units. In plotting (3.19), the quantity $\varPi _1$ is held constant at 5.8 as it falls within the typical spread of measured $\varPi _1$ values (3.3–8.3). We also plot the theoretical profile for the liquid velocity (black dashed line) according to Van Driest (Reference Van Driest1956), which in the limiting case of short distance from the wall is linear, whereas at far distances from the wall it is logarithmic. Notably, (3.19) can be further simplified in a few special cases. When the bubble always remains in the viscous sublayer while sliding, the solution becomes (shown as a red dotted line in figure 9)

(3.20)\begin{equation} U_{b}^{+}=\frac{\varPi_1}{\varPi_1+1}y^+=\frac{\varPi_1}{\varPi_1+1}U_{l}^+ . \end{equation}

It can be seen from the data that the measured bubble velocity follows this linear trend (i.e. equation 20) up until values of $y^+ \sim 10$ are reached. Once the bubble leaves the viscous sublayer and gets into the turbulent core, where $U_{l}^{+}=({1}/{\kappa })\ln y^+ + C^+$, then the bubble velocity may be expressed as (shown as a blue dotted line in figure 9)

(3.21)\begin{equation} U_{b}^{+}=U_{l}^{+} + \frac{\varPi_1}{\varPi_1+1}\frac{y_{{n}}^{+^{\varPi_1+1}}}{y^{+^{\varPi_1}}}-C^{+}\frac{y_{{n}}^{+^{\varPi_1}}}{y^{+^{\varPi_1}}}-\frac{1}{\kappa}\left(\frac{y_{{n}}^{+^{\varPi_1}}\ln y_{{n}}^+}{y^{+^{\varPi_1}}}+\frac{1}{\varPi_1} - \frac{y_{{n}}^{+^{\varPi_1}}}{\varPi_{1}y^{+^{\varPi_1}}} \right), \end{equation}

where $\kappa =0.41$, $C^+=5.0$ and $y_{n}^+$ is the location at which the log law and linear velocity profile intersect (${\sim }10.80$). The slip ratio (i.e. ratio of bubble velocity to liquid velocity) is a constant value less than unity when the bubble is in the sublayer (see (3.20)). Far from the wall, as $y^+$ approaches infinity, (3.21) becomes

(3.22)\begin{equation} U_{b}^{+}= U_{l}^{+}-\frac{1}{\kappa\varPi_1}, \end{equation}

and thus the slip ratio tends to unity. Overall, the measured bubble velocities can be explained by the analytical solutions proposed here, further suggesting that the high-pressure bubble dynamics can be properly accounted for without considering other forces related to bubble asymmetry, which are difficult to measure in practice and add significant complexity to bubble departure modelling. We acknowledge that this analysis assumes that the velocity profile is unaffected by the presence of bubbles, as mentioned before, and re-emphasize that this concern is alleviated in the limiting case of small bubbles and $U_{{b}}\sim U_{{l}}$ (i.e. low ${{Re}}_{b}$ and negligible slip), which is the case in high-pressure flow boiling conditions. The analytic bubble velocity expressions can also be used to calculate the vertical displacement of the bubble. For example, if the bubble remains in the viscous sublayer (i.e. $U_{l}^+=y^+$), then from (3.4) and (3.20), we get

(3.23)\begin{equation} x_{b}(t)=\frac{2}{3}\frac{\varPi_1}{\varPi_1+1}\frac{C_{RB}^{2}U_\tau^2}{\nu_{l}}t^{3/2}. \end{equation}

Otherwise, the bubble displacement can be computed using a specific eddy-diffusivity model and relatively simple numerical integration

(3.24)\begin{equation} x_{b}(t)=\int_{0}^{t}U_{b}^{+}U_\tau \,{\rm d}t. \end{equation}

The combination of (3.20), (3.21) and (3.24) (taking $y_{n}^+=10.8$) results in a completely analytical expression for the bubble trajectory during sliding. Note that, in both (3.23) and (3.24), we assumed that a bubble begins to slide as soon as it nucleates, an idea supported by the optical measurements in figure 4. To validate this approach, results from the numerical integration of the exact differential equation (3.13), analytical solution ((3.20), (3.21) and (3.24)) and experimentally measured bubble trajectories are shown in figures 10–12. Like the bubble growth predictions, variations in the bubble growth fitting parameter $C_{RB}$ can also explain most of the statistical spread in the bubble trajectory data. Unlike the bubble growth results, there is a clear dependence on the mass flux due to the derivations and analysis described above. Increasing pressure tends to slow down the bubble because of the decrease in $C_{RB}$ (e.g. see (3.23)). Finally, one can observe the excellent agreement between the differential equation solution and the simplified analytical solution, even though simplified drag coefficients were used, and buoyancy was neglected. A sensitivity study on the impact of drag and buoyancy was performed, and it was found, and shown in Appendix A, that the simplification of the drag coefficient had a larger, yet insignificant, impact on the model predictions. Figures 10–12 also show the magnitude of the forces acting on the bubble during sliding. At a fixed pressure, when $C_{RB}$ is effectively constant, increasing mass flux increases the magnitude of the drag force. This is of course expected, as viscous drag forces are approximately proportional to the velocity of the flow (noting that the slip ratio is approximately uniform), which will obviously increase with mass flux, because of both the larger mean fluid velocity and smaller viscous sublayer thickness. At a fixed mass flux, increasing pressure, i.e. smaller $C_{RB}$, tends to increase the role of drag forces, as suggested by (3.16).

Figure 8. Schematic of a force balance on a single, sliding vapour bubble.

Figure 9. Experimentally measured bubble velocities and distances from wall across all operating conditions, with theoretical predictions overlayed.

Figure 10. Distance bubble travels (a,c,e) and forces acting on bubble (b,d,f) while sliding at a mass flux of 500 kg m$^{-2}$ s$^{-1}$.

Figure 11. Distance bubble travels (a,c,e) and forces acting on bubble (b,d,f) while sliding at a mass flux of 1000 kg m$^{-2}$ s$^{-1}$.

Figure 12. Distance bubble travels (a,c,e) and forces acting on bubble (b,d,f) while sliding at the highest mass fluxes tested (2000 kg m$^{-2}$ s$^{-1}$ at 1 MPa and 1500 kg m$^{-2}$ s$^{-1}$ for 2–4 MPa).

Figure 13. Proposed departure criterion illustrated with departure data at 1 MPa, with a mass flux and heat flux of 1000 kg m$^{-2}$ s$^{-1}$ and 500 kW m$^{-2}$, respectively. Note that the scales for the graph and images are different.

In passing, we note that the bubble sliding phenomenon can be conceptualized and modelled as a small particle moving in a liquid medium, leading to simple and convenient equations of motion that are strongly supported by experimental measurements. The major conceptual difference between the two phenomena is the fact that the vapour bubble is growing with time. If the bubble remains at a fixed size, then the equation of motion is simply

(3.25)\begin{equation} \frac{{\rm d}U_{b}}{{\rm d}t}+\frac{9\nu_{l}}{2R_{b}^{2}}U_{b}=\frac{9\nu_{l}}{2R_{b}^{2}}U_{l}. \end{equation}

This first-order differential equation predicts that the particle velocity reaches the liquid velocity after ${\sim }1.1R_{b}^2/\nu _{l}$. In the case of the growing vapour bubble, there is always some non-zero slip that causes the bubble to never exactly reach the liquid velocity, unless the quantity $\varPi _1$ is very large, in which case the bubble velocity is essentially the liquid velocity for all practical purposes.

Figure 14. Departure diameter (a) and growth time (b) predictions as a function of pressure for various mass fluxes. Solid lines show predictions from (3.13), while the dashed lines show predictions from (3.29) and (3.30). For the sake of illustration, $C_{RB}=B/2$ and $q''=500$ kW m$^{-2}$.

3.3. Bubble departure

The observations of the HSV images and the analysis in the previous section support the hypothesis that bubble sliding begins as soon as the bubble nucleates. This interesting finding leaves us with an open question: How to define the bubble departure diameter and growth time? We may define the bubble departure diameter as the equivalent diameter that the bubble has when its optical footprint stops covering its original nucleation site, allowing for a new bubble to nucleate. This is defined to be consistent with the departure diameter and growth time definitions adopted in HFP models. This criterion is illustrated in figure 13, where the bubble is first covering the nucleation site and later moves past it after $\sim$0.2 ms with a diameter of 0.1 mm, which define the grow time and departure diameter, respectively. For a spherical bubble, the departure criterion can be represented mathematically as the time at which

(3.26)\begin{equation} x_{b}(t)=\int_{0}^{t}U_{b}(t)\,{\rm d}t=R_{b}(t). \end{equation}

Equation (3.26) can be readily solved analytically. In the case where the bubble departs from the nucleation site while inside the viscous sublayer (which will be shown to be the case the vast majority of the time), the growth time and departure diameter can be evaluated as

(3.27)\begin{equation} t_{g}=\frac{3}{2}\frac{\varPi_1+1}{\varPi_1}\frac{\nu_{l}}{U_\tau^2}, \end{equation}

and

(3.28)\begin{equation} D_{d}=C_{RB}\sqrt{\frac{6(\varPi_1+1)}{\varPi_1}\frac{\nu_{l}}{U_\tau^2}}, \end{equation}

respectively. In the limit where $\varPi _1\gg 1$, $U_{b}=U_{l}$

(3.29)\begin{gather} t_{g}=\frac{3}{2}\frac{\nu_{l}}{U_\tau^2}, \end{gather}

(3.30)\begin{gather} D_{d}=C_{RB}\sqrt{\frac{6\nu_{l}}{U_\tau^2}}, \end{gather}

which provide extremely simple expressions for the departure diameter and growth time. To demonstrate that these simple formulae can be sufficiently accurate, the departure diameter and growth time determined numerically using (3.26) and the complete force-balance model (i.e. (3.13)) are compared against (3.29) and (3.30). Results are shown in figure 14, where it is immediately clear that, for most pressures above 1 MPa, the simple analytical expressions for departure diameter and growth time compare extremely well with the more complex force-balance model. Note that (3.27) and (3.28) would lie in between the predictions of the force balance and (3.29) and (3.30) and are not shown in figure 14 for the sake of clarity.

It can also be seen that (3.29) and (3.30) consistently underpredict the results at lower pressures. Assuming the bubble remains in the viscous sublayer, this can be understood by considering the effect pressure has on the quantity $\varPi _1$ and, consequently, (3.27) and (3.28). Because decreasing pressure increases $C_{RB}$, $\varPi _1$ also decreases (see (3.17)). Since $\varPi _1$ is directly related to the vapour–liquid slip velocity ratio (3.20), decreasing $\varPi _1$ increases the relative slip between the two phases, making the assumption inherent to (3.29) and (3.30) (i.e. $U_{b}=U_{l}$) less accurate. It must also be confirmed that the bubble remains in the sublayer during sliding. Figure 14 shows that the bubble (see curve labelled $D_{d}^+$) remains in the viscous sublayer beyond a pressure of approximately 2 MPa, and is never far inside the turbulent core. Therefore, both assumptions (i.e. that the bubble remains in the viscous sublayer and that $U_{b}=U_{l}$) are more accurate at elevated pressures.

Based on (3.29) and (3.30), one can easily interpret how different operating conditions affect bubble departure. For example, it can be readily deduced that the departure diameter is approximately inversely proportional to the mass flux ($\propto G^{0.9}$) and directly proportional to $C_{RB}$, while the growth time is inversely proportional to $\propto G^{1.8}$ and does not at all depend on how fast the bubble grows. This observation can be deduced from the fact that the speed of the bubble is directly related to its size while in the viscous sublayer, which makes the time at which the bubble departs independent of its growth rate.

Interestingly, it can also be observed that all departure diameters shown in figure 14 lie below the Cole–Rohsenow predictions for pool boiling conditions (purple dash-dotted line figure 14). This is expected due to the presence of drag and added mass forces promoting bubble sliding. However, it must be noted that this finding should not be extrapolated to pressures lower than what was investigated in this study, as bubble deformation effects, due to their size, will create non-negligible surface tension forces (Favre et al. Reference Favre, Colin, Pujet and Mimouni2023) that make our force-balance model and analytical solution invalid.

Equation (3.26) can also be used to measure the departure diameter and growth time experimentally. This is achieved by our bubble tracking algorithm, where the intersection of the two curves in figure 13 is determined using linear interpolation. A comparison between average experimental measurements and theoretical predictions is shown in figure 15. The error bars in the departure diameter and growth time measurements come mostly from the spatial (${\sim }5\unicode{x2013}6\,\mathrm {\mu }{\rm m}$) and temporal ($33\unicode{x2013}67\,\mathrm {\mu }{\rm s}$) resolutions of the optical technique. In general, departure diameter predictions improve as pressure increases, as expected based on the model assumptions. It can also be seen that the growth time predictions are also reasonably accurate at higher pressures, but are not quite as precise as the departure diameter predictions. This happens because bubbles may not be precisely tracked the instant they are nucleated. Therefore, the experimentally measured growth time becomes less accurate the faster the bubbles depart. This is not the case for departure diameter, because the measurement is only taken when the detected bubble leaves the nucleation site area. In general, it is concluded that the model is in agreement with the optical measurements and their associated uncertainty. It is also clear that the measured departure diameters are indeed far smaller than what has been reported in past literature studies (see figure 1), as was expected. Even though there is no universal definition of departure diameter, our solution offers a physically sound view of the bubble departure process in high-pressure flow boiling compared with approaches previously used in the literature. Additionally, our analytical solution may offer much more tractable and computationally efficient ways to evaluate bubble departure and growth time compared with numerically solving complex equations of motion. In passing, we note that the assumption of departure diameter approaching zero as pressure increases has also been made in recent CFD flow boiling models (Kommajosyula Reference Kommajosyula2020), which is consistent with our mathematical models and experimental measurements.

Figure 15. Departure diameter (a) and growth time (b) predictions and measurements as a function of mass flux for various pressures.

Moreover, (3.29) and (3.30) can provide a simple and meaningful way to predict the effect of pressure and flow rate on the bubble departure process. As an example, we can use them to evaluate what would be the spatial and temporal resolutions needed to measure these quantities at pressurized water reactor conditions (${\sim }15\,{\rm MPa}$ and >2000 kg m$^{-2}$ s$^{-1}$). The departure diameter would be of the order of 100 nm and the growth time less than $10\,\mathrm {\mu }{\rm s}$. Such a small diameter cannot be resolved with the optical technique used in this study, which is limited to approximately $6\,\mathrm {\mu }{\rm m}$. It may be possible to achieve sub-$\mathrm {\mu }{\rm m}$ spatial resolution with optical microscopes, but such devices will have to be positioned a few millimetres away from the boiling surface, which makes the fabrication of a seal capable of withstanding $\sim$15 MPa increasingly difficult. In fact, flow boiling at a pressure of 7.6 MPa was attempted in the present study (shown in figure 16), but the uncertainty in the optical measurements was quite large. At a pressure of 7.6 MPa and mass flux of 1000 kg m$^{-2}$ s$^{-1}$, departing bubbles appear to be a few pixels or even less than one pixel. Based on the images shown in figure 16, the departure diameter appears to be approximately $2\unicode{x2013}3\,\mathrm {\mu }{\rm m}$. If we refer to figure 14, it can be seen that the model predicts a departure diameter of approximately $2\,\mathrm {\mu }{\rm m}$ at a pressure of 7.6 MPa and 1000 kg m$^{-2}$ s$^{-1}$ mass flux, which indicates that the present model and measurements are approximately in agreement, but will need to be validated with future measurements. Based on these findings and limitations, future studies that focus on bubble departure in high-pressure flow boiling applications will either have to use a more advanced optical set-up with higher spatial resolution or make use of light with shorter wavelengths (e.g. ultraviolet or X-rays, provided that these techniques have sufficient temporal resolution).

Figure 16. Bubble columns appearing from two nucleation sites in $10\,^{\circ }{\rm C}$ subcooled flow boiling at 7.6 MPa with a mass flux and heat flux of 1000 kg m$^{-2}$ s$^{-1}$ and 220 kW m$^{-2}$, respectively.