Cavitation is ubiquitous in liquid, and happens when the local pressure is below the saturated vapor pressure. As cavitation bubble collapse near a solid wall, the associated phenomena include instant high pressure, high velocity jets and high temperature, which closely relate with the cavitation erosion of the solid material surface. On the other side, the collapse of the cavitation bubble has been applied in environmental protection, ultrasonic therapy, lab on chip and material surface cleaning   .
The mechanism of the cavitation bubble collapse near solid wall is a fundamental issue for the above applications. However, as too many phenomena are involved, theoretical model of cavitation bubble collapse is difficult to be established, and for complex geometry boundary, the analytical solution is even impossible. Many publications are devoted to the investigation of the bubble collapse near the planar wall starting from the experimental work  . In recent years, many works have been developed to investigate the interaction between collapsing bubble and non-planar solid wall   , and numerical methods are becoming more and more important tools to study the bubble collapse near solid wall   . Commonly used numerical methods include the finite volume method (FVM)  , the finite element method (FEM)  , and the boundary element method (BEM)  . These macroscopic numerical modeling methods based on solving partial differential equations are limited in processing the multiphase flows and complex geometry boundary conditions. As describing the multiphase flow, macroscopic methods need the assistance of the schemes of the interface tracking or capturing, which will reduce the computational efficiency. For the complex geometry boundary, it is difficult and inefficient to implement by macroscopic methods.
Owing to the flexibility for complex geometry boundary and the simplicity of the algorithm, the lattice Boltzmann method (LBM) has been developed into a powerful tool for the flow, heat and mass transfer simulations relating with com- plex geometry boundary   . In LBM community, many multiphase models have been presented, which can be generally classified as the color-gradient model  , the pseudopotential model (or Shan-Chen model)   , the free-energy model  and the phase-field model. The pseudopotential model is widely and successfully used in the LBM multiphase community due to its conceptual simplicity and computation efficiency   . In the pseudopotential model, the fluid interactions are mimicked by an interparticle potential, from which a non-monotomic equation of state (EOS) can be obtained. As a result, the separation of fluid phases or components can be achieved automatically in this method, and the methods to track or capture the interfaces are not required.
In recent years, the pseudopotential model, as the top choice of the multiphase LB model, was introduced into the issue of cavitation by Sukop and Or  . Subsequently, several research efforts emerged to investigate the mechanism of cavitation. Chen  simulated the cavitation bubble growth using the modified pseudopotential LB model with the exact difference method (EDM) force scheme. The results in quiescent flows agree fairly well with the solution of Rayleigh- Plesset equation. Mishra  introduced a model of cavitation based on the pseudopotential LB model that allows for coupling between the hydrodynamics of a collapsing cavity and supported solute chemical species. Unfortunately, the above researches do not involve the interaction between bubble and solid wall. Until, the pseudopotential model was introduced into the researches on the mechanism investigations of cavitation bubble collapse near solid wall by Shan et al.   . However, these previous works involve just planar wall. In engineering cases, the rough walls are ubiquitous. The mechanism investigation of cavitation bubble collapse near rough wall is the basic issue for the efficient applications of cavitation bubble collapse.
One of the key challenges of the simulation of cavitation bubble collapse near rough wall is the stability of the LB model. For pseudopotential multiphase model, the stability is closely related to the thermodynamic consistency  , and great efforts have been made for this issue    . Recently, Li et al.    founded that there exists a suitable forcing scheme, which can meet the thermodynamic consistency requirement in an efficient way. In order to investigate the interaction mechanism between the collapsing bubble and complex geometry boundary, in the present work, the improved forcing scheme for the pseudopotential MRT LB model developed by Li et al.  is adopted to achieve the sufficient density ratio and stability of the numerical model.
2. Numerical Method
The pseudopotential LB model, also well known as Shan-Chen model, was developed by Shan and Chen in 1993  . In pseudopotential model, the fluid interactions are mimicked by an interparticle potential, which is now widely called pseudopotential. In original pseudopotential LB model, the Single-Relaxation- Time (SRT) collision operator was employed. In recent years, the MRT collision operator has been verified that it is superior to the SRT operator in terms of numerical stability. The MRT LB evolution equation can be given as follows:
where is the density distribution function, is its equilibrium distribution, is the time, is the spatial position, is the discrete velocity along the direction, is the time step, is the forcing term in the velocity space, and is an orthogonal transformation matrix. in Equation (1) is a diagonal matrix, and for D2Q9 lattice, is given by
Through the transformation matrix , and can be projected onto the moment space via and , and the collision step of MRT LB equation (Equation (1)) can be rewritten as 
where is the unit tensor, and is the forcing term in the moment space with . For D2Q9 lattice, can be given by
Here is the macroscopic density, is the macroscopic velocity which satisfies and is calculated by
where for two dimensional space is the force action on the fluid system. Then the streaming step of the MRT-LB equation can be formulated as
For the pseudopotential LB model in D2Q9 lattice case, the in Equation (5) is given by
where is the interparticle potential, is the interaction strength, and and are the weights. In the present work, the form of proposed by Yuan et al.  is adopted which can be formulated as
where is the pressure calculated by equation of state (EOS). The in Equation (8) loses the meaning of the interaction strength and is used to ensure that the whole term inside the square root is positive  . The Carnahan-Starl- ing (CS) EOS is adopted in the present work, which can be given by 
where and . Here and are the critical temperature and pressure, respectively.
in Equation (7) can be incorporated in evolution equation via with specific forcing scheme. Li et al.  proposed a MRT version forcing scheme to achieve thermodynamic consistency. For the D2Q9 lattice, Li’s forcing scheme can be given by
where is a parameter related with the particle interaction range  and is proved by Li et al. to have the function of adjusting the thermodynamic consistency  . Li’s forcing scheme has the following advantages: a) maintaining a uniform layout with the a general form of the LB forcing scheme; b) achieving thermodynamic consistency only by tuning one constant parameter; and c) fully retaining the LBM’s advantages of simple and efficient. In  it is found that the thermodynamic consistency is independent of kinematic viscosity for Li’s improved forcing scheme, and the surface tension is independent of the relaxation time . These features make it more convenient to investigate the physical mechanism of the multiphase flows.
Unless otherwise specified, the unit adopted in this paper is the lattice unit of LBM, i.e. the units of the length, time, mass, velocity, density and pressure are lu(lattice unit), ts(time step), mu(mass unit), mu lu−2, lu ts−1 and mu ts−2.
3. Cavitation Bubble Collapse near Rough Solid Wall
It is confirmed that the geometry is a crucial role on the process of the cavitation bubble collapse  . The interaction between collapse bubble and rough solid wall is a common issue abstracted from the cavitation applications involving the interactions between bubble and complex geometry boundary. In the present work, the rough solid wall is described by the wall with the geometries of periodic grooves with equal widths as shown in Figure 1. In this diagram, and are the depth and the width of groove, respectively. is the width of bulge, and is the period width of the periodically arranged geometry.
The computational domain is shown in Figure 2. R0, and are the initial radius, the pressure inside bubble and the ambient pressure, respectively. and are the distance from bubble center to groove bottom and bulge top, respectively. Based on the concept of the stand-off parameter  , two non- dimensional parameters are introduced as follows to describe the positional relations between bubble and rough solid wall
Figure 1. Schematic diagram of rough wall geometry.
Figure 2. Computational domain for single solid wall case.
Here, , express the stand-off distances between bubble and bulge top, groove bottom, respectively.
The pressure boundary conditions, implemented by the nonequilibrium extra- polation scheme  , are adopted in the directions of left and right. And, a constant pressure boundary condition is adopted in the top of computational domain by Zou-He scheme  .
A lattice is adopted in the simulations of this section. A vapor bubble with radius of is initially placed at the center of the domain. The density field is initialized as 
where is the center of the bubble, is the prescribed width of the phase interface and is set as 5 in present work, is the hyperbolic tangent function. The parameters of solid wall geometry are set as , , , and , respectively.
The parameters for the present MRT pseudopotential LB model are chosen as follows: , , , and . For CS EOS, , and . The temperature . In order to simulate the bubble collapse process, a positive pressure difference is achieved by artificially tuning the initial liquid density based on the equilibrium state. In this section, the pressure difference is .
4.1. Evolution of Density Field
The evolution of density field is shown in Figure 3. The deformation of bubble profile can be investigated from density field. The initial spherical bubble begins
Figure 3. Density field evolution of the collapsing bubble near rough solid wall.
to collapse motivated by the pressure difference from t = 0. The velocity of collapse is slow at the starting stage. Until t = 202, the bubble profile appears obvious deformation, i.e., the radius of curvature reduce towards the rough solid wall. Then, the bubble becomes an ellipsoid and the density just over the bubble increases. Along with the diffusion of the higher density area, the top of the bubble is flatten, and then concaved directing the solid wall. At the same time, the more dense and concentrated density area appears at the concave and accelerate the shift of the bubble profile. As a result, the jet is formed. The first collapse happens when the upper bubble wall clashes the bottom wall. The jet perforates the bottom wall and produces a crucial bubble. The shock wave generated by the first collapse propagates rapidly towards the rough solid wall and bounce back. The bubble experiences the complex and volatile deformation at the last stage under the jet which includes shock wave and bounce wave.
Figure 3 describes the whole process of bubble collapse with time t. Comparing with the flat solid wall case, the bubble collapse near rough solid wall appears the similar dynamics process  . The differences are mainly reflected in two aspects. The first one is the bounce of the shock wave in rough solid wall case is weaker than the one in the planar wall case. The strength and the reflection path of the shock wave are affected by the geometry of the rough solid wall. The second one is the geometry of the rough solid wall affects the bubble deformation, especially at the last stage. In order to be more intuitive, pressure field evolution of the collapsing bubble is investigated in the next section.
4.2. Evolution of Pressure and Velocity Field
Several representative moments of the pressure and velocity evolution are shown in Figure 4. At t = 528, one high pressure area emerges at the concave of bubble top. It is consistent with the analysis of the density distribution at the same moment. At the same time, the fluid velocity is increased towards the solid wall which formats the inceptive jet. The jet velocity at the concave is significantly higher than that of other areas in the moment. Along with the sagging of the upper bubble wall, the jet velocity is higher and higher. At t = 717, the first collapse occurs. The clash between the upper and the bottom walls induces shock wave as shown in t = 741. And at the same moment, the bubble is broken by the high velocity jet. Affected by the vortexes, which are induced by the jet and the bouncing back shock wave, the crucial bubble deforms and then collapses close to the rough solid wall. From Figure 4, we can find that the shock waves bounce intensively at the solid surfaces of the bulges top, and propagate unceasingly into the grooves with decreasing strength. The fluid jet follows the similar mode. It is inevitable that the shock wave and the jet will erode and impact the rough solid wall.
Figure 4. Pressure and velocity field evolution of the collapsing bubble near rough solid wall.
Figure 5 gives the finer detail display of the interaction between the pressure, velocity and the rough solid wall at t = 998. It is found that in the groove and near the bulge, which are off-center the collapse position, the vortexes are formed. Conversely, the laminar flows are formed in the groove just under the collapse position. The vortexes and the laminar flows induce the remarkable eroding effect on the solid surface of the bulge and the groove, respectively.
4.3. Analysis of Eroding Effect
In order to investigate the eroding process and effect induced by bubble collapse close to the rough solid wall, four test sections are set at the bulges top and the grooves side walls as shown in Figure 5. The coordinates of these four test sections are , , and , respectively. The pressure evolution and the velocity evolution of these four test sections, from 730 ts to 1000 ts, are displayed in Figure 6 and Figure 7.
Form Figure 6, we can find that the test section a experiences three pressure peaks during the time interval of observation. The closer to the groove entrance, the greater the pressure amplitude is. The first pressure peak with the greatest amplitude is the result of shock wave. The second is the water hummer effect of jet. However, as the test section is parallel with jet direction, the pressure amplitude is weaker than the first pressure peak. The third pressure peak is induced by the second collapse. The pressure evolution at the test section b displays similar pattern. Due to the impeding effect of the bulge, the lower pressure area is induced near the test section b. The test section c is just under the collapse position, so the three pressure peaks are greater than the test section a and b. The pressure pattern of the test section d is similar as that of b on account of the larger deviation
Figure 5. Pressure and velocity distribution at rough solid wall after the bubble collapse.
Figure 6. Pressure evolution at the test sections.
Figure 7. Velocity evolution at the test sections.
from the collapse station. But the pressure amplitudes are higher than the test section b because of the closer distance from the solid wall. The pressure evolution at the four test sections demonstrates that no matter the bulge solid wall or the groove entrance experience the oscillation of pressure.
The velocity evolution of the four test sections are shown in Figure 7. For the test section a, every test point experiences the similar process of the fluid velocity, which increases sharply and then gradually rises to the peak. The pattern of the velocity is same as the pressure, i.e., the closer to the groove entrance, the greater the velocity is. It demonstrated that the groove geometry impedes both of the pressure and the velocity at the groove side walls. Although the same near the entrance, the test section b shows a different pattern. Under the combined action of shock wave and jet, the vortexes are induced at b. The velocity evolution pattern of the test section c is symmetric in x direction. Due to the non-slid effect of solid wall, the fluid velocity at center of the bulge is lower than the sides, which are closer to the groove entrances. For the test section d, the velocity pattern is more diversified with lower fluid velocity than other test sections. The oscillation of velocity is more violent for farther distance from collapse position. The last velocity peaks are related with the collapse of the crucial bubble.
In summary, the eroding effect at the solid wall of the bulge part is weaker than that at the side walls, and is more obvious at the side walls of the grooves than other sections. From the pressure and velocity evolution patterns, it demonstrates that the fluid velocity accounts for the eroding effect at the side walls of the grooves, and the pressure accounts for the eroding effect at the bulge part of the rough solid wall.
For the modeling of the collapsing bubble near rough solid wall, an improved MRT pseudopotential LB model was adopted with the modified forcing scheme by Li et al. Then the bubble collapse near rough solid wall was simulated. The bubble collapse mechanism was investigated from the dynamic process including bubble profiles evolution, pressure and velocity field evolution. The eroding effects of shock wave and jet were analyzed in details. It is found that the process and the effect of bubble collapse are affect seriously by the geometry of solid boundary. In the present work, the interaction between the collapsing bubble and the geometries boundary structured with periodic grooves was investigated. We found that the fluid velocity is the major cause for the eroding effect at the side walls of the grooves, and the pressure is the major cause for the eroding effect at the bulge part of the rough solid wall.
The improved forcing scheme proposed by Li et al. provides a convenient and efficient approach to achieve thermodynamic consistency. By the modelling of the bubble collapse near the rough solid wall, it is demonstrated that the MRT pseudopotential LB model improved by the modified forcing scheme has enough stability to describe the whole process of bubble collapse and the interaction between the collapsing bubble and complex geometry boundary.
The present work is supported by the National Key Research and Development Program of China (Grant No. 2016YFC0401600), the Primary Research & Development Plan of Jiangsu Province, China (Grant No. BE2016056) and the Natural Science Foundation of Jiangsu Province, China (Grant No. SBK2014043338).