Influence of density and viscosity on deformation, breakage, and coalescence of bubbles in turbulence

We investigate the effect of density and viscosity differences on a swarm of large and deformable bubbles dispersed in a turbulent channel flow. For a given shear Reynolds number, ${\text{Re}}_{\tau }=300$, and a constant bubble volume fraction, $\mathrm{\Phi }\simeq 5.4\%$, we perform a campaign of direct numerical simulations of turbulence coupled with a phase-field method accounting for interfacial phenomena. For each simulation, we vary the Weber number (We, ratio of inertial to surface tension forces), the density ratio (${\rho }_r$, ratio of bubble density to carrier flow density) and the viscosity ratio (${\eta }_r$, ratio of bubble viscosity to carrier flow viscosity). Specifically, we consider two Weber numbers, $\text{We}=1.50$ and $\text{We}=3.00$, four density ratios, from ${\rho }_r=1$ down to ${\rho }_r=0.001,$ and five viscosity ratios, from ${\eta }_r=0.01$ up to ${\eta }_r=100$. Our results show that density differences have a negligible effect on breakage and coalescence phenomena, while a much stronger effect is observed when changing the viscosity of the two phases. Increasing the bubble viscosity with respect to the carrier fluid viscosity damps turbulence fluctuations, makes the bubble more rigid, and strongly prevents large deformations, thus reducing the number of breakage events. Local deformations of the interface, on the contrary, depend on both density and viscosity ratios: as the bubble density is increased, a larger number of small-scale deformations, small dimples and bumps, appear on the interface of the bubble. The opposite effect is observed for increasing bubble viscosities: the interface of the bubbles become smoother. We report that these effects are mostly visible for larger Weber numbers, where surface forces are weaker. Finally, we characterize the flow inside the bubbles; as the bubble density is increased, we observe, as expected, an increase in the turbulent kinetic energy (TKE) inside the bubble, while as the bubble viscosity is increased, we observe a mild reduction of the TKE inside the bubble and a strong suppression of turbulence.

Figure 1

Top view of four statistically stationary configurations ($t^{+}=4000$) for different combinations of density ratios (${\rho }_r=0.001$ and 1) and viscosity ratios (${\eta }_r=0.01,1,$ and 100). Panel (a) refers to $\text{We}=1.5$, while panel (b) refers to $\text{We}=3.0$. The subpanels are arranged in a plot using ${\rho }_r$ as $x$ coordinate and ${\eta }_r$ as $y$ coordinate. The effect of density can be appreciated in the sequence of panels on the middle row, while that of viscosity in the right column. The background of the plot shows the turbulent kinetic energy, $\mathrm{TKE}=(\rho /{\rho }_c)(u^{\prime 2}+v^{\prime 2}+w^{\prime 2})/2$ (white—low; black—high), computed on the central $x^{+}\text{−}{y}^{+}$ plane of the channel.

Figure 2

Time evolution of the number of bubbles, $N\left(t^{+}\right)$, normalized by its initial value $N_0$. Left column refers to $\text{We}=1.5$, while the right column to $\text{We}=3.0$. Top row: effect of density ratio, for ${\rho }_r=0.001,0.01,0.1,$ and 1 (with ${\eta }_r=1$); middle row: effect of viscosity ratio, for ${\eta }_r=0.01,0.1,1,10,$ and 100 (with ${\rho }_r=1$); bottom row: combined effect of density and viscosity, for the case with ${\rho }_r=0.1$, the cases ${\rho }_r=0.1, {\eta }_r=1$ and ${\rho }_r=1, {\eta }_r=0.1$ are reported for reference. On each line the left plot also includes the color code and a sketch with the definition of the property ratio considered (${\rho }_r, {\eta }_r$, or both ratios).

Figure 3

Time evolution of the normalized breakage rate, ${\dot{N}}_b\left(t^{+}\right)/N\left(t^{+}\right)$, and coalescence rate, ${\dot{N}}_c\left(t^{+}\right)/N\left(t^{+}\right)$. Left column refers to $\text{We}=1.5$, while right column to $\text{We}=3.0$. Top row: effect of density ratio, for ${\rho }_r=0.001,0.01,0.1,$ and 1 (with ${\eta }_r=1$); Middle row: effect of viscosity ratio, for ${\eta }_r=0.01,0.1,1,10,$ and 100 (with ${\rho }_r=1$); Bottom row: combined effect of density and viscosity ratios, for the case with ${\rho }_r=0.1, {\eta }_r=0.1$. Cases ${\rho }_r=0.1, {\eta }_r=1$ and ${\rho }_r=1, {\eta }_r=0.1$ are reported for reference. For each row of plots, the left plot also shows the color code and a sketch with the definition of the ratio considered (${\rho }_r, {\eta }_r$, or both).

Figure 4

Time evolution of the total interface area $A\left(t^{+}\right)$, normalized by its initial value $A_0$. Top row: effect of density, for ${\rho }_r=0.001,0.01,0.1,$ and 1 (with ${\eta }_r=1$); middle row: effect of viscosity, for ${\eta }_r=0.01,0.1,1,10,$ and 100 (with ${\rho }_r=1$); bottom row: combined effect of density and viscosity, for the case with ${\rho }_r=0.1, {\eta }_r=0.1$. Cases with ${\rho }_r=0.1, {\eta }_r=1$ and ${\rho }_r=1, {\eta }_r=0.1$ are reported for reference. These effects are shown for two different Weber numbers: (a, c, e) $\text{We}=1.5$ and (b, d, f) $\text{We}=3.0$. On each row the left plot also includes the color code and a sketch with the definition of the ratio considered (${\rho }_r, {\eta }_r$, or both ratios).

Figure 5

Top view of the mean curvature of the bubble surface, ${\mathcal{K}}^{+}$, for four different combinations of density ratios (${\rho }_r=0.001$ and 1) and viscosity ratios (${\eta }_r=0.01,1,$ and 100) once a statistically stationary configuration is reached ($t^{+}=4000$). Panel (a) refers to $\text{We}=1.5$ while panel (b) to $\text{We}=3.0$. The subpanels are arranged in a plot using ${\rho }_r$ as $x$ coordinate and ${\eta }_r$ as $y$ coordinate. The effect of density can be appreciated in the sequence of panels on the middle row, while that of viscosity in the right column. Bubble surface (isolevel $\varphi =0$) is colored according to the local value of the mean curvature (low—blue; high—red).

Figure 6

Probability density function of the mean curvature, ${\mathcal{K}}^{+}$. Left column refers to $\text{We}=1.5$, while right column to $\text{We}=3.0$. Effect of density ratio can be appreciated on the top row for ${\rho }_r=0.001,0.01,0.1,$ and 1 (with ${\eta }_r=1$). The effect of bubble viscosity can be observed in the middle row for ${\eta }_r=0.01,0.1,1,10,$ and 100 (with ${\rho }_r=1$). Finally, the combined effect of the density and viscosity ratio is shown on the bottom row for the case with ${\rho }_r=0.1, {\eta }_r=0.1$, with respect to the cases where a single effect is considered (with ${\rho }_r=0.1, {\eta }_r=1$ and ${\rho }_r=1, {\eta }_r=0.1$).

Figure 7

Wall-normal behavior of the streamwise mean velocity profiles. Left column refers to $\text{We}=1.5$, while the right column to $\text{We}=3.0$. Density ratios effects are shown on the top row for ${\rho }_r=0.001,0.01,0.1,1$. Viscosity ratio effects are shown on the middle row for ${\eta }_r=0.01,0.1,1,10,100$. Finally, the combined effect of the density and viscosity ratios is shown on the bottom row for the case ${\rho }_r=0.1$ and ${\eta }_r=0.1$, with respect to the cases where only one effect is considered. As a reference, the classical law of the wall, $u^{+}=z^{+}$ and $u^{+}=(1/k)log{z}^{+}+5$ (with $k=0.41$ the von Kármán constant) is also reported with a dashed line. For all cases, with respect to single-phase, we observe a minor increase of the mean velocity.

Figure 8

Scatter plot of the wall-normal location of each bubble over its size for the different cases considered. The two black dashed lines identify the condition for which the interface of the bubble intercepts the closer wall in the hypothesis of a perfectly spherical bubble. Smaller bubbles tend to disperse along the entire channel height can get rather close to one of the two walls while larger bubbles tend to accumulate at the center of the channel.

Figure 9

Contour map of the turbulent kinetic energy in a $x^{+}\text{−}{y}^{+}$ plane located at the channel center ($z^{+}=0$). Panel (a) refers to the case ${\rho }_r=0.01$ and ${\eta }_r=1$ while panel (b) refers to the case ${\rho }_r=1$ and ${\eta }_r=0.01$. Both panels refer to the lower surface tension case ($\text{We}=3.0$) and to the time instant $t^{+}=4000$ (statistically steady configuration). The interface of the bubbles is highlighted with a white line. For ${\rho }_r=0.01$, bubbles and characterized by a low and uniform value of the TKE while, for ${\eta }_r=0.01$, the TKE map is nonuniform and characterized by small scales fluctuations.

Figure 10

Mean value of the turbulent kinetic energy (TKE) inside the bubbles. In panel a, TKE is evaluated using the complete definition of specific TKE (i.e., including the prefactor $\rho /{\rho }_c$) while, in panel b, TKE is evaluated considering only the velocity contribution (i.e., not considering the prefactor $\rho /{\rho }_c$). For both panels, a dashed line ($\text{We}=1.5$) and a continuous line ($\text{We}=3.0$) are used to show the behavior of TKE as the density or viscosity ratios are changed. Each value of TKE is marked with a circle (empty for $\text{We}=1.5$ and filled for $\text{We}=3.0$), with a red-color scale for the nonmatched density cases and a blue-color scale for the nonmatched viscosity cases, while the black color is used for the reference case.

Figure 11

Possible cases considered for the algorithm: panel (a) corresponds to a translation, panel (b) to a breakage, and panel (c) to a coalescence. Red bubbles are at the current time step ($n$), while blue bubbles are at the next time step available ($n+1$). Semitransparent bubbles show the estimated position, ${\mathbf{x}}_{i,\text{est}}^{n+1}$. Arrows show the trajectory of the bubbles, $\mathrm{\Delta }T{\mathbf{u}}_i^n$.

Figure 12

Flow chart of the algorithm used to detect breakage and coalescence events in the post-processing of the simulations.

Figure 13

Coalescence and breakage rates obtained using three different grid resolutions: $N_x\times N_y\times N_z=256\times 128\times 257$ (triangles), $N_x\times N_y\times N_z=512\times 256\times 513$ (circles) and $N_x\times N_y\times N_z=1024\times 512\times 1025$ (squares). The results refer to the case $\text{We}=3.0, {\rho }_r=1.000,$ and ${\eta }_r=1.00$.