Numerical Study of Bubble Behavior under Gradient Flows during Subcooled Flow Boiling in Vertical Flow Channel

: In this study, we examined the condensing behavior of single and multiple bubbles of pure steam in a subcooled liquid phase using a fully compressible two-phase homogeneous mixture method that is solved by an implicit dual-time preconditioned technique. The interface between the liquid and vapor phases was determined by the advection equations using a compressive high-resolution interfacing capturing method. The spurious current reduced near the interface, a smoothing filter is applied to the progress curvature calculation. The sensitivity study carried out to predict the empirical constant by using Lee’s mass transfer model. A comparison of the numerical and experimental results highlighted that the proposed model accurately predicted the behavior of the definite condensing bubble. Furthermore, the single and multiple bubble condensation behaviors were investigated for different initial subcooled temperatures, and bubble diameters under various gradient flow, such as velocity gradient, temperature gradient, and velocity and temperature gradients. Subsequently, the effect of multiple bubbles flows in different bubble pattern forms, and their condensation was studied. The coalescence of bubbles depends on the subcooled temperature. Furthermore, the bubble diameter, the gap between the bubbles, and the flow rate of the bubbles were also observed.


Introduction
Heat transfer characteristics in bubble dynamics are essential phenomena in the design of nuclear reactors, boilers, and cooling system in electronic devices. Existence of vapor bubbles impacts the heat transfer phenomena, pressure drop, and flow instability. Heat transfer characteristics are highly pronounced under subcooled boiling flow conditions, and hence, the knowledge of bubble dynamics is crucial during the design of several industrial applications.
Nowadays, this research is focused on theoretical study, experimental analysis, and numerical simulation of the application, as mentioned above (bubble condensation in subcooled flow boiling). Theoretical studies were conducted [1] and [2] on stationary and dynamic bubble condensation. In the last two decades, the technological improvements in the measuring and monitoring instruments used in bubble experiments have caused an increase in the knowledge of bubble behavior. Significant research progress has been experimentally made to understand the bubble behavior in terms of bubble size history, shape, velocity, and collapse time. Furthermore, the interfacial heat transfer coefficient was investigated under high-pressure/temperature conditions [3][4][5] and low-pressure conditions [6][7][8]. It is difficult to obtain vital information on bubble behavior because their shape and interfacing area are too complicated for measurement. Furthermore, bubble growth and condensation occur in less than a microsecond; thus, they are challenging to investigate. Therefore, numerical simulation is regarded as the only method that can be employed to overcome such difficulties.
Condensation of a pure single steam bubble [9][10][11][12][13][14][15][16] and a mixture bubble [17,18] have been studied extensively. Several researchers computationally investigated the bubble nucleation, growth, departure, film-and flow boiling, and condensation [19][20][21]. The simulation performed [22] focused on bubble behaviors under low-and high-pressure, different bubble sizes, and different subcooled temperatures; and compared different velocities and bubble lives for a single bubble. Moreover, simulation of multi-bubble behavior (under different gaps between bubbles, different bubble sizes, and different velocities) and validation of the surface tension force using a square droplet were performed. However, only very few simulations have been conducted for the multi-bubble study [22][23][24][25]. Furthermore, multi-bubble behaviors under various gradient flow conditions have not been studied. Therefore, considering the perspectives mentioned above, the numerical simulation in this study focused on single-and multi-bubble condensation flow under velocity and temperature gradients. For simulating the gradient flow condensation, we used to apply different operating conditions, such as different bubble diameters, subcooled temperatures, velocities, and different multi-bubble pattern arrangements.
The in-house code for the condensation model is validated using a vapor bubble condensed in liquid [22]. A compressive interface-capturing scheme was applied in the simulation for interfacing the liquid and vapor phases. A detailed review of the advantages and disadvantages of this scheme was conducted [26][27][28][29]. The proposed model was further validated using experimental data and numerical results of cavitation. While numerically analyzing the flow with low Mach number, the preconditioning method improves the numerical stability and reduces the computational time [30][31][32][33][34]. For analysis the cavitation and flashing phenomena in the convergent divergent nozzle, the compressible homogeneous mixture model and a dual-time preconditioned approach is used for the calculation of the empirical constant range [35]. S Paramanantham et al. [22] extended the same computation method for bubble condensation using a vertical channel that was simulated. This prediction was successfully validated using the experimental data. Using the same method, Phan et al. [36] analyzed the air-steam mixture condensation flows in a vertical tube with different air-steam mixture flow conditions.

Geometrical Configuration
In this simulation, we considered the vertical rectangular channel geometrical models for the bubble rod condensation problem. The two-dimensional (2D) computational domain of the geometrical model was 10 × 20 mm, as shown in Figure 1a. In the vertical flow channel, the bubble rod was discharged from the bottom of the channel. We investigated the effect of bubble behavior and various phenomena under different flow conditions. We considered the applied velocity at the inlet, extrapolation pressure at the exit, and sides of the domain as a slip wall to be the boundary conditions in this simulation.

Governing Equations
In this numerical simulation, a fully compressible two-phase homogeneous mixture approach is used to model the direct contact condensation process of pure bubble condensation flows in subcooled water. For a two-dimensional axisymmetric and planer flow, the approach consists of a compressible Reynolds Averaged Navier-Stokes (RANS) system of equations. It is established on a mass conservation equation for the individual phase (water & steam), and a momentum and energy conservation equation for the mixture phase model consists of five conservative equations. A condensation mass transfer model is employed to simulate the process of mass and energy transfer across interfaces. For two-dimensional flows, dimensionless forms of these governing equations are written in generalized coordinates ( , ) as follows: where x and y denote the Cartesian coordinates; p denotes the pressure; u and v denote the velocity components in the x-and y-directions, respectively; T denotes the temperature; Y is denotes mass fraction; is the condensation rate; ξ and η denote the generalized coordinates; and are the vector of the primitive variables and Source term, respectively; and denote the flux terms; and denote the viscous terms; h, ht and hlv denote the mixture enthalpy, the total enthalpy and the latent heat of evaporation, respectively; g, ρ, µ, and α are the gravity, mixture density, mixture viscosity and void fraction, respectively; σ denotes surface tension coefficient and k denotes the surface curvature; and subscripts v denotes the steam and L denotes water, respectively. The factor denotes whether the governing equations are written for planar flow ( = 0) or axisymmetric flow ( = 1). τ is the pseudo-time and J denotes the Jacobian matrix: The mixture properties, i.e., density, viscosity, enthalpy, thermal conductivity, specific heat, and Prandtl number are given in S Paramanantham et al. [22]. The properties of water and steam were estimated by employing the IAPWS-97 (International Association for the Properties of Water and Steam) database [37].
The flux Jacobian, , and the preconditioning matrix, , are defined as: A study [22] provided the derivatives of mixture enthalpy, density, and other quantities in Equations (8) and (9).
The interfacial tension force was applied as a continuous surface force based on Brackbill et al. method [38]. The curvature is defined as: where is the unit normal vector to the interface: The governing equations are coupled to the system of interface advection equation for the volume fraction of the vapor phase to capture the interface between the phases (liquid-vapor). The interface advection equation for the volume fraction of steam is given in Ha et al. [39].
The simplified saturation model developed by Lee was used for the condensation and evaporation processes [40]. The mass transfer is considered to occur at constant pressure, and quasithermal-equilibrium state, which gives the relationship described in the Equation (12).
where * = / ∞ is the empirical scaling constant calculated using the mass transfer intensity factor.

Numerical Method
Implementation of the dual-time preconditioned system in the numerical method is established based on the references [26] and [31] . The discretization form of Equation (1) is given as follows: where I denote the identity matrix = ⁄ + 1.5 / ; denotes the first-order forward differencing; denotes the first-order backward differencing; ∆ , = , − , superscripts n+1 denotes the new physical-time level; the superscript m and m+1 denote the actual and the new pseudo-time-levels; Δt is the physical time-step; Δτ is the pseudo-time-step; and 1,m is the residual. For simplification, the superscript n+1, m is denoted as m. is the residual represented as follows: Time discretization is approximated by using the second-order backward difference method, physical-and pseudo-time derivatives are approximated with first-order backward differences. By using the MUSCL (Monotonic Upstream-Centered) upwind scheme, approximated for the spatial discretization and the inviscid terms. The second-order central differencing scheme applied to discretized for viscous and surface tension terms

Geometrical Model Validation
In this case, the liquid temperature (Tliquid) of the whole domain of the channel is defined from lower than the saturation temperature (Tliquid = Tsat − ΔT). To acquire the grid resolution using Figure  1a, the initial bubble diameter of 4 mm, and the subcooled temperature (ΔT = 10 K) and velocity (u = 0.3 m/s) are plotted in Figure 2. Four types of grid densities were considered: 50 × 100, 100 × 200, 150 × 300, and 200 ×400. Among the tested grid types, the 150 × 300 and 200 × 400 grids showed similarities in the bubble diameters and shapes, as shown in Figure 2a and Figure 2b, respectively, for the range of grid sizes showing satisfactory convergence. Thus, in this study, the bubbles were resolved using the 150 × 300 grid size. This grid arrangement was used to obtain the results presented in this study.

Mass Transfer Coefficient
The condensation constant (CRV) in the Lee phase change model [40] is essential for the prediction of a suitable value by comparing the bubble condensing behaviors while varying the CRV, as shown in Figure 4. A CRV value of 1750 is appropriate for this simulation, which is comparable to the experimental results reported by Kim and Park [6].   Figure 6 shows a comparison of the experimental results and the present numerical results for the bubble shapes. The bubble rising behavior is numerically calculated with a two-dimensional axisymmetric assumption. The bubble surface tension is calculated using the continuous surface force method. On using continuous surface force method, even surface force applied on the surface of the bubble. Therefore, even though the initial square droplet relaxed to the equilibrium circular shape after some time [22]. Due to these assumptions, a considerable difference between experimental and numerical results of a bubble shape.

Single-Bubble Rising Behavior
We performed an analysis of the single-bubble condensation behavior under different gradient conditions, such as velocity gradient, temperature gradient, and velocity and temperature gradients, as shown in Figure 1.

Velocity Gradient Flow
To study the effect of condensing bubble-moving trajectory under the velocity gradient flow, we used the computational 2D bubble released at the bottom of the vertical channel, and the flow direction was upward. A bubble is simulated in a subcooled liquid with a linear velocity gradient at the inlet, as shown in Figure 1b and test conditions are listed in Table 1. Figure 7a-d shows the simulated moving path of the bubble and its shape change under the velocity gradient with different scenarios. Figure 7a shows the bubble diameter, D of 2 mm and linear velocity gradient values, vmin and vmax of 0.1 m/s and 0.3 m/s, respectively. The moving bubble path change occurred in a significant amount, but it is negligible, and bubble condensed.   Moreover, increasing the subcooled temperature with the velocity gradient flow ensures that the moving bubble path is approximately on the bubble axis. The bubble moving path may exist in the liquid depending on the bubble shape and deformation of the bubble by interaction with the surrounding liquid. The influence of the surface tension is large enough in this small bubble size to restrict the changes in the shape. Thus, a small-sized bubble will rise vertically under the liquid velocity gradient.
Furthermore, noticeable changes in the bubble shape and moving path were observed for larger bubble sizes with diameters of 3 mm and 4 mm. Figure 7b and Figure 7c illustrates that the condensing bubble shape changes from spherical to ellipsoidal via a wing shape and finally retained its original shape. It was observed that in the velocity gradient flow, migrating bubbles did not move across the sides more than the bubble radius distance. If the subcooled temperature is increased, the bubble will move closer to a lower velocity region. As a result, bubble condensation accelerates the crosswise migrating bubble in the surrounding liquid. Decreasing the velocity gradient and increasing the subcooled temperature tend to restrict the migrating bubble. It moves approximately in a vertical line for a bubble diameter, D of 4 mm and linear velocity gradient values, vmin and vmax of 0.05 m/s and 0.15 m/s, respectively, as shown in Figure 7d. Commonly, the bubble shape is restrained by the buoyant force, surface tension, and the surrounding liquid's inertial force, and the condensation rate. The existence of the velocity gradient modifies the inertial force, which successively influences the bubble shape, condensation rate, and bubble trajectory. Figure 8 shows the comparison of bubble condensing life under with velocity gradient and uniform velocity flow in ΔT = 5 K. the condensation rate majorly affected in the velocity gradient compared with uniform flow. In the velocity gradient flow, inertial force is increasing the bubble life, i.e., the rate of condensation is reduced.

Temperature Gradient Flow
The effect of the local condensation rate and the bubble surface is simulated, as shown in Figure  1c and test condition are listed in Table 2, to examine the influence of the linear temperature gradient on liquid flow in the bubble condensation channel.  Figure 9a shows the predicted sizes of the bubble diameter. The bubble rate of condensation on the right side of the bubble interface area is more significant than the left side because of the local liquid temperature difference, which results in an asymmetrical shape of the migrating bubble. Nevertheless, this study shows that the effect of the local condensation rate does not affect bubble deformation because a negligible amount of deformation occurs. As shown in Figure 9b the initial normal velocity does not affect bubble deformation in the temperature gradient part of the simulation.

Velocity and Temperature Gradient Flow
The simultaneous effect of the temperature and velocity gradients is examined in this section. Figure 10 shows the condensation and deformation of the bubble in liquid along with the temperature and velocity gradients, as shown in Figure 1d and test conditions are lists in Table 3. The migrating bubble generated a sinuous path for the larger bubble and a negligible path for the smaller bubble. Although the velocity gradient dominates in the flow, the bubble moves to a higher temperature region due to the high condensation rate in the temperature gradient. The velocity and temperature gradients show that the lateral migrating bubble can significantly affect the bubble lifetime.    Compared to a single bubble condensation, the multi-bubble simulation is more complicated because it contains numerous physical phenomena that occurred, such as bubble coalescence, and it makes rapid and continuous shape changes. For examined these behaviors, and along with well as the effect of gradient under different bubble patterns such as two similar and dissimilar bubbles in a row with and without gradient flow, three similar bubbles in a row with and without gradient flow, and three similar bubbles in a triangular shape and delta shape with and without gradient flow, as shown in Figure 3. The following analyses were performed to investigate the multi-bubble rising behavior with mass transfer simulation: (1) Two-bubble rising behavior, (2) Three-bubble rising behavior.

Two-Bubbles Rising Behavior
In this simulation, two bubbles were placed at the centerline of the domain. Subsequently, the gap between the two bubbles was varied, and the bubble merge process was examined based on the bubble gaps (BGs). Figure 11 shows similar bubbles with different gaps between each bubble. When the gap is above 0.2 mm, the merging process will be affected, and condensation occurs separately. We observed that during the bubble merging and condensation processes in the liquid, the initial shape of the bubbles changed rapidly and continuously with time. For the case of no gap between the bubbles, the merging process started, and the bubble shape-changing occurred with very less time steps. Moreover, for the case with a BG of 0.2 mm, the bubbles merged after some time steps.
The bubble merging process break due to the increase in the certain gap between bubbles. When bubbles move vertically along the rising line, bubble coalescence did not occur because the lowest bubble condensed faster than the upper bubble. The upper condensing bubble releases the condensation heat towards the lower bubble, which increases the liquid temperature near the lower bubble. It will reduce the lower bubble condensation rate, but the random perturbation produced by the upper bubble can raise the condensation heat transfer of the lower bubble [22,23]. Local perturbation causes the lower bubble to condense faster than the upper bubble. Therefore, the condensation rate of the lower bubble is influence by liquid temperature and random perturbation. In the case of two dissimilar bubbles in a row without the gradient flow, the same phenomena are observed as in similar bubbles case.    Figure 12 shows the dissimilar bubble condensation at different bubble gap without gradient flow condition. Bubble condition. Increasing the bubble gap small enough, it breaks the bubble coalescence. Even though a small gap occurred between bubbles, the coalescence not occurred due to the lower bubble large than the upper bubble. However, the local perturbation causes the lower bubble to condense faster than the upper bubble. Figure 13 shows a comparison of the bubble equivalent diameter for bubble condensation behaviors with different BGs. Compared to the merged bubble process, the case without bubble merging may increase the condensation process because local saturation can accelerate condensation.   Figure 14 show the different subcooled temperature with BG = 0.2 mm for D = 3 mm. During the bubble merging, the bubble shape-changing occurred rapidly with different shapes that were observed. The merging process depends on the subcooling temperature. When the subcooling temperature increases, the bubble merging process will also be affected, and under the velocity gradient flow, the migrating bubble generated a negligible range of sinuous path because of the bubbles interlink and tends to move toward the vertical direction approximately.  Figure 15 shows that gradient velocity flow at different subcooling temperatures with a BG between 0.2 mm. When the gap between the two bubbles is small enough, it merges. The lower bubble moves faster than the upper bubble because of the buoyancy force and the two-bubble capillary reaction that occurs between them. The upper bubble flow creates a lower velocity region between the bubbles, and then the upper bubble pulls the lower bubble to merge, as shown in Figure 15a. Figure 15b shows that similar phenomena occurred because of the lower subcooled temperature (ΔT = 20 K) that those occurred with ΔT = 10 K. The upper bubble tends to pull the lower bubble toward itself; however, the lower bubble loses its mass before that. Thus, the capillary effect is broken, and the bubble merging process will not occur. Figure 16 shows the bubble under the temperature gradient. We applied different velocities and observed different characteristics. Decreasing the bubble velocity affects the bubble merging process and increases condensation in the temperature gradient flow. Figure 17 shows that under the velocity and temperature gradient flow, the complete merging process is restricted, even in the low-velocity gradient range. However, when compared to the case without the gradient flow, bubble condensation is affected by the gradient flow condition.  Figure 18 shows the simulation of three similar bubbles at an equal distance in a liquid to examine the condensation effects. Figure 18a illustrates the simulation of the bubble condensation process in liquid without the gradient flow. The bubble shapes were not affected, and the local perturbation improved the bubble condensation heat transfer. A considerable amount of bubble deformation occurred in the three bubbles; the lower bubble condensed, followed by the upper bubble, and these two bubbles slightly restricted the condensation of the middle bubble. Figure 18b show three bubbles rise in a row condensation under the gradient flow; the bubble merge from the lower bubble to the upper bubble occurred. The bottom bubble merged completely after coalescing with the upper bubble. The coalescing bubble shape bubble formed was similar to a bowling pin, and condensation occurred from the bottom of the bubble with continuously changing shape, as reported [22]. Figure 19 shows a comparison of bubble condensation, with gradient and without the gradient. It can be observed that the case without the gradient flow improved bubble condensation as compared to that with the gradient flow. Figure 19. Comparison of three-bubbles condensation in a row with and without the gradient.

Three-Bubbles Rising Behavior
Finally, we performed a simulation of more complicated multiple bubble patterns with triangular and delta arrangements to study the bubble behavior under velocity and temperature gradient with different gaps between the bubbles, as shown in Figure 20. Complicated shapes were formed during the bubble merging process, and the bubble shape transformed rapidly and continuously into an indistinguishable shape. The triangular and delta patterns showed a similar type of bubble merging shape characteristics because the initial pattern formed a small liquid drop in the merged indistinguishable shaped bubble. When the BG is sufficiently large, the merging process is restricted, and bubble condensation occurs separately. Initially, the lower temperature gradient region condensed faster than the higher region. However, the velocity gradient moves from the higher side to the lower side, i.e., from the right side to the left side, as shown in Figure 20.

Conclusions
In this study, we numerically investigated the condensation behaviors of single and multiple bubbles under different gradient conditions using a homogeneous mixture model. The solver is developed for solving the two-dimensional, and compressible flow. The mass transfer was modeled based on the Lee phase change model. At subcooled temperature conditions, the numerical results of bubble condensation characteristics were compared with two different experimental data. The numerical solutions are in fairly agreement with the experimental results. The condensation bubble behavior under the gradient flow was studied qualitatively under several conditions. The bubble size, BG, subcooled bubble temperature, and bubble pattern had significant effects on the bubble lifetime, bubble deformation, and bubble shapes. The following conclusions were drawn with the parameters for bubble condensation: 1. When the bubble moved in the velocity gradient, it migrated toward the lower velocity side, thereby decreasing the bubble diameter and reducing the sinuous path. Similarly, in the temperature gradient, the bubble moved toward the lower temperature side and increased the condensation on one side of the bubble. In the case of velocity and temperature gradient flow, the bubble initially moved toward the higher velocity side; however, it again transformed and moved to the lower velocity and higher temperature side. 2. A similar shape and size bubble merger occurred when the gap in the bubble was sufficiently small (BG = 0.2 mm). This occurred up to a certain subcooled temperature; however, it was interrupted beyond the subcooled temperature ΔT = 15 K. 3. Compared to the three bubbles in a row condensation, the case without gradient flow condensation, the case with velocity and temperature gradient flow condensation affects the life of the bubble's condensation, and this case causes the bubble merging process. 4. The bubble pattern affects the bubble shapes and condensation rates during the bubble condensation because of the attraction with the neighboring bubble. From these cases, the bubbles heat and mass transfer depend on the number of bubbles within a particular BG and subcooled temperature in the defined shape and random arrangements of bubble flow inside the bubble column. 5. When bubbles move vertically along the rising line, bubble coalescence did not occur because the lowest bubble condensed faster than the upper bubble due to the local perturbation were observed.