Simulation of cavitation induced by water hammer ＊

Jie Geng (耿介), Xiu-le Yuan (苑修樂), Dong Li (李冬), Guang-sheng Du (杜廣生)School of Energy and Power Engineering, Shandong University, Jinan 250061, China,E-mail: gj_8944@163.com

Introduction

The cavitation might happen when the local pressure is lower than the working fluid’s saturated vapor pressure, which is a discontinuous process. The cavitation widely exists in nature, such as, behind blades of the rotor propeller, under the waterfall and on the surface of rocks in high-speed rivers[1,2].

In the hydraulic system, the cavitation generally consists of four processes, the gas nuclear generation,the expansion, the compression, and the collapse process. The noise and the vibration, caused by the pressure pulsation with bubble collapsing, are harmful in most cases. The cavitation will seriously affect the reliability and durability[1]. Therefore, the prevention of the cavitation is an important issue in the industrial designs. The water hammer can create a pressure transient in the pipeline, which may lead to the cavitation[3]. The phenomenon is often induced by an instantaneous closure of the valve. Nowadays, with the increasing demands for the precise flow control, a high speed valve is very common. Consequently, the water hammer prevention becomes more and more important. Since the nineteenth century when Russian scientist Joukowski (1898) firstly introduced the theory of water hammer, the water hammer has been a research focus for more than a hundred years. The main numerical methods include: the MOC, the wave front method, the finite difference (FD) method, and the finite volume method (FVM). The MOC is the most popular technique for solving partial differential equations for the water hammer, with a high computational efficiency, and a programming simplicity.Acrivos[4]developed the method of characteristic(MOC) for heat and mass transfer problems.Tijsseling and Bergant[5]removed the need of grid development by introducing the meshless method. An Implicit Method of Characteristics is proposed by Afshar and Rohani[6]and any arbitrary combination of devices can be allowed with implicit method. Zhao and Ghidaoui[7]compared the application of FVM schemes and MOC schemes with space line interpolation and found that MOC produces the same results with first-order FV Godunov-scheme. Meniconi et al.[8]studied the interaction between valve action and an in-line devices using MOC. The predictions with the MOC are relatively reliable, widely used in energy plants, environment industry, agriculture automation,chemical industry, urban water supply and other fields.

The water hammer is a complex process, and the simulation is mostly limited to one-dimensional cases[9].The industrial design is restricted by an unreal dimension. With the progress of computer science and hardware, the computational cost is greatly reduced.Many approaches of three-dimensional cavitation simulation with the FVM were proposed in recent years. A general cavitation prediction method was developed by Singhal et al.[10]and introduced to improve the performance of turbo-machinery. Bakir et al.[11]investigated the cavitating behavior of an Inducer. Shang et al.[12]studied the cavitation around a high speed submarine. To the best of the authors’knowledge, the related references mainly focused on the cases of pumps and high speed airfoils[13-18], the cavitation induced by the water hammer has not yet been studied by using three-dimensional numerical methods.

In the following, joint simulation by one-dimensional water hammer software “Hammer” and open source FVM software “OpenFOAM” is introduced.Firstly, pressure changes on the boundary will be calculated by MOC with Hammer. Then, the threedimensional prediction of cavitation will be achieved and three-dimensional visualization will be realized with the help of FVM. By analysis of the causes of cavitation from a three-dimensional perspective, it can provide new ideas for the industrial design of the fluid equipment.

1. Model description

1.1 Equations of water hammer

The propagation velocity of the wave that would induce the water hammer[4]

wherevE is the bulk modulus of liquid, E is the elastic modulus of tube’s wall, D is the pipe diameter, and T is the wall thickness.

The basic differential equations of the water hammer[19,20]:

where H is the static pressure head, α is the inclined angle of the pipeline, λ is the resistance coefficient and s is the axial length.

1.2 Equations of turbulence and cavitation

The --kωSST model[21]is used to model the turbulence, which enjoys a robust convergence rate and a relatively accurate inverse pressure gradient and separated vortex. The cavitation model is based on the volume of fluid model (VOF), the interPhaseChangeFoam in OpenFOAM.

The k and ω equations in the - -kω SST model:

The phase transition is described by the Merkle mass transformation equation:

where m˙-is the transition rate from the liquid phase to the gas phase, and m ˙+is the transition rate from the gas phase to the liquid phase, Cc, Cv, t∞, U∞are constants depending on the mean velocity, and psatis the saturated vapor pressure.

1.3 Boundary condition and grid settings

The physical model, initially designed by Nicolaus Bernoulli, is a flow measurement system operated in the Fluid Dynamics Laboratory of Shandong University[13,22,23], as shown in Fig.1.

Fig.1 Physical model

Table 1 The material properties

High speed valves are used to ensure the accuracy of the flow control (the valve closing process involves nonlinear changes of the flow rate. Improving Note: d is the pipe diameter of the simulation, l is the length between the water tank and the valve, h stands for the pressure head, tΔ is the closing time for the valve,wavec is the velocity of the water hammer wave, andpipeλ is the resistance loss of the pipe.

the speed of the valve can reduce the impact of the nonlinear part of the measurement). In the actual process, the test bench will witness undesirable vibrations. The water hammer effect is significant. On the other hand, the ultrasonic flowmeter is very sensitive to the bubbles when there is cavitation. The potential damage caused by the pressure pulsation with bubble collapsing should also be considered. The cavitation control is the key in the design. The basic parameters of the structure are shown in Table 1.

In the OpenFOAM, the x plane at =0 mx is set as the total pressure for the fluid inlet. The x plane at =0.2 mx is set as the total pressure for the fluid outlet. Other boundaries are the wall.

In normal cases, the cavitation would happen at the valve end of the pipeline. However, this paper would specifically study the cavitation when the water hammer wave interacts with the flow resistance piece.The fluid delivering devices often connect the parts with the rubber gasket to strengthen the sealing. In some cases, the rubber gasket will extend into the pipe and generate disturbance in the flow, as shown in Fig.1. The stability of the equipment will be influenced in some extreme situation, such as that with cavitation. The boundary condition of the VOF is introduced based on the result of the MOC, and the grid setting is shown in Fig.2.

Fig.2 Normal section mesh of pipeline

The calculation region consists of a circular tube and a rubber gasket. Because the accuracy in the boundary layer will dramatically affect the simulation result around the gasket, the first layer of the boundary grid is set to 3×10-6m, and y+is 5. The direct solution is used instead of wall functions. The total grid number is 1.7×106, as shown in Fig.2. The calculation tolerance for P is 1×10-5, tolerances for U and k are 1×10-4.

1.4 Model verification

The flow rate fluctuation calculated by the FVM is verified by the counterpart based on the MOC. The method of characteristics, which has been applied for the water hammer simulation for years, is developed and tested[5,6]. The cost of the computing resource is small, and the accuracy of the calculation for the extreme pressure is reliable. Therefore, the MOC for the water hammer calculation is regarded as the verification basis.

The boundary conditions of the FVM simulation are based on the pressure transient calculated by the MOC, and the flow rates of the same section downstream obtained from the two simulation methods, on a normal section 0.5 m before the valve, are compared.The simulation results are shown in Fig.3.

Fig.3 Flow rate fluctuation based on MOC and FVM

Due to the introduction of the no slip wall and the cavitation in the FVM, different scale vortexes are generated. The results are slightly different from the counterpart in the MOC simulation, and the correlation between the data from the two methods is 0.87. If the two curves are time averaged over 5 ms, they will be perfectly overlapped and the difference will be within 1%. The results of the FVM are relatively smoother. Overally, the trend is consistent. Therefore,the joint simulation is reliable.

2. Results and discussions

The simulation is carried out by using jointly the one-dimensional MOC and the three-dimensional FVM. Firstly, the pressure transient of the water hammer is calculated by the MOC. Then, the cavitation process is predicted by the FVM simulation.

Based on the boundary conditions as shown in Table 1, an indirect water hammer is induced. The pressure fluctuation curve calculated by the Hammer is shown in Fig.4, with the highest pressure of 2.75×105Pa at 0.106 s and the lowest pressure of 3.30×103Pa at 0.115 s. The cycle length of the water hammer pulsation at about 0.1s is equal to the ratio of the pipe length to the wave propagation velocity. In a normal situation, the temperature is approximately, the saturated pressure is 2.30×103Pa, which is lower than the lowest point in the simulation. Theoretically, there will be no cavitation. However, the structure of the fluid segments is also an important factor,which can be considered by the three-dimensional FVM.

Fig.4 Boundary pressure calculated by MOC ( 0.5 m beforethe valve)

The boundary conditions in the OpenFOAM are set according to the results of the Hammer’s simulation. The FVM segment is 0.202 m long and the start point is 2 m before the valve. The gasket is 2 mm wide and the extending section is 2 mm long. For the pressure boundary, the information on the related nodes of the MOC is used, which covers the whole pipe of 5 m long. The averaged pressure of the entire fluid region over time is shown in Fig.5, with the highest point at 0.105 s and the lowest point at 0.109 s.It is reasonable that the pressure difference is greatly reduced due to the space averaging. However, it is hard to explain why the lowest pressure comes a half cycle earlier than that on the boundary. At the same time, the cavitation is observed (In Fig.5, the solid line represents the volume percentage of the liquid).

Fig.5 Spatial averaged pressure fluctuation and volume fraction of water ( αwater )

Due to the cavitation, the pressure rebounds a little. At the same time, the pressure declines on the boundary. The cavitation is induced by the local extremely low pressure, and the generation of bubbles would compensate for the low pressure. The cavitation is an indicator of the extreme situation.

When the pressure is lower than the saturation vapor pressure, the cavitation is induced. On one hand,the pressure is constrained by the boundary conditions.On the other hand, the cavitation can be influenced by the pressure drop caused by the flow field structure.The pressure drop caused by the local resistant part is proportional to the velocity squared. The averaged space velocity over time is shown in Fig.6.

Fig.6 Spatial mean velocity

When t ＜ 0 .100 s, the average velocity is about 2.5 m/s, which is the maximum over the time span.

The cavitation appears when the inverse velocity reaches the maximum at 0.109 s. The flow fluctuates regularly and changes directions every half cycle, and the fluctuation energy is gradually dissipated. The interaction between the high-speed motion and the flow resistance is an important factor for the pressure drop. The transient flow field is used to study the localized situation.

Figure 7(a), a combination of the pressure field and the velocity vector field, illustrates the flow around the gasket. Initially, the flow is in the positive direction of the -xaxis. At the moment of 0.109 s,the flow reverses its direction and a strong vortex is formed at the front edge of the gasket, as shown in the figure. The local pressure declines remarkably around the core of the vortex. Bubbles are generated at the same time, as shown in Fig.7(b). In this practical case of the water hammer, where the minimum boundary pressure is higher than the saturated vapor pressure condition, the cavitation only happens at this specified time in this specified area.

Fig.7 Cavitation at t =0.109 s

The cavitation is localized in time and space. To a certain extent, the temporal analysis of the velocity and the pressure can explain the localization. The inverse velocity comes to the maximum and the boundary pressure is relatively low at 0.109 s. However,the cavitation usually happens at the back edge of the resistance and the absolute velocity before the water hammer is even higher than that at 0.109 s. Comparison is needed to illustrate the more significant pressure drop at 0.109 s.

A pressure and velocity vector field at =t 0.900s is shown in Fig.8(a). It is obvious that vortexes are generated behind the gasket and propagate along the flow direction, which is remarkably different from the standing vortexes in Fig.7(a). The squared vorticity shown in Fig.8(b) can be expressed as

Fig.8 Flow information at t=0.900 s

Fig.9 Squared vorticity at t =0.109 s

Before t = 0.100s, the flow is relatively stable.The squared vorticity in Fig.8(b) is time-averaged over 0.05 s.

The surface integral of the squared vorticity fluctuates at the downstream of the gasket at t = 0.900s . By contrast, the same parameter concentrates locally at t = 0.109s demonstrated in Fig.9(a), and the concentration synchronizes with the cavitation shown in Fig.7(b). As illustrated in Fig.9(b),the surface integral of the squared vorticity reaches the maximum around x = 0.09 m and remains low in the range of x＜ 0 .09 m (the flow is along the negative direction of the x-axis). Therefore, the vorticity can be used to characterize the possibility of the cavitation. The transfer equation of the vorticity is

The cycle of the pressure transient in our case is relatively short due to the short length of the operating pipeline. The effects of the convection, the diffusion and the viscous loss are negligible in this process.Consequently, the pressure gradient would be dramatically influenced by the vorticity changes.Therefore, we may conclude that the standing vortexes and the cavitation are locally generated in Fig.7(a).

3. Conclusions

The MOC has been widely adopted over decades,since it was developed. However, it is limited to onedimension, and the variables of internal structures in the flow field cannot be adequately considered. The simulation is carried out by using jointly the one dimensional MOC and the three-dimensional FVM.The three-dimensional FVM enjoys advantages of three-dimensional visualization. A practical case of the water hammer, with the minimum boundary pressure higher than the saturated vapor pressure condition, is simulated. Based on the temporal and spatial analyses, the possible cavitation area around a gasket is predicted, at the front edge of the gasket0.09 m and at the moment of t = 0.109s. The boundary pressure is relatively low, not reaching the minimum. The reverse speed is the highest0.109=v

2.0 m/s. Due to the short cycle of the pressure transient caused by the water hammer, vortexes cannot be developed downstream. When the standing vortexes are triggered, a remarkable pressure drop induces the cavitation.

In conclusion, there are three causes of the cavitation in pipes, the low boundary pressure, the strong vorticity, and the fast pressure transient. This theory provides some insight for further studies of the cavitation induced by the water hammer. Secondly,the joint simulation is a good way to realize the visualization, which could help the industrial design.Thirdly, in the case of a pipe of 5 m long with an indirect water hammer, the structure has to be reinforced at the front edge of the gasket to compensate for the potential harm due to the cavitation.

[1] Washio S. Recent developments in cavitation mechanisms:A guide for scientists and engineers [M]. Sawston, UK:Woodhead Publishing, 2014.

[2] Brennen C. E. Cavitation and bubble dynamics [M].Cambridge, UK: Cambridge University Press, 2013.

[3] Ghidaoui M. S., Zhao M., Mcinnis D. A. et al. A review of water hammer theory and practice [J]. Applied Mechanics Reviews, 2005, 58(1): 49-76.

[4] Acrivos A. Method of characteristics technique: Application to heat and mass transfer problems [J]. Industrial and Engineering Chemistry, 1956, 48(4): 703-710.

[5] Tijssenling A. S., Bergant A. Meshless computation of water hammer [C]. Proceedings of the 2nd IAHR International Meeting of the Workgroup on Cavitation and Dynamic Problems in Hydraulic Machinery and Systems.Timisoara, Romania, 2007, 65-76.

[6] Afshar M. H., Rohani M. Water hammer simulation by implicit method of characteristic [J]. International Journal of Pressure Vessels and Piping, 2008, 85(12): 851-859.

[7] Zhao M., Ghidaoui M. S. Godunov-type solutions for water hammer flows [J]. Journal of Hydraulic Engineering, ASCE, 2004, 130(4): 341-348.

[8] Meniconi S., Brunone B., Ferrante M. In-line pipe device checking by short-period analysis of transient tests [J].Journal of Hydraulic Engineering, ASCE, 2010, 137(7):713-722.

[9] Zeng Q., Ma G. Y., Jiang D. F. et al. A review of calculation methods for pipeline water hammer [J]. Contemporary Chemical Industry, 2013, 42(8): 1189-1193.

[10] Singhal A. K., Athavale M., Li H. et al. Mathematical basis and validation of the full cavitation model [J].Journal of Fluids Engineering, 2002, 124(3): 617-624.

[11] Bakir F., Rey R., Gerber A. G. et al. Numerical and experimental investigations of the cavitating behavior of an inducer [J]. International Journal of Rotating Machinery, 2004, 10(1): 15-25.

[12] Shang Z., Emerson D. R., Gu X. Numerical investigation of cavitation around a high speed submarine using OpenFOAM with LES [J]. International Journal of Computational Methods, 2012, 9(9): 1-14.

[13] Zheng Y. Research on current rush to interception airmass with delivery pipeline system and water hammer with gas [D]. Doctoral Thesis, Nanjing, China: Hohai University, 2004(in Chinese).

[14] Li J., Chen H. S. Numerical simulation of micro bubble collapse near solid wall in Fluent environment [J].Tribology, 2008, 28(4): 311-5.

[15] Ding H., Visser F., Jiang Y. et al. Demonstration and validation of a 3D CFD simulation tool predicting pump performance and cavitation for industrial applications [J].Journal of Fluids Engineering, 2011, 133(1): 011101.

[16] Luo X., Jin B., Tsujimoto Y. A review of cavitation in hydraulic machinery [J]. Journal of Hydrodynamics, 2016,28(3): 335-358.

[17] Zhang H. M. Studies on numerical simulation of weakly compressible coupling water hammer and cavitation flow in hydroturbine [D]. Doctoral Thesis, Kunming, China:Kunming University of Science and Technology, 2013(in Chinese).

[18] Li W. F., Feng J. J., Luo X. J. et al. Numerical simulation of transient flows in Francis turbine runner based on dynamic mesh method [J]. Journal of Hydroelectric Engineering, 2015, 34(7): 64-73.

[19] Bergant A., Simpson A. R., Tijssenling A. S. Water hammer with column separation: A historical review [J].Journal of Fluids and Structures, 2006, 22(2): 135-171.

[20] Mambretti S. Water hammer simulations [M]. Southampton, UK: Wit Press, 2013.

[21] Menter F. R. Improved two-equation -kω turbulence models for aerodynamic flows [R]. NASA STI/Recon Technical Report N, 1992, 93.

[22] GUO L. L. Research on the key questions of the influence on the performance of a water flowrate calibration facility[D]. Doctoral Thesis, Shandong, China: Shandong University, 2014(in Chinese).

[23] Guo L. L., Geng J., Shi S. et al. Computing research of the water hammer pressure in the process of the variable speed closure of valve based on UDF method [J]. Journal of Shandong University (Natural Edition), 2014, 49(3):27-30(in Chinese).