Abstract

The two-dimensional (2D) axial-symmetric model is applied to investigate the transient cavitating flows in the reservoir-pipeline-valve (RPV) system. Firstly, the MacCormack scheme is used to solve the governing equations, and compared to the numerical results of the one-dimensional (1D) model, the 2D head peaks and time-dependent evolutions predicted by five-region turbulence model in the frozen form are in better agreement with the experimental results, and the comparisons show that the maximum head relative errors of the 2D model are generally smaller than those of the 1D model. Then, further numerical simulations are carried out to investigate the performances of different turbulence models incorporated in the 2D model. The comparisons between the numerical results and the experimental ones show the head curves of the two-region turbulence model are similar to those of the five-region turbulence model, which indicates that the transient cavitating flows are insensitive to the magnitude and distribution of the eddy viscosity in the core region. In addition, the sensitivity to the quasi-steady form and the frozen form of the five-region turbulence model is implemented in the 2D model; the numerical results predicted by the two forms both agree better with the experimental results, as the non-dimensional parameter P increases.

INTRODUCTION

Water hammer causes not only pressure rise but also pressure drop in pipes. When the pressure drops below the vapor pressure of pure liquid, cavitation will occur. The occurrence, development and collapse of the cavitation bubble are processes of phase changes, which have the characteristics of complex physical mechanisms and difficult cognition. Therefore, it is important to analyze the problem of transient vaporous cavitation.

A lot of experiments and theoretical research have been carried out and many mathematical models of transient cavitating flows have been proposed in recent years. The Discrete Vapor Cavity Model (DVCM) with constant wave-speed proposed by Streeter (1983) assumed that the whole cross section of the pipe is full of vapor when the pressure of the section is lower than the vapor pressure and its pressure holds at vapor pressure until the collapse of the vapor. The evolution of the cavitation cavity is calculated by mass conservation of the liquid column. The model shows the basic characteristic of the transient cavitating flow and is easily programmed, so it is widely applied in the existing commercial water hammer software (Flowmaster, Hammer, etc., Marcinkiewicz et al. 2008). However, the DVCM may predict unrealistic pressure peaks in the case of collapse of multi-cavities (Simpson & Bergant 1995). The Discrete Gas Cavity Model (DGCM) (Simpson & Bergant 1994), which is similar to the DVCM model, has the hypothesis that there is always a bubble (i.e. gas nucleus) at the computational node and the size of the bubble changes with the change of ideal gas (isothermal or adiabatic). The calculated results predicted by the DGCM model are in good agreement with the experimental results and less affected by the input parameters than those of the DVCM model. The Generalized Interface Vaporous Cavitation Model (GIVCM) carried out by Simpson (1986) included water hammer equations for the liquid phase, two-phase equations for a distributed vaporous cavitation region, shock equations for condensation of the liquid-vapor mixture back to liquid, and equations for a discrete vapor cavity. Considering the two-phase vapor-liquid mixture as a single-phase mixture with average physical characteristics, the model accords better with the physical process but the algorithm is complex and time-consuming. In order to avoid the unrealistic pressure peaks, Adamkowski & Lewandowski (2012) proposed a new DVCM model, which assumes that the cavity volumes calculated in several computational nodes are shifted to the main cross section without breaking the conservation law. Taking the cavitation transport model into consideration, Sumam et al. (2010) adopted the MacCormack numerical scheme and a new pressure relaxation method to overcome the numerical instability caused in the phase transition process. The results predicted by the two-phase homogeneous equilibrium vaporous cavitation model with frequency dependent friction proposed by Shu (2003) were compared with those of the DVCM model, and the effectiveness of the new model was tested and verified. Great achievements have been obtained for transient cavitating flows, but there are significant differences in pressure amplitudes and phases (Szymkiewicz 2002) between the simulated and experimental results. This is related to the steady friction term or quasi-steady friction term used in these models, which cannot describe the complex characteristics of the shear stress and underestimate the energy dissipation observed. To solve such problems, a series of one-dimensional (1D) unsteady friction models are developed from the point of unsteady wall shear stresses. For a fast closure of the downstream-end valve in a simple reservoir-pipeline-valve (RPV) system, Urbanowicz et al. (2012) came to the conclusion that the combination of the 1D unsteady convolution-based (CB) (Zielke 1968) and DGCM model brought the results much closer to the experimental results. Bergant et al. (2005) introduced the Momentum Correction Factor (MCF) in the momentum equation and verified that the results of the DGCM-CB-MCF model were more consistent with the experiment. Storli & Nielsen (2011) found that the calculated results by using the Brunone unsteady friction model with a constant empirical coefficient k matched only at specific nodes and suggested that a position-dependent and time-dependent k should be employed. As another efficient method, the more accurate two-dimensional (2D) unsteady friction model which takes the velocity distribution of the cross section into consideration attracted more and more attention. Vardy & Hwang (1991) established the quasi 2D five-region model for both laminar and turbulent flows and verified the effectiveness of CB expressions in transient laminar flows. The numerical method is accurate and numerically stable, but it has the disadvantage of low efficiency in solving the sparse matrix in the coupling linear equation. For the same space and temporal discretization, the different numerical schemes proposed by Zhao & Ghidaoui (2003) and Korbar et al. (2014) had the same precision, but with less computation time. To reduce the computational cost, Ohmi et al. (1985) obtained the pressure and average velocity from the 1D continuity equation and the velocity distribution was derived from the 2D momentum equation. A similar procedure was used by Jang et al. (2016). A mixing length algebraic turbulence model was employed in smooth pipes (Silva-Araya 1997) and rough pipes (Silva-Araya & Chaudhry 2001). The two-region turbulence model proposed by Pezzinga (1999) was applied to simulate the transient non-cavitating flows in smooth and rough pipe networks. Zhao & Ghidaoui (2006) developed the quasi 2D model and pointed out that this model was not suitable for transient flow with high initial Reynolds number. Riasi et al. (2013) put forward the quasi 2D model and carried out the dimensionless analysis. The Baldwin-Lomax turbulence model, which was performed in both quasi-steady and frozen forms, was investigated by Wahba (2009) in transient non-cavitating flows. Pezzinga & Santoro (2017) studied two different cavitation models, both in 1D and 2D forms. Moreover, the 2D computational fluid dynamics (CFD) method was adopted by Wang et al. (2016) to analyze the vapor distribution, pressure and velocity fields in pipelines and the results demonstrated that the method could simulate the pressure variations and visualize the associated physical processes. Based on the full cavitation model, the commercial code CFD-ACE+ is used by Pinho et al. (2014) to investigate the corresponding characteristics of the transient cavitating events.

Based on the above facts, further research is carried out in the present work for transient cavitating flows based on the 2D model proposed by Pezzinga & Cannizzaro (2014). The comparisons between the 1D model and 2D model are related to the classical reservoir-pipe-valve system. The effectiveness and superiority of the 2D model is comprehensively tested and verified against a series of experimental data. In order to judge the reliability of the turbulence models and their assumptions in the transient cavitating 2D model, a further sensitivity investigation of different eddy viscosity models (five-region turbulence model and two-region turbulence model) on the accuracy is carried out. Moreover, both the quasi-steady and frozen forms for implementing the five-region turbulence model are examined in order to assess their range of applicability related to the non-dimensional parameter P.

MATHEMATICAL MODEL

If it is hypothesized that there is no slip between the liquid and the vapor in the transient cavitating flow, the quasi-2D model can be expressed by the continuity and momentum equations of the vapor-liquid mixture. And the 2D flow field with axial symmetry is assumed. Vardy & Hwang (1991) came to the conclusion that the maximum values of the radial velocities are of the order of 10–20 μm/s both in transient laminar and turbulent flow. The radial velocity component v and its derivatives can be neglected in the momentum equation and the associated governing equations in rectangular coordinates can be written in the following forms:  
formula
(1)
 
formula
(2)
where = density of liquid-vapor mixture; = total cross-sectional area of the pipe; = mean velocity; = time; = distance along pipe; = velocity components in the longitudinal direction; = gravitational acceleration; = pressure head; = distance from the axis; and = shear stress.
The density of the liquid-vapor mixture can be expressed by:  
formula
(3)
where = vapor void fraction, = liquid density; and = vapor density.
It can be noted that the vapor density can be neglected in comparison with the liquid density. When is very small, Equation (3) can be simplified as:  
formula
(4)
Substituting Equation (4) into Equation (1). Neglecting the convective term and with the product being a function of pressure p, the continuity equation can be rewritten as:  
formula
(5)
where a is the wave-speed for the liquid-vapor mixture, ; and p is the pressure.
Introduce the auxiliary variable and define it as:  
formula
(6)
The continuity equation of the liquid-vapor mixture can be written in the following form:  
formula
(7)
The calculation formulas of variables p and are as follows:  
formula
(8)
 
formula
(9)
where = vapor pressure. More details on the logarithmic cavitation model for the calculation of the pressure and void fraction can be found in Pezzinga & Cannizzaro (2014).
The mean velocity is calculated by the following formula:  
formula
(10)
where = radius of the pipe.
In the 1D model, the following 1D momentum equation is used instead of Equation (2), and is associated with Equation (7)  
formula
(11)
where = Darcy-Weisbach friction factor; and = diameter.

NUMERICAL SCHEME

For the 2D model, the pipe is discretized into cylinders with constant area in the radial direction and reaches with constant length in the longitudinal direction. The 1D model is a special case of the 2D model . The adopted time step is . For each discrete cylinder, axial velocity u is defined at each radial grid center; shear stress is defined on the inside and outside of the grid surface, as shown in Figure 1.

Figure 1

Discrete cylindrical grid of computational domain.

Figure 1

Discrete cylindrical grid of computational domain.

The MacCormack two-step scheme (Mirzaei-Shiri et al. 2017) is adopted. For the x-derivatives, forward finite differences are used in the predictor step:  
formula
(12)
 
formula
(13)
and backward finite differences in the corrector step:  
formula
(14)
 
formula
(15)

Indices ‘i’ and ‘j’ refer, respectively, to node numbers indicating longitudinal direction x and radial direction r; index ‘t’ refers to time t; indices ‘pre’ and ‘cor’ refer to the predictor and corrector steps. The values at time t are calculated as the average of the predictor and corrector values.

The steady state velocity profile should be used as the initial condition for the calculation of unsteady flow. For fully developed steady turbulent flow, the shear stress achieves its largest value at the pipe wall and decreases linearly to zero at the pipe centerline. Hence, the shear stress term can be expressed in the following form:  
formula
(16)
where is defined as the total kinematic viscosity of each cylinder cross section; and is the friction velocity. The algebraic expression of the eddy viscosity is evaluated from the five-region turbulence model proposed by Kita et al. (1980) for oscillating flows. An initial assumption for the value of is needed to start the computation and a shooting procedure for the friction factor is adopted to achieve the desired flow rate.

The numerical simulations are concerned with the fast valve closure in the reservoir-pipe-valve system. At the upstream point, the value of the pressure head is given while the local velocity is obtained by time-marching the momentum equation, and the x-derivative of the head is evaluated forward both in the predictor and the corrector step. At the downstream point, the pressure head during the valve closure is interpolated from a set of tabular experimental results at regular increments of time and evaluated by time-marching the continuity equation, and the x-derivative of the mean velocity is evaluated backward both in the predictor and the corrector step. When the valve is completely closed, the mean velocity is zero and the pressure head is evaluated backward both in the predictor and the corrector step. Moreover, a symmetry boundary condition is applied to the local velocity at the pipeline centerline and a no-slip condition is enforced at the pipe wall.

RESULTS AND DISCUSSION

In order to validate the proposed mathematical model and numerical scheme, the numerical results were compared with experimental tests carried out by Simpson (1986). Layout of the laboratory apparatus is shown in Figure 2. The laboratory apparatus is composed of a smooth copper straight sloping pipe (internal diameter 19.05 mm, length 36 m and thickness 1.6 mm) connecting two pressurized tanks. A specified pressure in each of the tanks is controlled by a computerized pressure control system. The ends of the pipe have an elevation difference of 1.0 m. Three piezoelectric transducers are located at the downstream valve, at the upstream 1/4 point (9.0 m from the upstream end), and at the downstream 1/4 point (9.0 m from the downstream valve), respectively. The wave speed is 1,280 m/s.

Figure 2

Schematic of laboratory apparatus.

Figure 2

Schematic of laboratory apparatus.

Nine different experimental tests were studied and a summary of all corresponding parameters considered is listed in Table 1. All simulated parameters are consistent with the experimental ones.

Table 1

Summary of corresponding parameters for different tests

  Initial velocity (m/s) Hydraulic grade line of the upstream tank (m) Friction factor Valve closure time (ms) 
R1 0.239 24.30 0.0325 31 
R2 0.332 23.41 0.0315 22 
R3 0.401 23.38 0.029 17 
R4 0.466 23.16 0.028 28 
R5 0.507 23.00 0.028 25 
R6 0.596 23.24 0.027 29 
R7 0.696 22.91 0.026 26 
R8 0.938 22.54 0.024 35 
R9 1.125 21.74 0.023 43 
  Initial velocity (m/s) Hydraulic grade line of the upstream tank (m) Friction factor Valve closure time (ms) 
R1 0.239 24.30 0.0325 31 
R2 0.332 23.41 0.0315 22 
R3 0.401 23.38 0.029 17 
R4 0.466 23.16 0.028 28 
R5 0.507 23.00 0.028 25 
R6 0.596 23.24 0.027 29 
R7 0.696 22.91 0.026 26 
R8 0.938 22.54 0.024 35 
R9 1.125 21.74 0.023 43 

Selection of discrete parameters

When several cylinders exist, there are many methods to determine the radial thickness. Bratland (1986) implied that the equal area can give more detail close to the pipe wall than the equal thickness . Since the radial disturbances are caused primarily by the no-slip condition at the wall and since diffusion processes are much slower than wave processes in low Mach number flows, the equal area has a significant advantage. So all results presented herein are obtained with cylinders of equal area. In addition to the selection of the thicknesses of various cylinders, the total number of them should also be selected. First of all, the grid independence verification of the 1D model is carried out for each test to determine Nx. Eight different numbers of pipe reaches (Nx = 8, 16, 32, 64, 128, 256, 512 and 1024, respectively) were considered and the parameter of the maximum head at the valve was used to assess the performance of the model. The results for R2 are shown in Figure 3(a). The maximum head at the valve increased gradually with the increase of Nx. When Nx > 32, the maximum head at the valve almost remained unchanged, therefore the optimal value of Nx is 32. Then, the grid independence verification of the 2D model is carried out to determine Nx and Nr. Four different numbers of discrete cylinders (Nr = 21, 30, 50 and 100, respectively) at a range of different Nx were considered and the parameter of the maximum head at the valve was used to assess to performance of the model. The results for R2 are shown in Figure 3(b). For each specific Nr value, the maximum head at the valve almost remained unchanged when Nx > 256. For the particular problem considered here, the use of significantly fewer cylinders Nr = 50 can not only save computer resources, but also guarantee the stability of the results. Therefore the value of Nx is 256 and the value of Nr is 50. So for the cavitating flow, the optimum ratio of the number of reaches (Nx) and the number of discrete cylinders (Nr) is found to be about 5.0.

Figure 3

Grid independent verification for R2, (a) maximum head at the valve for variable Nx (1D model), (b) maximum head at the valve for variable Nr (2D model).

Figure 3

Grid independent verification for R2, (a) maximum head at the valve for variable Nx (1D model), (b) maximum head at the valve for variable Nr (2D model).

Model verification

Preliminarily, the mathematical model and numerical scheme were tested for a non-cavitating test (R1), in which the transient pressures stay above the vapor pressure of the liquid. The pressure curves after the rapid closing of the valve at each pipeline monitoring point calculated by the 1D and 2D models are shown in Figure 4. It can be observed from the figures that there exists a small attenuation of the maximum head predicted by the 1D model. The first four peak pressures at the valve are 55.48 m, 55.15 m, 54.82 m and 54.50 m, respectively. And the first four peak pressures at the valve for the 2D model are 55.36 m, 54.33 m, 52.91 m and 51.46 m, respectively. Based on the above-mentioned quantitative analyses, the results predicted by the 2D model are closer to the experimental results because the velocity profile in the cross section is taken into account in the 2D model, which allows the wall shear stress to be more accurately evaluated. It can be noted that, also for upstream 1/4 point and downstream 1/4 point, the improvements made by the 2D model are similar to those previously shown.

Figure 4

Comparison between the computed heads by 2D model, 1D model and the experiment for R1 at different locations, (a) head at the valve, (b) head at the upstream 1/4 point, (c) head at the downstream 1/4 point.

Figure 4

Comparison between the computed heads by 2D model, 1D model and the experiment for R1 at different locations, (a) head at the valve, (b) head at the upstream 1/4 point, (c) head at the downstream 1/4 point.

The measured head oscillations compared with the results of the 2D model and 1D model are shown in Figure 5 for R2 at all monitoring points. The phenomenon of this test is characterized initially by a Joukowsky pressure rise caused by the rapid closure of the valve, and then the pressure reduced to vapor pressure and the subsequent pressure peak was caused by the collapse of the cavity. It can be noted that the 2D model is much better than the 1D model in reproducing the experimental results at all monitoring points. Compared to the experiment, there is a phase lag phenomenon of the pressure predicted by the 1D model and the discrepancies are magnified with the increase of time. The dissipation caused by friction could be responsible for the shaping and shifting of the pressure trace predicted by the 1D model. On the contrary, the 2D model reproduces the experimental oscillations with very good accuracy throughout the experiment. From the quantitative point of view, both models can predict that the pressure peak caused by the collapse of the first cavity at the valve, , is greater than the first peak pressure due to Joukowsky pressure rise. is 91.00 m in the experiment, 104.09 m predicted by the 2D model and 106.34 m by the 1D model.

Figure 5

Comparison between the computed heads by 2D model, 1D model and the experiment for R2 at different locations, (a) head at the valve, (b) head at the upstream 1/4 point, (c) head at the downstream 1/4 point.

Figure 5

Comparison between the computed heads by 2D model, 1D model and the experiment for R2 at different locations, (a) head at the valve, (b) head at the upstream 1/4 point, (c) head at the downstream 1/4 point.

An analogous comparison is shown in Figure 6 for R3. Also, in this case both the two models can reproduce well the above-mentioned phenomenon and predict that the pressure peak caused by the collapse of the first cavity at the valve, , is greater than the first peak pressure due to Joukowsky pressure rise. is 95.15 m in the experiment, 101.07 m predicted by the 2D model and 102.01 m by the 1D model. For the phase of the pressure trace, there is a good agreement between the calculated pressure curves of the 1D model and the experiment. But the 2D model has a tendency to anticipate the head oscillations in advance with respect to the experimental ones, which had also been observed and attributed to the neglecting of the convective terms and the release of dissolved gas in Pezzinga & Cannizzaro (2014).

Figure 6

Comparison between the computed heads by 2D model, 1D model and the experiment for R3 at different locations, (a) head at the valve, (b) head at the upstream 1/4 point, (c) head at the downstream 1/4 point.

Figure 6

Comparison between the computed heads by 2D model, 1D model and the experiment for R3 at different locations, (a) head at the valve, (b) head at the upstream 1/4 point, (c) head at the downstream 1/4 point.

Similar comparisons for the rest of the tests were also carried out and analyzed. All the comparisons show that the 2D model has significant improvements in both magnitude and phase shift of pressure head traces compared with the 1D model in most cases. Due to the different evaluation of dissipation, the 1D model allows a good reproduction of the instant of the pressure peaks only for the first and second oscillations, whereas the results of the 2D model are in good agreement with the experimental results in most cases. The oscillation phases computed by the 1D model are more consistent with the experimental ones, only for R3.

A summary of the comparisons between the results of the numerical models and the experimental data is shown in the bar diagrams of Figure 7, where the relative error (in percent) of the 1D model and the 2D model on maximum heads at the valve are reported. It can be observed clearly that the relative errors of the 2D model are generally smaller than those of the 1D model except for R7 and R8. The 2D model is always in better agreement with the experimental results than the 1D model because the 2D model takes the velocity profile into account. However, the reasons for the poor performance of the 2D model for R7 and R8 are worth studying in depth. Various simplifying hypotheses in the 2D model can affect the predicted results. Firstly, the dissipation caused by friction is not the only factor responsible for the shaping and shifting of the pressure trace in the transient flow. The effects of inertia caused by the unsteady velocity profiles, which is not considered, may be far greater than the wall shear stress effect on the pressure trace (Storli & Nielsen 2013). Then, according to Simpson (1986), the vapor fraction at the valve and the time of vapor existence increase with the increase of the initial velocity. The hypothesis of axial-symmetry is likely to lose efficacy in the presence of conspicuous vapor cavities, which would lead to the above-mentioned large relative errors. More detailed experimental results and more advanced models are needed to deepen these aspects.

Figure 7

Relative errors in percent between the computed maximum heads at the valve by 2D model and 1D model and experimental ones.

Figure 7

Relative errors in percent between the computed maximum heads at the valve by 2D model and 1D model and experimental ones.

There are two important timescales in transient flows, namely the timescale of the pressure wave and the timescale of the generated shear pulse due to passage of the pressure wave. Ghidaoui et al. (2002) introduced the dimensionless parameter P, which is the ratio of the two mentioned timescales and defined as:  
formula
(17)

A comprehensive study of the behavior of axial velocity profiles in the transient non-cavitating flow was numerically investigated by Riasi et al. (2009). A further investigation into the behavior of the velocity profiles in the transient cavitating flow was carried out in the paper. The predicted profiles of velocity divided by the initial mean velocity at the upstream 1/4 point at different instants by the 2D model for R2 (P = 22.99) and R9 (P = 7.94) are shown in Figure 8. When the positive pressure wave passes from the upstream 1/4 point at time t1, there is a uniform shift in the velocity profile in the core region and a strong gradient and a vortex ring (reverse flow) in the near-wall region because of the no-slip condition at the pipe wall, which produces a great wall shear stress. The generated vortex ring is subsequently diffused toward the pipe center and the reverse velocity moves to the core region. When the reflected positive pressure wave passes from the upstream 1/4 point at time t2, the whole cross section has a negative flow, which presents an S shape distribution and another vortex is generated at the pipe wall. The same condition occurs when the negative pressure wave comes from the valve at time t3 and also when the reflected positive wave passes from the upstream 1/4 point at time t4. The process of velocity shift, vortex generation and vortex diffusion will continue to occur as long as there are pressure waves in the pipe. It can be noted that the velocity profiles are different from the corresponding ones in the 1D condition with equal discharge and are very similar to those in the transient non-cavitating flow.

Figure 8

Axial velocity distributions at upstream 1/4 point computed by 2D model, (a) R2, (b) R9.

Figure 8

Axial velocity distributions at upstream 1/4 point computed by 2D model, (a) R2, (b) R9.

Five-region vs. two-region eddy viscosity distribution simulations

Turbulence models are needed to calculate the shear stress related to energy dissipation in the 2D model. Algebraic turbulence models, especially the two-region and the five-region eddy viscosity models, which relate the eddy viscosity to the mean velocity through some algebraic functions, have been widely used. The five-region model divides the flow cross-sectional region into a viscous sublayer, a buffer zone, and a core region. The buffer zone and core region are further sub-divided to reflect different viscosity characteristics. The two-region model divides the flow cross-sectional region into a viscous sublayer and a turbulent region. Ghidaoui et al. (2002) investigated the two-region and the five-region eddy viscosity models regarding the accuracy of the turbulent non-cavitating transient flows. And further systematic studies to evaluate the effects of these assumptions on cavitating transient flows were carried out in the paper.

Figure 9 depicts the computed head trace results obtained from the two-region and the five-region models and the measured head trace results for R2. Obviously, the pressure peaks as well as their shapes predicted by the two models almost coincide and they are in good agreement with the experimental results. The conclusions of the rest of the tests are consistent with those of R2, which shows that the turbulence models are insensitive to the form of the spatial distribution of the eddy viscosity within the core region. This is because the turbulent stresses are negligible compared to both the inertia and the pressure gradient in the core region. In future investigations, it is suggested that more efforts should be concentrated on the establishment of accurate models for the near-wall region.

Figure 9

Comparison between the computed heads by five-region and two-region turbulence model and the experiment for R2 at different locations, (a) head at the valve, (b) head at the upstream 1/4 point, (c) head at the downstream 1/4 point.

Figure 9

Comparison between the computed heads by five-region and two-region turbulence model and the experiment for R2 at different locations, (a) head at the valve, (b) head at the upstream 1/4 point, (c) head at the downstream 1/4 point.

Quasi-steady vs. frozen form simulations

The five-region turbulence model is developed from steady flows and assumptions are needed when applying the model to transient flows. For the quasi-steady form, the turbulence model is updated at each time step, providing values for the instantaneous eddy viscosity as a function of the instantaneous velocity profile. For the frozen form, the eddy viscosity values at any spatial location and any time step are kept frozen and equal to their initial steady state values at the same spatial location. The influence of the quasi-steady form on the accuracy of the transient non-cavitating flows was studied by Ghidaoui et al. (2002). In order to evaluate the influence of the quasi-steady and frozen forms regarding transient cavitating flows, a further intensive research of both forms on cavitating flows was carried out in the paper.

Numerical simulations for both forms of five-region turbulence model are performed first for R2 (P = 22.99). The predicted pressure plots at the monitoring points are compared with the corresponding experimental results in Figure 10. The results obtained by both forms are almost on top of each other and in good agreement with the experimental ones. This is because the shear pulse is mainly within the viscous sublayer and has no influence on the turbulence intensity and structure (He & Jackson 2000). That is to say, the flow characteristics are essentially identical to their steady state characteristics after the wave passes the particular location. Therefore, the frozen eddy form and the quasi-steady form are equally acceptable. Greenblatt & Moss (1999) obtained a similar conclusion for the temporally accelerating flow.

Figure 10

Comparison between the computed heads by five-region turbulence model (quasi-steady form and frozen form) and the experiment for R2 at different locations, (a) head at the valve, (b) head at the upstream 1/4 point, (c) head at the downstream 1/4 point.

Figure 10

Comparison between the computed heads by five-region turbulence model (quasi-steady form and frozen form) and the experiment for R2 at different locations, (a) head at the valve, (b) head at the upstream 1/4 point, (c) head at the downstream 1/4 point.

Analogous representations are performed in Figure 11 for R9 (P = 7.94). In this case, the results obtained by both forms are close to each other and have a slight difference from the experimental ones. It seems that the two forms work much better for R2 than R9 and this is partly due to the fact that P for the former test is greater than that for the latter test. The conclusion in the cavitating flows is similar to that summarized by Ghidaoui et al. (2002) in non-cavitating flows. This is because when the value of P is of the order of 1, the wave-induced shear pulse diffuses out of the viscous sublayer into the buffer zone and core region. As a result, the velocity profile does not present a logarithmic-like shape and the turbulent production and dissipation are different from their initial steady state characteristics in the region outside the viscous sublayer. Therefore, the structure and intensity of turbulence experience significant changes during the transient event, which makes the frozen and quasi-steady form questionable. Unfortunately, there are no available data that could be used to judge on the better form in such a complex phenomenon. Therefore, there is an urgent need for a systematic experimental study so as to determine the more accurate forms.

Figure 11

Comparison between the computed heads by five-region turbulence model (quasi-steady form and frozen form) and the experiment for R9 at different locations, (a) head at the valve, (b) head at the upstream 1/4 point, (c) head at the downstream 1/4 point.

Figure 11

Comparison between the computed heads by five-region turbulence model (quasi-steady form and frozen form) and the experiment for R9 at different locations, (a) head at the valve, (b) head at the upstream 1/4 point, (c) head at the downstream 1/4 point.

CONCLUSIONS

Based on the governing equations of the vapor-liquid homogeneous mixture, 1D and 2D models for transient cavitating flow were constructed, respectively. By comparing the experimental pressure curves of the quick closing of the downstream valve in a simple RPV system, the main findings of the paper are summarized below:

  • (1)

    For the transient cavitating flow considered here, the optimal ratio of the number of grid reaches (Nx) and the number of discrete cylinders (Nr) is found to be about 5.0 in the 2D model.

  • (2)

    The results predicted by the 2D model are in better agreement with the experimental pressure head traces including magnitude and phase than those predicted by the 1D model in most tests. Compared with the experimental results, the pressure curve predicted by the 1D model under most conditions has a phase lag phenomenon and the discrepancies of phase differences are magnified with the increase of time.

  • (3)

    The relative errors on maximum heads predicted by the 2D model are generally smaller than those of the 1D model. In cases of conspicuous vapor cavities, the relative errors of the 2D model are relatively larger.

  • (4)

    The pressure curves predicted by five-region and two-region turbulence models are similar, which indicates that the transient cavitating flows are insensitive to the magnitude and distribution of the eddy viscosity in the core region. Therefore, future efforts in turbulence modeling of transient flow should be concentrated on the establishment of accurate models for the near-wall region.

  • (5)

    If the non-dimensional parameter P is much greater than 1, the turbulence distribution rapidly changes and the velocity profile has enough time to form the fully-developed velocity. So the quasi-steady and frozen forms are both acceptable. If the non-dimensional parameter P is of the order of 1, the quasi-steady and frozen forms become less accurate.

ACKNOWLEDGEMENTS

This work was supported by the National Natural Science Foundation of China (Grant Nos 51779257, 51479196, 51179192), and the Program for New Century Excellent Talents in University (NCET) (Grant No. NETC-10-0784).

REFERENCES

REFERENCES
Adamkowski
,
A.
&
Lewandowski
,
M.
2012
Investigation of hydraulic transients in a pipeline with column separation
.
Journal of Hydraulic Engineering
138
(
11
),
935
944
.
Bergant
,
A.
,
Karadžić
,
U.
,
Vítkovský
,
J.
,
Vušanović
,
I.
&
Simpson
,
A. R.
2005
A discrete gas-cavity model that considers the frictional effects of unsteady pipe flow
.
Journal of Mechanical Engineering
5111
(
11
),
692
710
.
Bratland
,
O.
1986
Frequency-dependent friction and radial kinetic energy variation in transient pipe flow
. In:
Proc. 5th Int. Conf. on Pressure Surges
,
BHRA Fluid Engrg
,
Germany
.
Ghidaoui
,
M. S.
,
Mansour
,
S. G. S.
&
Zhao
,
M.
2002
Applicability of quasisteady and axisymmetric turbulence models in water hammer
.
Journal of Hydraulic Engineering
128
(
10
),
917
924
.
Greenblatt
,
D.
&
Moss
,
E. A.
1999
Pipe-flow relaminarization by temporal acceleration
.
Physics of Fluids
11
(
11
),
3478
3481
.
He
,
S.
&
Jackson
,
J. D.
2000
A study of turbulence under conditions of transient flow in a pipe
.
Journal of Fluid Mechanics
408
,
1
38
.
Jang
,
T. U.
,
Wu
,
Y.
,
Xu
,
Y.
,
Newman
,
J.
&
Sun
,
Q.
2016
Efficient quasi-two-dimensional water hammer model on a characteristic grid
.
Journal of Hydraulic Engineering
142
(
12
),
060160191
060160196
.
Kita
,
Y.
,
Adachi
,
Y.
&
Hirose
,
K.
1980
Periodically oscillating turbulent flow in a pipe
.
JSME International Journal
23
(
179
),
656
664
.
Korbar
,
R.
,
Virag
,
Z.
&
Šavar
,
M.
2014
Truncated method of characteristics for quasi-two-dimensional water hammer model
.
Journal of Hydraulic Engineering
140
(
6
),
040140131
040140137
.
Mirzaei-Shiri
,
R.
,
Ghader
,
S.
,
Mazraeh-Farahani
,
M.
&
Aliakbari-Bidokhti
,
A.
2017
Numerical solution of the shallow water equations using fourth-order compact MacCormack scheme
.
Journal of the Earth & Space Physics
43
(
1
),
209
228
.
Ohmi
,
M.
,
Kyomen
,
S.
&
Usui
,
T.
1985
Numerical analysis of transient turbulent flow in a liquid line
.
Bulletin of JSME
28
(
239
),
799
806
.
Pezzinga
,
G.
&
Cannizzaro
,
D.
2014
Analysis of transient vaporous cavitation in pipes by a distributed 2D model
.
Journal of Hydraulic Engineering
140
(
6
),
040140191
0401401910
.
Pezzinga
,
G.
1999
Quasi 2D model for unsteady flow in pipe networks
.
Journal of Hydraulic Engineering
125
(
7
),
676
685
.
Pinho
,
J.
,
Lema
,
M.
,
Rambaud
,
P.
&
Steelant
,
J.
2014
Multiphase investigation of water hammer phenomenon using the full cavitation model
.
Journal of Propulsion & Power
30
(
1
),
105
113
.
Riasi
,
A.
,
Nourbakhsh
,
A.
&
Raisee
,
M.
2009
Unsteady velocity profiles in laminar and turbulent water hammer flows
.
Journal of Fluids Engineering
131
(
12
),
121202
.
Riasi
,
A.
,
Nourbakhsh
,
A.
&
Raisee
,
M.
2013
Energy dissipation in unsteady turbulent pipe flows caused by water hammer
.
Computers & Fluids
73
(
6
),
124
133
.
Shu
,
J. J.
2003
Modelling vaporous cavitation on fluid transients
.
International Journal of Pressure Vessels & Piping
80
(
3
),
187
195
.
Silva-Araya
,
W. F.
1997
Computation of energy dissipation in transient flow
.
Journal of Hydraulic Engineering
123
(
2
),
108
115
.
Silva-Araya
,
W. F.
&
Chaudhry
,
M. H.
2001
Unsteady friction in rough pipes
.
Journal of Hydraulic Engineering
127
(
7
),
607
618
.
Simpson
,
A. R.
1986
Large Water Hammer Pressure due to Column Separation in Sloping Pipes
.
University of Michigan Dept. of Civil Engineering
,
Ann Arbor, USA
.
Simpson
,
A. R.
&
Bergant
,
A.
1994
Numerical comparison of pipe-column-separation models
.
Journal of Hydraulic Engineering
120
(
3
),
361
377
.
Simpson
,
A. R.
&
Bergant
,
A.
1995
Developments in pipeline column separation experimentation
.
Journal of Hydraulic Research
33
(
3
),
429
430
.
Storli
,
P. T.
&
Nielsen
,
T. K.
2011
Transient friction in pressurized pipes. II: two-coefficient instantaneous acceleration–based model
.
Journal of Hydraulic Engineering
137
(
6
),
679
695
.
Streeter
,
V. L.
1983
Transient cavitating pipe flow
.
American Society of Civil Engineers
109
(
11
),
1407
1423
.
Sumam
,
K. S.
,
Thampi
,
S. G.
&
Sajikumar
,
N.
2010
An alternate approach for modelling of transient vaporous cavitation
.
International Journal for Numerical Methods in Fluids
63
(
5
),
564
583
.
Szymkiewicz
,
R.
2002
Developments in unsteady pipe flow friction modelling
.
Journal of Hydraulic Research
40
(
5
),
647
656
.
Urbanowicz
,
K.
,
Zarzycki
,
Z.
&
Kudźma
,
S.
2012
Universal weighting function in modeling transient cavitating pipe flow
.
Journal of Theoretical & Applied Mechanics
50
(
4
),
889
902
.
Vardy
,
A. E.
&
Hwang
,
K.-L.
1991
A characteristics model of transient friction in pipes
.
Journal of Hydraulic Research
11
(
5
),
391
411
.
Wang
,
H.
,
Zhou
,
L.
,
Liu
,
D.
,
Karney
,
B.
,
Wang
,
P.
,
Xia
,
L.
,
Ma
,
J.
&
Xu
,
C.
2016
CFD approach for column separation in water pipelines
.
Journal of Hydraulic Engineering
142
(
10
),
04016036
.
Zhao
,
M.
&
Ghidaoui
,
M. S.
2003
Efficient quasi-two-dimensional model for water hammer problems
.
Journal of Hydraulic Engineering
129
(
12
),
1007
1013
.
Zhao
,
M.
&
Ghidaoui
,
M. S.
2006
Investigation of turbulence behavior in pipe transient using a k-ε model
.
Journal of Hydraulic Research
44
(
5
),
682
692
.
Zielke
,
W.
1968
Frequency dependent friction in transient pipe flow
.
Journal of Fluids Engineering
90
(
1
),
109
115
.