Investigation of Pressure Oscillation and Cavitation Characteristics for Submerged Self-Resonating Waterjet

: The cavitation phenomenon of the self-resonating waterjet for the modulation of erosion characteristics is investigated in this paper. A three-dimensional computational ﬂuid dynamics (CFD) model was developed to analyze the unsteady characteristics of the self-resonating jet. The numerical model employs the mixture two-phase model, coupling the realizable turbulence model and Schnerr–Sauer cavitation model. Collected data from experimental tests were used to validate the model. Results of numerical simulations and experimental data frequency bands obtained by the Fast Fourier transform (FFT) method were in very good agreement. For better understanding the physical phenomena, the velocity, the pressure distributions, and the cavitation characteristics were investigated. The obtained results show that the sudden change of the ﬂow velocity at the outlet of the nozzle leads to the forms of the low-pressure zone. When the pressure at the low-pressure zone is lower than the vapor pressure, the cavitation occurs. The ﬂow ﬁeld structure of the waterjet can be directly perceived through simulation, which can provide theoretical support for realizing the modulation of the erosion characteristics, optimizing nozzle structure.


Introduction
Submerged self-resonating waterjet is a simple means with strong self-oscillating cavitation characteristics [1,2], which have been widely applied in the industries of cutting of solid materials, cleaning and burry removal of complicated mechanical products [3][4][5].
Over the past few years, deep-sea resources have attracted worldwide attention, and selfresonating waterjet also has broad applications in the exploitation of deep-sea mineral resources. The deep-sea environment is so complex, and the mineral resources have a changeable outlook when facing different erosion materials and cavitation behaviors, the confining pressure and the velocity of the flow decision on the oscillation characteristics of self-resonating waterjet. Although the modulation mechanisms of oscillation characteristics are a key technology in a complex deep-sea environment, it still not perfect.
Thus far, many experimental studies on submerged self-resonating waterjet have been done with the driving pump pressure [6], structure and size of the nozzle [7], cavitation erosion [8], etc., which, evaluated by the losing mass of the specimen, have the disadvantage of long erosion time [9]. Additionally, by contrast images, one or two cycles from the highspeed visualization in the article [10], the research of the cavitation cloud is defective for the influence of the unsteady pressure.
Moreover, in recent years, many researchers have focused on the mathematical modeling of oscillation characteristics of submerged self-resonating waterjets. Belyakov [11] studied the critical cavitation number and the vortex shedding by various liquid pressures and several jet exit velocities under the submerged conditions. They suggested the empirical formula of the force ratio at different cavitation numbers. Alehossein [12] established a numerical finite-difference solution for simulating the process of cavitation bubbles form, collapse, rebound and travel along the jet. They combined the computational fluid dynamics results and the modified Euler and Runge-Kutta-Fehlber, and successfully solved more wide range results of Rayleigh-Plesset formula solutions. In this research [13], we had pointed out that the resonance frequency of the nozzle played an important role in modulating the oscillation characteristics. According to the "triad resonance" theory, the frequencies bands of the self-resonating waterjet had acquired.
Furthermore, with the development of CFD, numerical studies for submerged selfresonating waterjet have seen flourishing development, for prediction cavitation erosion [14], bubble dynamics [15], aided design [16], etc. However, the systematically specialized research on submerged self-resonating waterjet velocity, pressure, and cavitation distribution characteristics is not enough.
For the investigation of the modulation mechanisms of the submerged self-resonating waterjet, the present study has started from the resonance frequency of the self-oscillating jet, based on mixture multiphase models, coupling realizable turbulence model and Schnerr-Sauer cavitation model, a two-phase flow model of submerged self-resonating waterjet is developed. Define and adjust model parameters, cater to the "triad resonance" theorythe fundamental frequency, the nature frequency, the resonance frequency are almost equal [17], and comparison of frequency distribution characteristics of numerical simulations against experimental data of the pressure oscillation is done. Based on this, the structural characteristics of the flow field at the nozzle exit are systematically studied.

Multi-Fluid Model and Numerical Method
The self-resonating waterjet is a mixture of liquid and gas. In this manuscript, the velocity of the objective fluid is about 170 m/s, which is far away from the Mach number 0.3, so the compressibility is ignored. The mixture was modeled as incompressible flow with the assumption that both phases are incompressible, and the no-slip condition is assumed between them.

Governing Equations
The self-resonating waterjet mixture flow with two separate phases, which are determined as a linear function of the vapor volume fraction and the gas and liquid properties [18][19][20]. The mixture flow properties of the mixture density ρ m is defined as follows, where ρ l is the liquid density, α v is the vapor volume fraction, and ρ v is the vapor density. The governing equations consist of the continuity and momentum equations, and the flow field can be modeled by Navier-Stokes equations expressed as (2)-(5) [21][22][23], where t, g, V, P, µ m , and F are time, gravitational acceleration, velocity vector, static pressure, viscosity, and body force vector, respectively.

Turbulence Model
In engineering applications, RKE is one of the most popular models with the advantages of high computational efficiency and simplicity [24,25]. With reference to [26][27][28], the researchers adopted a new dissipation rate model and realizable eddy viscosity formula, and RKE improved the standard k − ε model. The RKE includes two equations, one is a turbulent kinetic energy equation written as Equation (5), the other is a turbulent kinetic energy dissipation rate equation, Equation (6), in the simulation the velocity was use on the water in Equations (5) and (6).
where k represents turbulent kinetic energy as Equation (7), and ε represents turbulent kinetic energy dissipation rate as Equation (8), In the equation, turbulent viscosity µ t is Equation (9), and the turbulence kinetic energy is p k . The research focuses on the submerged self-resonating waterjet, without losing generality, so that we choose the turbulence parameters provided in the model, and the values of the constants were set as C ε1 = 1.9, C ε2 = 1.9, ε = 0.75,

Schnerr-Sauer Cavitation Model (SS)
Based on the Rayleigh-Plesset equation, the Schnerr-Sauer cavitation model is the function of the radius of the bubble and the number of the bubbles in the unit volume [29][30][31], the mass transfer between liquid and vapor was described by the source terms expressed as follows, → V is the fluid velocity, the source terms R p − R c is the mass transfer between the vapor and liquid phases in cavitation, where 3 is the radius of the bubbles, which was assumed to be a function of the vapor volume fraction, n b = 1 × 10 13 m 3 is the bubble number density.

Numerical Setup
The schematic of the computational domain used was shown in Figure 1a, using axisymmetric simulations. As the size of the nozzle was small, we simplified the fluid field to the maximum diameter of 160 mm, maximum height of 80 mm, and the stand-off distance is 6 mm, which can fully reflect the pressure oscillation and cavitation characteristics at the outlet of the waterjet nozzle. The boundary conditions were shown in Figure 1b, the inlet and outlet boundaries were defined as pressure, the adjacent side of the wall was defined as fluid, and the fluid properties were shown in Table 1. We employed the HyperMesh commercials software [32] to acquire the structured grid of the simulation domain, which was shown in Figure 2. The grid was refined at the outlet of the nozzle, the number of the grids elements was 14,465,000. During the simulation the time step was defined as 5 × 10 −6 s. The computer we used was with two Intel(R) CPU E5-2683 v3, which had a main frequency of 2 GHz, and the number processors was defined as 28. Symmetric boundary conditions were used, half of the computational domain was considered for the numerical simulation. The SIMPLE algorithm for unsteady simulations was employed to solve the coupling between pressure and velocity [33]. At the nozzle exit, the fluid oscillated violently, so that the transient term was discretized using a first-order upwind method in the governing equations to stabilize the numerical simulation. All the governing equations could be found in the software users' manual ANSYS/Fluent (2017).
is the radius of the bubbles, which was assumed to be a function of the vapor volume fraction, 13 3 b n =1 10 m × is the bubble number density.

Numerical Setup
The schematic of the computational domain used was shown in Figure 1a, using axisymmetric simulations. As the size of the nozzle was small, we simplified the fluid field to the maximum diameter of 160 mm, maximum height of 80 mm, and the stand-off distance is 6 mm, which can fully reflect the pressure oscillation and cavitation characteristics at the outlet of the waterjet nozzle. The boundary conditions were shown in Figure 1b, the inlet and outlet boundaries were defined as pressure, the adjacent side of the wall was defined as fluid, and the fluid properties were shown in Table 1. We employed the Hy-perMesh commercials software [32] to acquire the structured grid of the simulation domain, which was shown in Figure 2. The grid was refined at the outlet of the nozzle, the number of the grids elements was 14,465,000. During the simulation the time step was defined as 5 × 10 −6 s. The computer we used was with two Intel(R) CPU E5-2683 v3, which had a main frequency of 2 GHz, and the number processors was defined as 28. Symmetric boundary conditions were used, half of the computational domain was considered for the numerical simulation. The SIMPLE algorithm for unsteady simulations was employed to solve the coupling between pressure and velocity [33]. At the nozzle exit, the fluid oscillated violently, so that the transient term was discretized using a first-order upwind method in the governing equations to stabilize the numerical simulation. All the governing equations could be found in the software users' manual ANSYS/Fluent (2017).

Experimental Setup
The experimental setup was shown in Figure 3, which mainly consists of a water tank, a plunger pump, which supplies a maximum operating pressure of 60 MPa and a maximum flow rate of 60 L/min, a flowmeter, a pressure gauge, a relief valve, a high-pressure tank, a nozzle, and a hit plate were installed inside the high-pressure tank. The Piezotronics ICP model 102B03 pressure sensor, with a resonant frequency of 500 kHz, has a non-linearity of less than 1%FS. In our experiment, we used two sensors, one was arranged on the pipeline outside the high-pressure tank to record the inlet pressure P in , while the other was placed at the drainage of the high-pressure tank to record the confining pressure P a . An RHS-20 hydrophone with the flat response 120 kHz and ±3 dB, with oscillations up to 200 kHz, as well as a LMS SCADAS multi-functional data acquisition system.  Self-resonating nozzle geometry schematic description is shown in Figure 4, the structure parameters are shown in Table 2. The nozzle consists of three parts, an upstream domain contraction (Ds/D), a downstream domain contraction (D/d1), and a resonant chamber with a length of L with a diameter of D.   In the experiment, the temperature of water is about 10 • C, the mean value of inlet pressure was P in = 15.51 MPa, the mean value of measured confining pressure was P a = 1.503 MPa, the mean value of measured flow rate of the pump was Q = 19.98 L/min, the distance between the hit plate and the nozzle was 6 mm, and the sampling frequency was f a = 204,800 Hz.
Self-resonating nozzle geometry schematic description is shown in Figure 4, the structure parameters are shown in Table 2. The nozzle consists of three parts, an upstream domain contraction (Ds/D), a downstream domain contraction (D/d1), and a resonant chamber with a length of L with a diameter of D.  Self-resonating nozzle geometry schematic description is shown in Figure 4, the structure parameters are shown in Table 2. The nozzle consists of three parts, an upstream domain contraction (Ds/D), a downstream domain contraction (D/d1), and a resonant chamber with a length of L with a diameter of D.    For the sudden changes of flow velocity at the downstream domain contraction, pressure waves will be formed when high-speed flow passes through the nozzle. As a result, fluid resonance formed at the upstream and downstream contractions. According to the fluid resonance, the peak resonance will be achieved when the self-excited frequency of the water jet f s was matched with the acoustic natural frequency f n of the nozzle [13]. The fundamental frequency of the water jet structuring f s corresponds with the Strouhal number S d expressed as (13), and the average flow is defined as (14). The acoustic natural frequency of the chamber f n is determined by the length of the chamber L, which is expressed as (15). When (D S /D) > 1, the "mode parameter" K N is given by (16).
where d1 is the nozzle orifice diameter, V out is the jet velocity, m/s, N is mode number of the organ-pipe chamber of the nozzle, and the sound wave velocity c = 1450 m/s is given by the data from [34,35]. Joint Formulas (13)-(17), combined parameters in Table 2 and experiment parameters, and the acoustic natural frequency of the nozzle is about 14 kHz.

Validation
The CFD simulation model which we built up was validated by comparing with the experimental data, and applying the Fast Fourier Transform (FFT) analysis method on the collected data in order to identify main frequencies. Figure 5 is the time domain and spectrum of frequencies obtained through the FFT method of the experiment pressure signal, which was taken from the inlet pressure P in , which was shown in Figure 1. The frequencies contain broadband spectrum with a series of peaks, f e (1) = 2894 Hz ≈ 2.89 kHz, the sub-harmonic frequency is f e (2) = 5992 Hz ≈ 5.99 kHz, the fundamental frequency is f e (3) = 1.196 × 10 4 Hz ≈ 11.96 kHz, and the super-harmonic frequency is f e (4) = 1.788 × 10 4 Hz ≈ 17.88 kHz. The most important external factor affecting the structures of the waterjet is the hydrodynamic characteristics at the outlet of the nozzle. For this reason, during the simulation, the pressure oscillation characteristic at different locations of the nozzle was monitored. The location of the monitoring points in the simulation were shown in Figure 6, case 1 was any point selected from the domain 1, case 2 was any point selected from the domain 2, and case 3 was any point selected from the domain 3. Among them, case 1 and case 2 were used for monitoring cavitation frequencies, and case 3 was used for monitoring the acoustic natural frequency.
Appl. Sci. 2021, 11, x FOR PEER REVIEW The most important external factor affecting the structures of the waterje drodynamic characteristics at the outlet of the nozzle. For this reason, during tion, the pressure oscillation characteristic at different locations of the nozzle tored. The location of the monitoring points in the simulation were shown i case 1 was any point selected from the domain 1, case 2 was any point selecte domain 2, and case 3 was any point selected from the domain 3. Among them, case 2 were used for monitoring cavitation frequencies, and case 3 was used fo ing the acoustic natural frequency.  The most important external factor affecting the structures of the waterj drodynamic characteristics at the outlet of the nozzle. For this reason, during tion, the pressure oscillation characteristic at different locations of the nozzle tored. The location of the monitoring points in the simulation were shown case 1 was any point selected from the domain 1, case 2 was any point select domain 2, and case 3 was any point selected from the domain 3. Among them case 2 were used for monitoring cavitation frequencies, and case 3 was used f ing the acoustic natural frequency.  Additionally, through the FFT method, the simulation frequencies containing broadband spectrums with a series of peaks is shown in Figure 7 As we know, the self-resonating water jet has the characteristics of the pulse and cavitation, which is a typical nonstationary, time-varying signal [36]. The cavitation also has the modulating function on the mechanism [37]. The modulation signals model can be expressed as Equation (18)   As we know, the self-resonating water jet has the characteristics of the pulse and cavitation, which is a typical nonstationary, time-varying signal [36]. The cavitation also has the modulating function on the mechanism [37]. The modulation signals model can be expressed as Equation (18) [38], the signals of different positions at the outlet of the nozzle can be considered as a modulation model as Equation (18), the frequencies bands with the experimental and simulation are compared in Table 3, and the results are in a good agreement. Error evaluation was provided according to different frequency bands of the experimental and simulation results. Define the absolute error E (a) (i) as Equation (19), where f s (i) are the simulation frequencies and f e (i) are experimental frequencies, then we can take the absolute error at different frequencies bands from the results that are shown in Table 3. The absolute error E (a) (i) ranges from 1.3% (very good) to 17.2% (satisfactory). For the reason that the acoustic natural frequency of the nozzle is about 14 kHz, E (a) (4) = 16.7% is the most suitable error for evaluate the simulation model.
Based on the above verification, the result had been shown that the numerical simulation adequately predicts the behavior of the self-resonating waterjet.

Simulations Results and Discussions
The characteristics and cavitation of the flow field have a significant influence on the impact erosions of self-resonating waterjets. However, measurement of velocity distributions, pressure distributions, and cavitation intensity in the flow field was hard to implement.

Velocity Distribution
For the sudden change of the nozzle geometry as shown in Figure 4, the high-speed flow passes through the nozzle, which would lead to the rapid fluid changes, and forms velocity vortexes in the direction of the flow. Figure 8 is the vector field of the velocity distribution at the outlet of the nozzle, to the knowledge that it consists of the first and the second shear zone. Figure 8a illustrates the first shear zone, the interaction between the shear layer of the jet and the inner surface of the nozzle exciting the jet shear layer and forms the high energy velocity vortices at the throat of the nozzle; Figure 8b illustrates the second shear zone, when the high-speed jet interacts with the relatively stationary liquid, the formation of the vortexes is caused. Figure 9 is the streamline of the velocity distributions, where the vortexes can be seen clearly in the first shear zone, and much bigger vortexes in the second shear zone, and the stationary vortex adjacent to the hit plate. second shear zone. Figure 8a illustrates the first shear zone, the interaction between the shear layer of the jet and the inner surface of the nozzle exciting the jet shear layer and forms the high energy velocity vortices at the throat of the nozzle; Figure 8b illustrates the second shear zone, when the high-speed jet interacts with the relatively stationary liquid, the formation of the vortexes is caused. Figure 9 is the streamline of the velocity distributions, where the vortexes can be seen clearly in the first shear zone, and much bigger vortexes in the second shear zone, and the stationary vortex adjacent to the hit plate.

Pressure Oscillation Characteristics
For the sudden change of the flow velocity, pressure fluctuation would be formed at the downstream area contraction. In the center of the velocity vortices, low-pressure domains were created, the distributions of the pressure are shown in Figure 10, a very wide range of pressures is −1.371 to 16.85 MPa. From the local zoom, it can be seen that the local low-pressure areas are periodically shedding from the throat of the nozzle, and moving down the axis. The dynamic process of pressure could be observed, which means that the distribution of the bubbles along the jet path is responsible for the cavitation. However, whether strong cavitation occurs or not depends on the relations between the local pressure and the vapor pressure.

Pressure Oscillation Characteristics
For the sudden change of the flow velocity, pressure fluctuation would be formed at the downstream area contraction. In the center of the velocity vortices, low-pressure domains were created, the distributions of the pressure are shown in Figure 10, a very wide range of pressures is −1.371 to 16.85 MPa. From the local zoom, it can be seen that the local low-pressure areas are periodically shedding from the throat of the nozzle, and moving down the axis. The dynamic process of pressure could be observed, which means that the distribution of the bubbles along the jet path is responsible for the cavitation. However, whether strong cavitation occurs or not depends on the relations between the local pressure and the vapor pressure. mains were created, the distributions of the pressure are shown in Figure 10, a very wide range of pressures is −1.371 to 16.85 MPa. From the local zoom, it can be seen that the local low-pressure areas are periodically shedding from the throat of the nozzle, and moving down the axis. The dynamic process of pressure could be observed, which means that the distribution of the bubbles along the jet path is responsible for the cavitation. However, whether strong cavitation occurs or not depends on the relations between the local pressure and the vapor pressure.

Cavitation Characteristics
When the local pressure reduced down to the vapor pressure, cavitation generally occurs, consequently, the gas-filled or gas and vapor-filled bubbles are formed [39,40]. Usually, jet-induced cavitation occurs in the jet shear layers primarily, at the same time the cores of vortex are present in the jet flow, and the cavitation number is usually defined as Equation (20) for submerged waterjets [5,41]. The occurrence of cavitation, following the difference between the inlet pressure, the confining pressure and the vapor pressure of the fluid, had been assessed through the cavitation number, the vapor pressure set as P c = 3169 Pa in the simulation, and combining with experimental parameters. In this study, the cavitation number is σ = 0.085, The vapor pressure was defined as P c = 3169. Pa was an absolute pressure, so that we need to change it to gauge vapor pressure as P g = P c −atmospheric pressure = (3169 − 1.013 × 10 5 ) Pa =−97131 Pa. As mentioned before, the condition of cavitation is that the simulation flow field pressure is lower than the gauge vapor pressure. The errors ranges reult was shown in Figure 10, the minimum pressure of the simulation flow field was distribution at the throat of the nozzle, the value is −1.371 × 10 6 , which is lower than the gauge vapor pressure P g , and it can be concluded that strong cavitation occurs at the throat of the nozzle. However, in other places of the simulation domain, the dynamic process of local pressures m did not observe obvious cavitation phenomenon, for the reason that the pressure was higher than the gauge vapor pressure P g . The volume fraction of the vapor phase, as shown in Figure 11, ranges from 0 to 0.676. It can be seen that the strongest cavitation occurs at the start point of the straight section of the nozzle throat, which is clearly shown in the local zoom of Figure 11. As we know, erosion characteristics have a strong correlation with the cavitation intensity, and we can through the cavitation characteristics to realize the control and modulation of the erosion characteristics. tion phenomenon, for the reason that the pressure was higher than the gauge vapor pressure g P . The volume fraction of the vapor phase, as shown in Figure 11, ranges from 0 to 0.676. It can be seen that the strongest cavitation occurs at the start point of the straight section of the nozzle throat, which is clearly shown in the local zoom of Figure 11. As we know, erosion characteristics have a strong correlation with the cavitation intensity, and we can through the cavitation characteristics to realize the control and modulation of the erosion characteristics. Figure 11. Vapor phase volume fraction contours at the outlet of the nozzle. Figure 11. Vapor phase volume fraction contours at the outlet of the nozzle.

Conclusions
In this article, the combination of the numerical and experimental methods is employed to investigate the cavitation and erosion of submerged self-resonating waterjet. The multiphase, turbulence and cavitation models were employed, and we adopted the simple algorithm method to research the fluid domain characteristics of the waterjet. The performance of the selected simulation models was evaluated by the experiments data, and E (a) (4) = 16.7% is the most suitable error for evaluating the simulation modal, which had good agreement in the frequencies bands. Important conclusions are briefly summarized as follows: (1) The sudden change of the flow velocity for the geometry of the nozzle formation the velocity vortexes consists of two shear zones. The first shear zone is the formation by the interaction between the shear layer of the jet and the inner surface of the nozzle; the second shear zone is because of the high-speed jet inflow into the relatively stationary liquid, and the strength of the vortex is much lower than the first shear zone. (2) Velocity vector backwash produces a low-pressure region, and at the strongest position of our research, the negative pressures inside the conduits can drop down to −1.3 MPa. At the same time, the distribution of the bubbles along with the jet shedding from the throat of the nozzle periodic, the cavitation depends on the relation between the pressure and the vapor pressure. (3) The results show that the vapor phase at the throat of the nozzle has the highest of the volume fraction, where the strongest cavitation in our research occurs. This shows that when the velocity vortexes appear and the pressure is lower than vapor pressure, the cavitation occurs.
The results indicate the capability of the numerical model in predicting the cavitation characteristics of the waterjet, and erosion characteristics had a strong correlation with the cavitation intensity, then we can realize the control and modulation of the erosion characteristics.
For future studies, we will do more research on different working conditions to find out a better case in high amplitude when the frequency bands are similar.