Due to the complex phase change and heat transfer processes, the mechanisms of cavitation bubble collapse near a rigid boundary are well recognized to be complicated. Based on a modified large-density ratio multi-relaxation-time pseudo-potential lattice Boltzmann model, a single and a dual cavitation bubble collapse process near a rigid boundary with large-density and various viscosity ratios are simulated in the present study. Effects of density ratio, viscosity ratio, initial pressure difference, and distance between the cavitation bubble and wall on the cavitation process are studied. Furthermore, the evolution of maximum pressure, micro-jet velocity, lifetime, deformation index, and the first introduced total kinetic energy of cavitation bubbles are analyzed in the development of cavitation. Simulations show that the interaction mode of the bubbles and the distance between the rigid boundary and the lower bubble are key factors in determining the effect of aeration reduction. The study also shows that the proposed lattice Boltzmann pseudo-potential model is a robust and effective tool for studying the collapse of near-wall cavitation bubbles and has potential to predict the interaction of cavitation bubbles in the presence of complex boundaries.

Near-wall cavitation is a widespread phenomenon in both natural and engineering processes, such as hydraulic engineering, navigation, and humoral circulation. In growth and collapse process, there is a rapid phase transformation between gas and liquid phases, accompanied by a large increase in temperature and the generation of micro-jet and shock wave, which results in noise, vibration, and cavitation erosion. Although a substantial quantity of research has been conducted in the last few decades,1–3 some mechanisms are still not clear, including the evolution of the temperature field,4 deformation of cavitation bubbles,5 and threshold parameters of bubble formation and collapse near the complex boundary.6 However, cavitation in the near-wall region is more complicated and the cavitation bubble interface changes very sharply. It is difficult to theoretically study cavitation bubble collapse in the near-wall region.

To study the evolution behavior of bubbles while interacting with other bubbles or complex boundaries, a series of experiments were conducted. How the cavitation bubble deforms and the causes of formation of micro-jet and shock wave were studied.7,8 Blake9 studied the interaction mode of two cavitation bubbles near a rigid boundary. Later, bubbles collapsing and interacting with each other under complex boundary conditions were also studied experimentally.10,11 Through extensive observations, Philipp12 concluded that the mechanism of cavitation damage is a hydro-combined action of micro-jet, shock wave, and dramatic temperature change. All the experimental studies mainly focused on how the micro-jet, shock wave, and high temperature form, and how the bubbles deform and interact with each other.13–16 

With the development of computer technology, computational fluid dynamics (CFD) has become an important tool to study the evolution process of cavitation bubbles. When compared to experiments, the numerical simulation could test the evolution of cavitation bubbles under various condition combinations or even extreme ones, which are difficult to be achieved in experiments with little expense. In addition, the development of all physical fields (e.g., pressure, velocity, and density fields) could be checked in higher temporal and spatial resolution. Traditional CFD methods to study the evolution of cavitation bubbles mainly depend on the Navier–Stokes equation and Rayleigh–Plesset equation. Depending on different assumptions regarding the compressibility of cavitation bubbles and different interface capture methods, numerical simulation methods of macroscopic cavitation bubbles can be divided into three categories to study the causes of micro-jet, shock wave, high temperature, and cavitation noise.17–19 The first method assumes that the gas in the bubble and the surrounding liquid has the same fluid component. The gas–liquid phase change during bubble growth and collapse is then based on the equation of state (EOS) relating the density, pressure, and temperature of this fluid, and the compressibility of the fluid is also taken into account.20,21 According to the research of Ghahramani,21 a small time step and fine grids are needed to capture the liquid–gas interface. The second method is based on the Navier–Stokes equations together with a mass transfer equation to obtain the mass exchange between the gas and liquid phases. To improve its accuracy, some empirical constants are usually introduced into the Rayleigh–Plesset equation of this kind of model to obtain the gas–liquid mass transfer rate.22,23 The third method, also called the discrete bubble model (DBM), uses the Euler equation to simulate the continuum flow and a Lagrange method to capture the changes in shape and the motion of cavitation bubbles.19,24 This method also needs a fine grid, and can be used to simulate complex gas–liquid interface changes, such as the movement and collapse of a group of cavitation bubbles.

As a mesoscopic method, the lattice Boltzmann method (LBM), which is based on Boltzmann kinetic molecular dynamics, has become a robust and efficient tool to simulate the multiphase phenomenon. All the established LBM multiphase models can be categorized into color gradient model,25 pseudo-potential model,26,27 free energy model,28,29 interface tracking model,30 and entropy model.31 The pseudo-potential model, proposed by Shan and Chen,26,27 uses the interaction force between particles to form the interface automatically, and this special treatment makes interface tracking easier. With a non-ideal gas EOS introduced in this model, the pseudo-potential model is capable of simulating the multiphase phenomenon with phase change processes, such as liquid boiling32 and crystallization.33 After the first application of the pseudo-potential model in studying the growth of cavitation bubbles by Sukop,34 LBM gradually occupies a competitive place in the field of the numerical study of cavitation. When the Bhatnagar–Gross–Krook (BGK) collision operator35–38 is used to simulate the growth and collapse of cavitation bubbles, the calculated density ratio is 30–70, while an multi-relaxation-time (MRT) collision operator39,40 is used, the calculated density ratio can be increased to 720, range of density ratio values found in actual cavitation phenomena. In addition, the relevant research on cavitation bubble collapse is summarized in Table I.

TABLE I.

Parameters of the relevant research on cavitation bubble collapse.

ResearchersCollision operatorForce schemeT/Tcρl/ρgνg/νl
Shan37  BGK EDM 0.689 40 
Mao38  BGK EDM 0.68 40 
Shan39  MRT Li45  0.5 720 
Xue40  MRT Li45  0.8 12.8 
ResearchersCollision operatorForce schemeT/Tcρl/ρgνg/νl
Shan37  BGK EDM 0.689 40 
Mao38  BGK EDM 0.68 40 
Shan39  MRT Li45  0.5 720 
Xue40  MRT Li45  0.8 12.8 

Experiments and numerical results showed that the shock wave and micro-jet generated in the cavitation process is related to the gas–liquid viscosity ratio.41,42 However, because of the poor numerical stability of the LBM pseudo-potential model with a BGK collision operator with a large liquid/gas viscosity ratio, the effect of the viscosity ratio is not involved in recent numerical simulations of the cavitation process with LBM,35–40 and the viscosity ratio of gas–liquid was often set as 1 in these research studies. Considering the LBM multiphase flow model combining with the heat transfer equation can simulate both the flow field and temperature field,43 LBM seems to have more advantages to study the mechanism of the cavitation process. Based on a modified MRT pseudo-potential model, the collapse process of cavitation bubbles in the near-wall region is studied in the present study. The liquid–gas density ratio is about 750, and the gas–liquid viscosity ratio can be adjusted independently. The influence of different parameters, e.g., the viscosity of gas and liquid, the initial pressure difference, and the distance between the bubble and the wall, on the cavitation process is discussed. The maximum micro-jet velocity, the pressure field, and the lifetime of the bubble are studied. To explain the formation of the shock wave and micro-jet during the collapse of the cavitation bubble, the total kinetic energy (TKE) of cavitation bubble is first introduced in the present study, and the evolution of TKE during the collapse process is studied. Later, the interaction between two equal-size bubbles in the near-wall region is studied. This paper is organized as follows: a modified multi-relaxation time pseudo-potential model, which can simulate a large liquid–gas density ratio and different gas–liquid viscosities, is described in Sec. II. The evolution process of a single cavitation bubble and the interaction of two cavitation bubbles in the near-wall region are studied, and the results are shown in Sec. III. Finally, our conclusions are given in Sec. IV.

Compared to the BGK collision operator, the pseudo-potential model with an M RT collision operator has greater numerical stability and lower spurious currents when simulating some high-density ratio multiphase phenomenon.44,45 The particle distribution functions with external force terms can be expressed as 46 

fix+eiΔt,t+Δt=fix,tjΛ¯ijfjfjeq|x,tjIij12Λ¯ijSjΔt|x,t,
(1)

where fi is the particle distribution function, fjeq is the equilibrium distribution, x is the spatial position, ei is the discrete velocity of ith direction, Δt is the time step, Iij is the unit matrix, Λ¯ = M−1ΛM is the collision matrix, and the diagonal matrix Λ can be described as

Λ=diag(τρ1,τe1,τζ1,τj1,τq1,τj1,τq1,τν1,τν1),
(2)

where the relaxation time τν is correlated with the kinematic viscosity ν=1/cs2(τν0.5)Δt, cs is the lattice sound speed, M is the transformation matrix, and M−1 is its inverse matrix. Using M, the equilibrium moments meq can be obtained by projecting the equilibrium distribution feq onto the moment space, and for the D2Q9 LBM model, meq is expressed as44 

meq=Mfeq=ρ1,2+3ux2+uy2,13ux2+uy2,ux,ux,uy,uy,uxuy,uxuy.
(3)

Li45,46 proposed a modified external forcing scheme S, which can achieve thermodynamic consistency, given by the following expression:

S=ρ0,6uxFx+uyFy+12εFm2ψ2Δt(τe0.5),6uxFx+uyFy12εFm2ψ2Δt(τζ0.5),Fx,Fx,Fy,Fy,2(uxFxuyFy),uxFx+uyFy,
(4)

where the value of the parameter ε can be adjusted to achieve thermodynamic consistency. This modified forcing term enlarges the range of density ratio and increases the numerical stability of the model. ρ is the macroscope density and u is the macroscope velocity, which can be obtained as

ρ=ifi,ρu=ifiei+12FΔt,
(5)

where F = Fm + G + ⋯ is the total force acting on the fluid particles, Fm is the interaction force between particles, which can be obtained as44 

Fm=Gψxiωαψx+eiΔt,
(6)

where G is the interaction strength between two particles, ψ is the inter-particle potential, and ωαs are the weights. For the D2Q9 LBM model,47,ω1−4 = 1/3 and ω5−8 = 1/12. ψ can be calculated with a non-ideal gas EOS, which is introduced as48 

ψρ=2(pEOSρcs2)Gc2,
(7)

where pEOS is the pressure calculated by the non-ideal gas EOS.

Finally, the collision process can be described as

mRHS=mΛ¯mmeq+ΔtIΛ¯2S,
(8)

and the stream process can be expressed as

fix+eiΔt,t+Δt=fi*x,t=M1mRHS.
(9)

The computational domain is a 401 × 401 lu2 square domain as shown in Fig. 1. A Zou–He pressure inlet condition49 is applied at the upper boundary, non-equilibrium extrapolation conditions are applied at both the left and right boundaries, and the no-slip boundary conditions46 in Eq. (11) are applied at the bottom boundary. A spherical bubble with radius Rini = 50 lu is initialized in the computational domain, d is the distance between the bubble center and rigid boundary, pv is pressure inside the bubble, and p is the ambient and the inlet boundary pressures,

f2=f4,f5=f70.5f1f30.25Δt(Fx+Fy),f6=f8+0.5f1f3+0.25Δt(FxFy).
(10)
FIG. 1.

Schematic of the computational domain for the evolution of a single cavitation bubble.

FIG. 1.

Schematic of the computational domain for the evolution of a single cavitation bubble.

Close modal

The Carnahan–Starling (C–S) EOS is applied in the present study,48 

pEOS=ρRT1+bρ/4+(bρ/4)2(bρ/4)3(1bρ/4)2aρ2,
(11)

where a=0.4693R2Tc2/pc and b = 0.1873RTc/pc are the two parameters of the C-S EOS. The parameter a can affect the interface thickness w; considering both numerical stability and validity of the results, the initial interface thickness is suggested to be chosen as w = 4–5 lattice units.50 Therefore, the parameters in C-S EOS can be chosen as follows: a = 0.25, b = 4, and R = 1.

The density field is initialized as follows:

ρx,y=ρl+ρg2+ρlρg2× tanh2(xx0)2+(yy0)2Riniw,
(12)

where (x0, y0) is the center of the bubble and the initial interface width is set as 5 lu. In the present study, the interaction strength is set as G = −1, and the relaxation parameters are chosen as follows: τρ1=τj1=1.0,τe1=τζ1=1.1,τq1=1.1. Considering the density ratio of the actual cavitation phenomena in water as about ρl/ρg = 750, we choose the temperature T = 0.5Tc with the corresponding liquid–gas density as ρl/ρg ≈ 720. The viscosity relaxation time τν can be obtained by the following linear interpolation:51 

τν=τg+ρρgρlρg(τlτg),
(13)

where τg is the relaxation time of gas, τl is the relaxation time of liquid, ρg and ρl are the gas density and liquid density, respectively. Equation (13) provides a tunable viscosity ratio method for the present study. ε is chosen as 0.112 to maintain the numerical stability and thermodynamic consistency.

Some normalized parameters are introduced in this paper to describe the effect of the rigid boundary, including the deformation parameter Ra = h/l and the distance ratio of radius to distance between the center of the bubble and wall γ = d/Rini. In addition, the dimensionless time step t′ = t/tmax and dimensionless maximum pressure p′ = pmaxp are also introduced in our study. The units of measurement discussed in this paper include mass unit (mu), time step (tu), velocity unit (lu · tu−1), and pressure unit (mu · lu−1tu−2).

The radius of evolution of cavitation bubble in the near-wall region can be predicted by the Rayleigh–Plesset equation,36 

lndRBṘB2+RBR̈BṘB222νlRBṘB+σρlRB=pvpρl,
(14)

where RB is the bubble radius, νl is the liquid viscosity, and σ is the surface tension. The evolution of a bubble with d = 1.6Rini is simulated, at an initial pressure difference of 0.002 53 mu · lu−1tu−2, with the corresponding liquid viscosity νl = 0.0125 and surface σ = 0.0105. The comparison between the normalized result of LBM and Rayleigh–Plesset prediction is shown in Fig. 2. It is found that the simulation result of LBM agrees quite well with the Rayleigh–Plesset prediction.

FIG. 2.

Comparison between the normalized simulation result of LBM and Rayleigh–Plesset prediction.

FIG. 2.

Comparison between the normalized simulation result of LBM and Rayleigh–Plesset prediction.

Close modal

Philipp12 referred that when d > 2.2Rini, the effect of the wall on the evolution of the cavitation bubble can be neglected. The influence of different parameters on the evolution characteristics of the cavitation bubble is studied in detail in this part.

For the LBM pseudo-potential model with non-ideal EOS, the pressure and surface tension are related to the temperature and the liquid/gas density ratio, which makes it difficult to compare the maximum pressure, maximum micro-jet velocity, and lifetime of the cavitation bubble directly under different density ratios. Therefore, a dimensionless initial pressure Δp′ = Δpmaxpe is introduced to study the effect of the density ratio on cavitation bubble collapse, where Δpmax is the maximum pressure difference that can maintain the numerical stability during the cavitation bubble collapse process under different density ratios, and Δpe is the pressure difference of the cavitation bubble in the equilibrium state. Four different density ratios ρl/ρg = 12.8, 40, 135.6, and 720 are chosen to study the effect of the liquid/gas density ratio, and the change in Δp′ with density ratio is shown in Fig. 3. It is found that Δp′ decreases with the increase in the density ratio. The surface tension of the cavitation bubble increases with the increase in the density ratio of the LBM pseudo-potential model, leading to the increase in the equilibrium pressure difference Δpe. Due to the density gradient on the interface, the numerical stability decreases with the increase in the density ratio, which also makes the dimensionless initial pressure Δp′ decrease with the increase in the density ratio.

FIG. 3.

The effect of the liquid/gas density ratio on the dimensionless pressure Δp′.

FIG. 3.

The effect of the liquid/gas density ratio on the dimensionless pressure Δp′.

Close modal

The fluid viscosity can be changed with the relaxation time in LBM. In order to study the influence of fluid viscosity on the evolution process, different relaxation times are chosen for gas and liquid. The initial pressure difference Δp = ppv is set as 0.002 53 mu · lu−1tu−2 and dimensionless distance γ is 2.0. The lifetime tmax, maximum micro-jet velocity umax, and maximum pressure pmax during the evolution process vary with gas relaxation times τg and liquid relaxation times τl, which are shown in Fig. 4.

FIG. 4.

The effect of (a) τg and (b) τl on the lifetime of the bubble tmax, maximum microjet velocity umax and maximum pressure pmax with Δp = 0.002 53 mu · lu−1tu−2.

FIG. 4.

The effect of (a) τg and (b) τl on the lifetime of the bubble tmax, maximum microjet velocity umax and maximum pressure pmax with Δp = 0.002 53 mu · lu−1tu−2.

Close modal

First, we set the liquid relaxation time τl as 0.5375, and the gas relaxation time τg is set from 0.5375 to 1.0625, with the corresponding gas viscosity νg ranging from 0.0125 to 0.1875, to discuss the effect of gas viscosity. As shown in Fig. 4(a), the numerical results indicate that the lifetime tmax, maximum velocity umax, and maximum pressure pmax of different gas relaxation times τg are consistent with each other. Then, we set the gas relaxation time τg as 1.0625, and the liquid relaxation times τl ranging from 0.5375 to 0.7, with the corresponding liquid viscosities νl ranging from 0.0125 to 0.0667. The results show that the lifetime of the bubble decreases from 769 to 714 with νl decreasing, while umax and pmax decrease with νl increasing. Maximum pressure pmax increases about 11 times with τl change from 0.021 mu · lu−1tu−2 to 0.266 mu · lu−1tu−2, while νl decreases from 0.0667 to 0.0125. Meanwhile, the maximum micro-jet velocity increases 110% from 0.448 lu · tu−1 to 0.943 lu · tu−1 with the reduction in νl from 0.0667 to 0.0125. The numerical results are consistent with the conclusion of Popinet41 and Minsier,42 and lower liquid viscosity is preferred to obtain a higher cavitation pressure and larger micro-jet velocity.

To investigate the influence of initial pressure difference between the bubble and ambient liquid on the evolution process, five different initial pressure differences are chosen from 0.0013 mu · lu−1tu−2 to 0.0068 mu · lu−1tu−2 in our study.39Figure 5 shows that the deformation index changes with dimensionless time under different Δp, the maximum Ra stays in the range between 1.1 and 1.2, and the deformation index Ra decreases with the increase in the initial pressure difference Δp. Figure 6 shows the lifetime tmax, maximum micro-jet umax, and maximum pressure pmax change with an initial pressure difference Δp. The maximum pressure pmax and maximum micro-jet umax velocity increase with the increase in Δp. When Δp increases from 0.0013 mu · lu−1tu−2 to 0.0068 mu · lu−1tu−2, umax increases 5.19 times from 0.47 lu · tu−1 to 2.44 lu · tu−1, and pmax increases 88.57 times from 0.014 mu · lu−1tu−2 to 1.24 mu · lu−1tu−2. On the contrary, it is found that the collapse time of the cavitation bubble decreases with the increase in the initial pressure difference, when Δp increases from 0.0013 mu · lu−1tu−2 to 0.0068 mu · lu−1tu−2 and tmax decreases from 925 tu to 487 tu.

FIG. 5.

The evolution of the deformation parameter with dimensionless time with different initial pressure differences Δp.

FIG. 5.

The evolution of the deformation parameter with dimensionless time with different initial pressure differences Δp.

Close modal
FIG. 6.

The effect of the initial pressure difference Δp on lifetime tmax, maximum micro-jet velocity umax, and maximum pressure pmax.

FIG. 6.

The effect of the initial pressure difference Δp on lifetime tmax, maximum micro-jet velocity umax, and maximum pressure pmax.

Close modal

Five dimensionless distances γ are chosen to study their relationship with the lifetime tmax. Figure 7(a) shows that when the bubble stays in the near-wall region, where d < 2.2Rini, tmax increases linearly when γ decreases. The rigid wall affects the evolution of the pressure field between the wall and cavitation bubble, making the deformation index Ra larger, as shown in Fig. 7(b).

FIG. 7.

The effect of the dimensionless distance γ between the bubble center and rigid boundary on (a) the lifetime of bubble and (b) deformation parameters change.

FIG. 7.

The effect of the dimensionless distance γ between the bubble center and rigid boundary on (a) the lifetime of bubble and (b) deformation parameters change.

Close modal

The evolution process of the cavitation bubble is presented in this section where the initial pressure difference Δp = 0.0038 mu · lu−1tu−2, the initial radius Rini = 50 lu, and the corresponding gas–liquid viscosity ratio is 15. The evolution process of the cavitation bubble obtained by our simulation agrees well with the experimental result,12 as shown in Fig. 8.

FIG. 8.

Collapse process of simulation results: (a) t = 530 tu, (b) t = 550 tu, (c) t = 580 tu, (d) t = 590 tu, (e) t = 600 tu, and (f) t = 610 tu, and the experimental results12 with Δp = 0.0038 mu · lu−1tu−2.

FIG. 8.

Collapse process of simulation results: (a) t = 530 tu, (b) t = 550 tu, (c) t = 580 tu, (d) t = 590 tu, (e) t = 600 tu, and (f) t = 610 tu, and the experimental results12 with Δp = 0.0038 mu · lu−1tu−2.

Close modal

Due to the influence of the rigid boundary, the vertical compression rate of the cavitation bubble is larger than the horizontal one with evolution time from t = 0 to t = 470 tu. At this collapse stage, the velocity on the interface is mainly the radial velocity along the radius of the centroid of the bubble, making the cavitation bubble change from spherical to ellipsoid, with the deformation parameter Ra varying from 1 to 1.104. Then, a dent occurred at the top of the interface of the bubble, the shrinkage rate in the vertical direction is higher than that in the horizontal direction. The height of the cavitation bubble decreases rapidly, and Ra changes from 1.104 to 0.

The flow field around the cavitation bubble is shown in Fig. 9. A high-pressure zone is formed, while a dent occurs at the top of the interface at t = 500 tu. In addition, the pressure of the area below the cavitation bubble is affected by the wall with a low-pressure zone formed here. With the dramatic deformation of the cavitation bubble in the collapse stage, the area of the high-pressure region and the maximum pressure pmax increase rapidly, while the pressure difference between the high-pressure region and low-pressure region increases rapidly, which accelerates the deformation of the cavitation bubbles and eventually leads to the collapse of the cavitation bubbles. In addition, the pressure at the collapse point increases sharply to 0.529 mu · lu−1tu−2, and the shock wave is formed. In addition, the maximum velocity and maximum pressure variation with time during the cavitation bubble collapse are shown in Fig. 10. The maximum pressure and the maximum velocity increase slowly in the early collapse process. However, with the violent deformation of the cavitation bubble during t = 500−610 tu, the maximum velocity increases rapidly to 1.36 lu · tu−1, and the maximum pressure increased to 0.529 mu · lu−1tu−2.

FIG. 9.

The flow field changes at the collapse stage: (a) t = 100 tu, (b) t = 400 tu, (c) t = 500 tu, (d) t = 550 tu, (e) t = 580 tu, and (f) t = 610 tu at Δp = 0.0038 mu · lu−1tu−2.

FIG. 9.

The flow field changes at the collapse stage: (a) t = 100 tu, (b) t = 400 tu, (c) t = 500 tu, (d) t = 550 tu, (e) t = 580 tu, and (f) t = 610 tu at Δp = 0.0038 mu · lu−1tu−2.

Close modal
FIG. 10.

Deformation parameter Ra, maximum pressure pmax, maximum velocity umax changes with time at Δp = 0.0038 mu · lu−1tu−2.

FIG. 10.

Deformation parameter Ra, maximum pressure pmax, maximum velocity umax changes with time at Δp = 0.0038 mu · lu−1tu−2.

Close modal

The evolution of the cavitation bubble is a process with energy accumulation and release, and the TKE of the cavitation bubble is introduced to describe the collapse process from the perspective of energy. The surrounding liquid works on the cavitation bubble through pressure; as the volume of the cavitation bubble shrinks, the distance between gas molecules decreases, which makes the molecular potential energy and molecular total kinetic energy increase. In addition, the energy release to the surrounding liquid in a short time with the phase change happens and shock waves form on the macroscope. The definition of TKE of the cavitation bubble is described as follows:51 

TKE=ρv|u|2.
(15)

The density cutoff value ρv ≤ (ρg + ρl)/2 is used to identify the gas phase. The evolution of the cavitation bubble can be divided into two stages from the perspective of energy, energy accumulation, and release. The energy accumulation with the cavitation bubble shrinks in most of the collapse time. However, the TKE releases sharply in a short time with the gas phase transfer into liquid. As shown in Fig. 11, the time of accumulation stage of TKE is from t = 0 to t = 585 tu, reaching a maximum value of 1.041 mu · lu−1tu−2. Later, from t = 585 tu to t = 611 tu, TKE releases radically in a short time, decreasing from 1.041 mu · lu−1tu−2 to 0.041 mu · lu−1tu−2.

FIG. 11.

TKE varies during the collapse process with Δp = 0.0038 mu · lu−1tu−2.

FIG. 11.

TKE varies during the collapse process with Δp = 0.0038 mu · lu−1tu−2.

Close modal

In most natural phenomena, cavitation occurs in large numbers; bubbles not only interact with the boundary, but also interact with each other.9 Previous studies indicate the interaction between multiple cavitation bubbles, and the wall may change the direction of the micro-jet,52 the propagation process of the shock wave,53 and the wall-effect range.52 

The interaction of two cavitation bubbles near a rigid boundary is studied in the present study. The computational domain is a 401 × 601 lu2 rectangular area, as shown in Fig. 12, and Zou–He pressure inlet boundary conditions49 are applied in the top of the domain; both the left and the right side are non-equilibrium extrapolation boundaries, and the no-slip boundary46 is used in the bottom. The density ratio is set as ρl/ρg ≈ 720, and the relaxation time of the gas phase is τg = 1.0625 with the liquid phase set as τg = 0.5375, making the ratio of the viscosity νg/νl = 15. All these make the gas–liquid density ratio and the viscosity ratio close to the actual cavitation phenomena. In our simulation, the initial pressure difference between the bubble and ambient liquid is Δp = 0.0025 mu · lu−1tu−2. Three cases are simulated, and the corresponding parameters are shown in Table II. Two normalized parameters, γ1 = d/Rini2 and γ2 = d2/Rini1, are introduced to describe the effect of the rigid boundary, where Rini1 and Rini2 are the radii of the lower bubble and upper bubble, respectively.

FIG. 12.

Schematic of the computational domain for the evolution of the interaction of two cavitation bubbles.

FIG. 12.

Schematic of the computational domain for the evolution of the interaction of two cavitation bubbles.

Close modal
TABLE II.

Parameters for different interaction cavitation cases.

Rini1 (lu)Rini2 (lu)d (lu)d2 (lu)
Case 1 50 50 311 159 
Case 2 50 50 219 86 
Case 3 50 50 180 50 
Rini1 (lu)Rini2 (lu)d (lu)d2 (lu)
Case 1 50 50 311 159 
Case 2 50 50 219 86 
Case 3 50 50 180 50 

In case 1, γ2 is set as 3.18, which means that the lower bubble stays out of the boundary-effect region, and the simulation results are shown in Fig. 12(a). The result shows that the interaction between cavitation bubbles makes the cavitation bubbles collapse toward each other; it seems there is a “wall” between the upper and lower cavitation bubbles. The upper cavitation bubble collapses a little faster than the lower one, and the differences of maximum deformation parameters Ra and the lifetime of two bubbles are small. The evolution process means that the interaction between two bubbles is stronger than the interaction with the rigid boundary in this case. Compared with the experimental result,9 the result of LBM seems more reasonable than that with the boundary integral method9 in the second column of Fig. 13(a).

FIG. 13.

Comparison of the cavitation bubble shape obtained by LBM results (first row), boundary integral method results (second row), and experimental results (the third row)9 with (a) case 1, (b) case 2, and (c) case 3.

FIG. 13.

Comparison of the cavitation bubble shape obtained by LBM results (first row), boundary integral method results (second row), and experimental results (the third row)9 with (a) case 1, (b) case 2, and (c) case 3.

Close modal

When γ2 is set as 1.72, the lower cavitation bubble is located in the wall-effect region. The simulation result shows that the deformation process of the upper cavitation bubble is unaffected by the wall, and its maximum deformation parameter Ra is close to 1.1. On the contrary, the lower cavitation bubble is affected by both the wall and the upper cavitation bubbles, and it becomes an ellipsoidal bubble with maximum Ra larger than 2.0. The lifetime of the lower bubble is longer than the upper one. The numerical result is shown in Fig. 13(b), and the result of LBM simulation agrees well with those from the boundary integral method and experiment.9 

When γ2 = 1, the lower bubble is in direct contact with the rigid boundary. The numerical result of case 3 is shown in Fig. 13(c). The deformation parameter Ra of the upper cavitation bubble is close to 1.1. Due to the effect of the rigid boundary and upper bubble, the lower cavitation bubbles mainly shrink along the horizontal direction, and the vertical compression is very small.

The deformation parameters of two bubbles varying during the process are shown in Fig. 14. The time-deformation parameter curves of the upper bubble in three cases almost overlap with each other, and the maximum Ra has only a mild difference. However, the time-deformation parameter curves of the lower bubble are different, when the lower cavitation bubble stays outside of the near-wall region; the time-deformation parameter curve is consistent with the upper bubble. When the lower bubble stays in the near-wall region, the shrinkage rate in the horizontal direction is much greater than that in the vertical direction. The smaller the γ2 is, the slower the vertical shrinkage rate is, thus gradually increasing the Ra value.

FIG. 14.

Deformation parameter Ra of two bubbles changes with time in different cases.

FIG. 14.

Deformation parameter Ra of two bubbles changes with time in different cases.

Close modal

In case 1, a dent occurs at the top of the interface of the upper bubble at t = 600 tu, and the high-pressure region is formed, while the pressure between the lower bubble and the rigid boundary is affected by the wall, and a low-pressure region is formed between two bubbles. At t = 700 tu, a dent also occurs at the bottom of the lower bubble, and a high-pressure region is formed between the lower bubble and the wall. It interacts with the high-pressure region of the upper cavitation bubble, making the upper bubble and lower bubble move toward each other. The influence of the wall makes the area of the high-pressure region of the lower bubble change slowly. However, the high-pressure region above the upper bubble increases rapidly, which accelerates the collapse of the upper bubble, as shown in Fig. 15. During the collapse process of two cavitation bubbles, the formation of the microjet can be observed, and the two cavitation bubbles collapse toward each other.

FIG. 15.

The evolution of the flow field in case 1 at (a) t = 200 tu, (b) t = 500 tu, (c) t = 650 tu, (d) t = 700 tu, (e) t = 730 tu, and (f) t = 747 tu.

FIG. 15.

The evolution of the flow field in case 1 at (a) t = 200 tu, (b) t = 500 tu, (c) t = 650 tu, (d) t = 700 tu, (e) t = 730 tu, and (f) t = 747 tu.

Close modal

Compared with case 1, the lower bubble in case 2 is affected by the rigid boundary and the upper cavitation bubbles. The low-pressure regions generated at both the top and bottom of the lower bubble make the lower bubble compress in the horizontal direction, and the lower bubble becomes an ellipsoid type. The evolution of the flow field near the upper bubble is similar to case 1. When t = 650 tu, a dent occurs at the top of the upper bubble, and a local high-pressure region appears at the same time. The maximum pressure value of the high-pressure zone increases rapidly, making the upper cavitation bubble collapse in a short time. The micro-jet is formed and flows in the direction toward the lower bubble, as shown in Fig. 16.

FIG. 16.

The evolution of the flow field in case 2 at (a) t = 200 tu, (b) t = 500 tu, (c) t = 650 tu, (d) t = 700 tu, (e) t = 750 tu, and (f) t = 781 tu.

FIG. 16.

The evolution of the flow field in case 2 at (a) t = 200 tu, (b) t = 500 tu, (c) t = 650 tu, (d) t = 700 tu, (e) t = 750 tu, and (f) t = 781 tu.

Close modal

The TKE of both upper and the lower bubbles in case 1 increases slowly in the early stage of the collapse process. Because the two bubbles stay outside of the wall-affect region, the curve of TKE in the early collapse stage overlaps with each other. However, with a dent occurring at the top of the upper bubble, the TKE of the upper bubble grows faster. The maximum value of TKE reached near 0.6 mu · lu−1tu−2 before the release process, and the energy releases rapidly with the bubble collapse later. The development of the high-pressure field near the lower bubble is affected by the wall; this decreases the TKE growth rate of the lower bubble. When the upper bubble collapses, the TKE of the lower bubble is 0.4 mu · lu−1tu−2, as shown in Fig. 17(a). In case 2, when the lower bubble is affected by the wall and the upper bubble, no high-pressure region is formed around it. When the upper bubble collapses, the lower bubble still stays in the collapse process. Before TKE of the upper bubble releases, it reaches the maximum value of 0.58 mu · lu−1tu−2, while the TKE of the lower bubble is only 0.29 mu · lu−1tu−2, as shown in Fig. 17(b).

FIG. 17.

The evolution of TKE in (a) case 1 and (b) case 2.

FIG. 17.

The evolution of TKE in (a) case 1 and (b) case 2.

Close modal

In this paper, an MRT pseudo-potential model with a high-density ratio and a high viscosity ratio is used to simulate the collapse process of a single bubble and a dual bubble collapse near a rigid boundary, with density ratio ρl/ρg = 720. The simulation agrees well with experiments and theoretical analysis, which means that the proposed LBM pseudo-potential model is a robust and effective tool for the study of the collapse of near-wall cavitation bubbles. Effects of the liquid/gas density ratio, viscosity, initial pressure difference, and bubble distance to the wall on the evolution of the cavitation bubble are studied. The deformation of cavitation bubbles, micro-jet velocity, and the maximum collapse pressure are analyzed. The TKE is first introduced to describe the formation of the shock wave from the perspective of energy. Simulation results indicate that

  1. The gas viscosity νg does not affect the maximum pressure and maximum micro-jet velocity in the cavitation process, while the liquid viscosity νl does. The maximum velocity and pressure decrease with an increase in the liquid viscosity νl, while the lifetime of the bubble increases. The initial pressure difference Δp mainly affects the time of the compression stage, but it has a little effect on the time of the collapse stage during cavitation bubble evolution. With an increase in Δp, the maximum velocity and pressure increase, while the lifetime of the cavitation bubble decreases.

  2. In the near-wall region, where d < 2.2Rini, the lifetime of the cavitation bubble decreases with the increase in the distance ration γ. The closer the cavitation bubble to the wall, the larger the deformation parameter Ra, and the longer the lifetime of the bubble.

  3. The evolution of the cavitation bubble can be divided into two stages with TKE introduced in the present study: energy accumulation and energy release. The shock wave is formed with TKE release in a short time.

  4. When two cavitation bubbles interact with each other in the near-wall region, the deformation curves and energy curves of the upper cavitation bubbles of different cases are similar to each other. The evolution of the lower cavitation bubble is related to the distance between the lower cavitation bubble and the wall. When γ2 > 2.2, the evolution of the lower cavitation bubbles is not affected by the wall, but the high-pressure of the lower cavitation bubble is affected by the rigid boundary, making the lower cavitation bubble collapse toward the upper one, and its lifetime is longer than that of the upper one. When γ2 < 2.2, the evolution of the lower cavitation bubble is affected by both the wall and the upper cavitation bubble, and low-pressure regions are formed at the top and bottom of the lower cavitation bubble. In addition, the lifetime of the lower cavitation bubble becomes much longer than that of the upper one.

This work was supported by the National Science Fund for Distinguished Young Scholars (Grant No. 51625901) and the National Natural Science Foundation of China (Grant No. 51879176).

1.
M.
Plesset
and
R.
Charpman
, “
Collapse of an initially spherical vapour cavity in the neighbourhood of a solid boundary
,”
J. Fluid Mech.
47
,
283
290
(
1971
).
2.
J.
Blake
,
B.
Taib
, and
G.
Doherty
, “
Transient cavities near boundaries
,”
J. Fluid Mech.
170
,
479
497
(
1986
).
3.
A.
Aganin
,
M.
Ilgamov
,
L.
Kosolapova
, and
G.
Malakhov
, “
Dynamics of a cavitation bubble near a solid wall
,”
Thermophys. Aeromech.
23
,
211
220
(
2016
).
4.
Y.
Shen
,
K.
Yasui
,
Z.
Sun
,
B.
Mei
,
M.
You
, and
T.
Zhu
, “
Study on the spatial distribution of the liquid temperature near a cavitation bubble wall
,”
Ultrason. Sonochem.
29
,
394
400
(
2016
).
5.
Z.
Izadifar
,
P.
Babyn
, and
D.
Chapman
, “
Ultrasound cavitation/microbubble detection and medical applications
,”
J. Med, Biol. Eng.
39
,
259
(
2018
).
6.
S.
Mishra
,
T.
Lacy
, and
S.
Kundu
, “
Effect of surface tension and geometry on cavitation in soft solids
,”
Int. J. Nonlinear Mech.
98
,
23
31
(
2018
).
7.
W.
Lauterborn
and
H.
Bolle
, “
Experimental investigations of cavitation-bubble collapse in the neighbourhood of a solid boundary
,”
J. Fluid Mech.
72
,
391
399
(
1975
).
8.
J.
Blake
and
D.
Gibson
, “
Cavitation bubbles near boundaries
,”
Annu. Rev. Fluid Mech.
19
,
99
123
(
1987
).
9.
J.
Blake
,
P.
Robinson
,
A.
Shima
, and
Y.
Tomita
, “
Interaction of two cavitation bubbles with a rigid boundary
,”
J. Fluid Mech.
255
,
707
721
(
1993
).
10.
T.
Ogasawara
,
N.
Tsubota
,
H.
Seki
,
Y.
Shigaki
, and
H.
Takahira
, “
Experimental and numerical investigations of the bubble collapse at the center between rigid walls
,”
J. Phys.: Conf. Ser.
656
,
012031
(
2015
).
11.
N.
Bremond
,
M.
Arora
,
S.
Dammer
, and
D.
Lohse
, “
Interaction of cavitation bubbles on a wall
,”
Phys. Fluids
18
,
121505
(
2006
).
12.
A.
Philipp
and
W.
Lauterborn
, “
Cavitation erosion by single laser-produced bubbles
,”
J. Fluid Mech.
361
,
75
116
(
1998
).
13.
O.
Lindau
and
W.
Lauterborn
, “
Cinematographic observation of the collapse and rebound of a laser-produced cavitation bubble near a wall
,”
J. Fluid Mech.
479
,
327
348
(
2003
).
14.
L.
Calvisi
,
O.
Lindau
, and
J.
Blake
, “
Shape stability and violent collapse of microbubbles in acoustic traveling waves
,”
Phys. Fluids
19
,
047101
(
2007
).
15.
T.
Brlansky
and
L.
Calvisi
, “
Shape mode stability of lipid-coated ultrasound contrast agent microbubbles
,”
J. Acoust. Soc. Am.
139
,
2093
(
2016
).
16.
S.
Fujikawa
and
H.
Takahira
, “
Dynamics of two nonspherical cavitation bubbles in liquids
,”
J. Fluid Dyn. Res.
4
,
179
194
(
1988
).
17.
F.
Kunz
,
A.
Boger
,
R.
Stinebring
,
S.
Chyczewski
,
W.
Lindau
, and
J.
Gibeling
, “
A preconditioned Navier–Stokes method for two-phase flows with application to cavitation prediction
,”
Comput. Fluids
29
,
849
875
(
2000
).
18.
A.
Osterman
,
M.
Dular
, and
B.
Sirok
, “
Numerical simulation of a near-wall bubble collapse in an ultrasonic field
,”
J. Fluid Sci. Technol.
4
,
210
221
(
2009
).
19.
S.
Yakubov
,
T.
Maquil
, and
T.
Rung
, “
Experience using pressure-based CFD methods for Euler–Euler simulations of cavitating flows
,”
Comput. Fluids
111
,
91
104
(
2015
).
20.
E.
Goncalves
,
M.
Champagnac
, and
R.
Patella
, “
Comparison of numerical solvers for cavitating flows
,”
Int. J. Comput. Fluid D
24
,
201
216
(
2010
).
21.
E.
Ghahramani
,
M.
Arabnejad
, and
R.
Bensow
, “
A comparative study between numerical methods in simulation of cavitating bubbles
,”
Int. J. Multiphase Flow
111
,
339
359
(
2019
).
22.
G.
Chahine
,
A.
Kapahi
,
J. K.
Choi
, and
C. T.
Hsiao
, “
Modeling of surface cleaning by cavitation bubble dynamics and collapse
,”
Ultrason. Sonochem.
29
,
528
549
(
2016
).
23.
C. T.
Hsiao
,
J. S.
Ma
, and
G.
Chahine
, “
Multiscale tow-phase flow modeling of sheet and cloud cavitation
,”
Int. J. Multiphase Flow.
90
,
102
117
(
2017
).
24.
S.
Yakubov
,
B.
Cankurt
,
M.
Abdel-Maksoud
, and
T.
Rung
, “
Hybrid MPI/OpenMP parallelization of an Euler–Lagrange approach to cavitation modeling
,”
Comput. Fluids
80
,
365
371
(
2013
).
25.
D.
Grunau
,
S.
Chen
, and
K.
Eggert
, “
A lattice Boltzmann model for multiphase fluid flows
,”
Phys. Fluids
5
,
2557
2562
(
1993
).
26.
X.
Shan
and
H.
Chen
, “
Lattice Boltzmann model for simulating flows with multiple phases and components
,”
Phys. Rev. E
47
,
1815
1819
(
1993
).
27.
X.
Shan
and
H.
Chen
, “
Simulation of nonideal gases and liquid-gas phase transitions by the lattice Boltzmann equation
,”
Phys. Rev. E
49
,
2941
2948
(
1994
).
28.
M.
Swift
,
W.
Osborn
, and
J.
Yeomans
, “
Lattice Boltzmann simulation of nonideal fluids
,”
Phys. Rev. Lett.
75
,
830
833
(
1995
).
29.
M.
Swift
,
E.
Orlandini
,
W.
Osborn
, and
J.
Yeomans
, “
Lattice Boltzmann simulations of liquid-gas and binary fluid systems
,”
Phys. Rev. E
54
,
5041
5052
(
1996
).
30.
X.
He
,
S.
Chen
, and
R.
Zhang
, “
A lattice Boltzmann scheme for incompressible multiphase flow and its application in simulation of Rayleigh–Taylor instability
,”
J. Comput. Phys.
152
,
642
663
(
1999
).
31.
M.
Mazloomi
,
S.
Chikatamarla
, and
I.
Karlin
, “
Entropic lattice Boltzmann method for multiphase flows
,”
Phys. Rev. Lett.
114
,
174502
(
2015
).
32.
R.
Zhang
and
H.
Chen
, “
Lattice Boltzmann method for simulations of liquid-vapor thermal flows
,”
Phys. Rev. E
67
,
066711
(
2003
).
33.
L.
Chen
,
Q.
Kang
,
Y.
Mu
,
Y. L.
He
, and
W. Q.
Tao
, “
A critical review of the pseudopotential multiphase lattice Boltzmann model: Methods and applications
,”
Int. J. Heat Mass Transfer
76
,
210
236
(
2014
).
34.
M.
Sukop
and
D.
Or
, “
Lattice Boltzmann method for homogeneous and heterogeneous cavitation
,”
Phys. Rev. E
71
,
046703
(
2005
).
35.
S.
Mishra
,
P.
Deymier
,
K.
Muralidharan
,
G.
Frantziskonis
,
S.
Pannala
, and
S.
Simunovic
, “
Modeling the coupling of reaction kinetics and hydrodynamics in a collapsing cavity
,”
Ultrason. Sonochem.
17
,
258
265
(
2010
).
36.
X.
Chen
,
C.
Zhong
, and
X.
Yuan
, “
Lattice Boltzmann simulation of cavitating bubble growth with large density ratio
,”
Comput. Math. Appl.
61
,
3577
3584
(
2011
).
37.
M.
Shan
,
C.
Zhu
,
X.
Zhou
,
C.
Yin
, and
Q.
Han
, “
Investigation of cavitation bubble collapse near rigid boundary by lattice Boltzmann method
,”
J. Hydrodyn.
28
,
442
450
(
2016
).
38.
Y.
Mao
,
Y.
Peng
, and
J.
Zhang
, “
Study of cavitation bubble collapse near a wall by the modified lattice Boltzmann method
,”
Water
10
,
1439
(
2018
).
39.
M.
Shan
,
C.
Zhu
,
C.
Yao
,
C.
Yin
, and
X.
Jiang
, “
Pseudopotential multi-relaxation-time lattice Boltzmann model for cavitation bubble collapse with high density ratio
,”
Chin. Phys. B
25
,
104701
(
2016
).
40.
H.
Xue
,
F.
Shan
, and
X.
Guo
, “
Cavitation bubble collapse near a curved wall by the multiple-relaxation-time Shan-Chen lattice Boltzmann model
,”
Chin. Phys. Lett.
34
,
084301
(
2017
).
41.
S.
Popinet
and
S.
Zaleski
, “
Bubble collapse near a solid boundary: A numerical study of the influence of viscosity
,”
J. Fluid Mech.
464
,
137
163
(
2002
).
42.
V.
Minsier
,
J.
De Wilde
, and
J.
Proost
, “
Simulation of the effect of viscosity on jet penetration into a single cavitating bubble
,”
J. Appl. Phys.
106
,
084906
(
2009
).
43.
Q.
Li
,
Q.
Kang
,
M.
Francois
,
Y.
He
, and
K.
Luo
, “
Lattice Boltzmann modeling of boiling heat transfer: The boiling curve and the effects of wettability
,”
Int. J. Heat Mass Transfer
85
,
787
796
(
2015
).
44.
Z.
Yu
and
L.
Fan
, “
Multi-relaxation-time interaction-potential-based lattice Boltzmann model for two-phase flow
,”
Phys. Rev. E
82
,
046708
(
2010
).
45.
Q.
Li
,
K.
Luo
,
Q.
Kang
,
Y.
He
,
Q.
Chen
, and
Q.
Liu
, “
Lattice Boltzmann methods for multiphase flow and phase-change heat transfer
,”
Prog. Energy Combust. Sci.
52
,
62
105
(
2016
).
46.
Q.
Li
,
K.
Luo
,
Q.
Kang
, and
Q.
Chen
, “
Contact angles in the pseudo-potential lattice Boltzmann modeling of wetting
,”
Phys. Rev. E
90
,
053301
(
2014
).
47.
X.
Shan
, “
Pressure tensor calculation in a class of nonideal gas lattice Boltzmann models
,”
Phys. Rev. E
77
,
066702
(
2008
).
48.
P.
Yuan
and
L.
Schaefer
, “
Equations of state in a lattice Boltzmann model
,”
Phys. Fluids
18
,
042101
(
2006
).
49.
Q.
Zou
and
X.
He
, “
On pressure and velocity boundary conditions for the lattice Boltzmann BGK model
,”
Phys. Fluids
9
,
1591
1598
(
1997
).
50.
Q.
Li
,
K.
Luo
, and
X.
Li
, “
Lattice Boltzmann modeling of multiphase flows at large density ratio with an improved pseudopotential model
,”
Phys. Rev. E
87
,
053301
(
2013
).
51.
H.
Dalgamoni
and
X.
Yong
, “
Axisymmetric lattice Boltzmann simulation of droplet impact on solid surfaces
,”
Phys. Rev. E
98
,
013102
(
2018
).
52.
L. W.
Chew
,
E.
Klaseboer
, and
S. W.
Ohl
, “
Interaction of two oscillating bubbles near a rigid boundary
,”
Exp. Therm. Fluid Sci.
44
,
108
113
(
2013
).
53.
R.
Betney
,
B.
Tully
,
A.
Hawker
, and
Y.
Ventikos
, “
Computational modelling of the interaction of shock waves with multiple gas-filled bubbles in a liquid
,”
Phys. Fluids
27
,
036101
(
2015
).