Analysis of pressure surges for water filling in deep stormwater storage tunnels with entrapped air-pocket using a VOF model

The deep stormwater storage tunnel (DSST) system is a new type of municipal engineering dealing with large-scale urban water circulation problems. It relieves pressures of underground pipe networks through planned and targeted water storage and drainage, in order to reduce urban local waterlogging and sewage confluence pollution. DSST systems have been practically applied in large cities such as Chicago (Cataño-Lopera et al. 2014) and Shanghai (Wu et al. 2016), and have achieved relatively good drainage and storage effects. However, during the storage process with rapid water filling from shafts, it is unrealistic to completely discharge air in the DSST and thus air-pocket will be entrapped. With the flow pattern of the tunnel changing from an open channel flow to a pressurized flow, the pressurized rainwater impacting on the entrapped air-pocket will cause abnormal pressure surges, and even endanger the safety of the DSST system (Vasconcelos & Wright 2005).

The pressure surges with the entrapped air-pocket are subordinate to typical air–water two-phase flow issues in pressurized pipelines. For pipelines with entrapped air-pocket, most of the research works mainly focus on constructing a mathematical model, including the elastic water model (EWM) and rigid column model (RCM). Among them, the EWM is widely used in research works of pressure surges in closed pipelines (Guarga et al. 1996) or air released in a single hole (Zhou et al. 2002). For a closed pipe containing air-pocket in the end or air releasing in an exhaust hole, the air-pocket impacted by pressurized water with valve opening will produce several times direct water hammer pressure (Zhou et al. 2011a, 2011b, 2018), and the elasticity of water cannot be ignored in this event. In contrast, the wave velocity of the air phase is much smaller than that of the liquid phase, and the elasticity of air can be ignored (Lee 2006). Except for the system dominated by elastic deformation of water, in the pipeline with entrapped air-pocket dominated by inertial force, the RCM is more convenient and the calculation accuracy is guaranteed (Guarga et al. 1996). In addition, the applicability of the RCM is also affected by air-pocket parameters and inlet pressure (Lee & Martin 1999). In fact, for rapid water filling impacting entrapped air-pocket at the end of a pipe, the pressure simulation result of the RCM is consistent with that of the experiment (Huang & Zhu 2021). However, the assumption that the water–air interface is perpendicular to the conduit axis in the RCM is difficult to apply in the complex alternating process of open and full flow (Hou et al. 2014), especially under the condition of multiple inflows of the DSST system.

In addition to the mathematical models, the volume of fluid (VOF) method is a computational fluid dynamics (CFD) two-phase flow model solving the stratified flow of incompatible liquids (Hirt & Nichols 1981). It can describe the transient deformation of the water–air interface, which is also applicable to pipelines with entrapped air-pocket (Hargreaves et al. 2007), especially for alternating open channel and pressurized flow. Some experiments show that, for geysers caused by air discharge through different orifices to the shafts, the VOF model simulation results were very close to the experiments, especially for internal flow characteristics and transient movement of the water–air interface (Zhou et al. 2011a, 2011b; Cong et al. 2017; Qian et al. 2020). The other experiments reflect that there is an obvious heat exchange process in the air–water interaction in the elbow pipe of the horizontal–vertical pipeline, and obvious white fog and local temperature rise can be observed in the experiment and the VOF simulation separately (Zhou et al. 2018). For a similar pipeline with water filling simulated by VOF, the maximum transient pressure in the release process can reach 70 m (Martins et al. 2017). In addition, the VOF model can also be used to calculate the downstream orifice wall deformation with the water column impacting air-pocket, coupled with solid analysis (Maddahian et al. 2020). Although the VOF model can accurately locate the air-pocket position and describe the process of air-pocket release, it lacks applicability for the multi-inlet and multi-branch DSST system for a long time of simulation.

The existing experimental research of the DSST system shows that air-pocket entrapped is unavoidable due to the complicated situation of open channel flow reflected waves in the water filling process with multi-shaft inflow (Vasconcelos & Wright 2005, 2009). However, due to the unclear law of air-pocket retention, most studies still tend to predetermine specific air-pocket positions, especially to study pressure surges or geysers of the manhole with unilateral inflow or tunnel diameter changes (Wright et al. 2011; Vasconcelos & Leite 2012). For the multiple inflows model, inflow conditions have been evaluated for the statistical size and number of entrapped air-pockets in a four-shaft system (Vasconcelos & Wright 2017). In order to study pressure surges of multi-inflow water impacting entrapped air-pocket, parameter conditions of air-pocket must be given, and the multi-inflow DSST system models must be constructed, of which the simplest one is the two shafts and one tunnel system. For this simple model, the CFD three-dimensional (3D) numerical simulation based on the VOF model can undoubtedly monitor the pressure surges and describe the water–air interface transient changes. The present work aims to study the mechanism for pressure surges of DSST systems with entrapped air-pocket during the water filling process. Specifically, this study focuses on (1) establishing a 3D numerical model of the two shafts and one tunnel system and comparing simulation results of the water filling process with test; (2) predicting pressure surges of different initial air-pocket length/diameter ratios; (3) simulating pressure surges of different initial air-pocket volume fractions; and (4) comparing pressure surges with different initial air-pocket center positions.

In order to reduce the impact of the water filling drop of shafts on the energy equations, the bottom of two shafts is defined as flow inlet, and an extension section is added at the bottom to ensure the full development of inflow. The top of two shafts is set as pressure outlet, the wall of shafts and tunnel is treated as non-slip velocity condition. The initial water–air interface position is defined by the method of selected area to patch, while the gas phase (air) model is defined as an ideal gas and the liquid phase (water) model is defined as a compressible liquid. It is assumed that there is no heat exchange in the whole process. The simulation of water filling process needs to consider the influence of gravity, and the acceleration is defined as 9.81 m/s², pointing to Z-axis direction. The solver adopts the pressure benchmark PISO algorithm for transient numerical simulation, with a time step of 0.02 s and a convergence standard of 10⁻⁴. In order to obtain the maximum transient pressure value in the rising process of the shaft water level, the time is set as 300 s.

Using CFD models to simulate pressure surges under different initial air-pocket parameters can be realizable, but for experiment, it is extremely difficult to control the shape and position of air-pocket especially for multiple or suspended bubbles. For the same CFD physical model, describing the pressure surges of the water filling process in empty tunnel or pressurized pipe with entrapped air-pocket, the control equations and multiphase model settings are consistent. Therefore, if the simulation of empty tunnel water filling process using the CFD model is credible, it can also be reliable for the process of water impacting entrapped air-pocket. The comparison of relative head H_s at the bottom of lower tunnel entrance in the original model between experimental data and CFD simulation results is shown in Figure 3(a). It can be clearly seen that the transient pressures estimated by CFD compare very well with those obtained experimentally, and the relative difference in maximum pressure head is less than 1%. The good agreement indicates that the CFD model can be used to characterize the maximum pressure of empty tunnel water filling process in the simplified system. On this basis, taking the maximum relative head at the bottom of the tunnel entrance H_m as the parameter variable, the grid dispersion error is evaluated according to the GCI convergence index of the Richardson extrapolation method recommended by the American Society of Mechanical Engineers. According to the requirements of GCI grid sensitivity analysis (Celik et al. 2008), four grid schemes with different average unit grid volume parameters are set. The refinement coefficients of the two adjacent schemes are greater than 1.3. Based on the simulation results of the grid convergence index of the parameter variable, as shown in Figure 3(b), grid 3 is selected and the total number is about 3.89 × 10⁵.

Under the constant initial volume fraction of entrapped air-pocket, D = D_t = 10 m is the most dangerous working condition. Compared with the others, the frequency is faster and the transient energy is larger. This is due to the radial deformation of air-pocket is affected by buoyancy, and the deformation rate and inertial force increases. Through experimental observation and daily experience, it can be found that the air-pocket in the tunnel always tends to appear in the form of a larger length/diameter ratio after a period of standing. Compared with shafts full flow time (177 s), the air-pocket radial deformation time is much shorter. Therefore, under a long time of water filling process, it is difficult for the air-pocket to stay in the form of blocking water flow. Since the air-pocket occupying the tunnel cross-section is the most dangerous working condition, the structural design and water filling control obtained by the risk assessment are the most conservative. Therefore, D = D_t = 10 m is used in the following research to define the initial diameter of air-pocket.

Taking the maximum pressure and cycle time T_s as the references, the comparison of 0.2–5.0% air-pocket with different initial air volume fractions is obtained as shown in Figure 5(b). It can be seen that with the increase of φ_a, the maximum pressure and cycle time T_s all show a gradual increasing trend. Among them, the maximum pressure increased by 29.0% in total, and its growth trend gradually slowed down and finally approached about 49 m. Also, the total increase of cycle time is 282.0%, except for φ_a = 0.5%, the growth rate under other working conditions of T_s is relatively close. In general, the influence of air volume fraction on the amplitude and frequency of pressure surges is obvious, especially the frequency, which is basically linear with the air volume fraction. However, in the range of φ_a = 2.0–5.0%, the influence of air volume fraction on the maximum pressure is small, while the pressures may be more affected by the static pressure of shaft water level. It can be speculated that with the further increase of air volume fraction and gradually occupying most of the tunnel, the increase of φ_a may lead to further reduce of the maximum pressure. But due to the shape of tunnel structure and the form of water filling process, entrapped air-pocket with large volume fraction could not exist with air discharging from shafts in the DSST system.

To explain the phenomenon, the maximum pressures at different air-pocket centers and the pipeline center (P1000 m) are compared as shown in Figure 6(b), where Pothers is defined as the monitoring point located at the air-pocket center. It can be seen that, except for P250 m and P1000 m, the maximum pressures at the center of other air-pockets is higher than that of the tunnel center. The maximum pressure of the air-pocket center at different initial positions shows a downward trend from the middle to two sides, while the pressure of air-pocket near higher shaft will be slightly higher. This is because that the air-pockets close to both sides of tunnel are subjected to different pressures on both sides, and the deformation and compression of it with the movement of water flow are less. At the same time, there will be some air overflow in the shafts. Also, the maximum pressure point appears at 1,250 m instead of the center of the tunnel. Due to the influence of tunnel slope, although the initial inflow of two shafts is equal, the final outflow is different, and a large amount of water flows from the high-level shaft to the low-level shaft. At the point slightly higher than the center of the tunnel, the sum of static pressure and dynamic pressure on both sides of the entrapped air-pocket is relatively close. At this time, the compression deformation of the air-pocket is the largest, and the maximum transient pressure is also the largest.

From the unmarked air–water two-phase distribution diagram in the figure and subsequent, the change of air-pocket during water filling process can be divided into the following three stages. (1) Deformation stage on both sides (0–5.0 s): due to the upward deformation of air-pocket under the influence of buoyancy, the water–air interface on both sides gradually changes from perpendicular to the tunnel axis to the water flow incoming direction. Also, the water cross-section changes along the gradient of water–air boundary line, the water inflows on both sides near the center of air-pocket continue to concentrate, the velocity increases, and finally converges to one point. (2) Air-pocket separating stage (5.0–10.65 s): The water inflows on both sides converge to a point, which converts the kinetic energy of water into potential energy, resulting in surge at the center point. The water level rises rapidly, and finally reaches the top of tunnel, separating the air-pocket into two parts. (3) Axial deformation stage: The two separated air-pockets gradually move away from the center. With the subsequent pressure surges, the air-pocket deformation is mainly through the axial direction (due to that the length diameter ratio of the tunnel is too large 200:1, the subsequent air–water two-phase distribution is not given in the paper).

It can be seen from Figure 4(b) that the air-pocket length/diameter ratio has an influence on maximum pressure, and the largest maximum pressure occurs when D = D_t = 10 m. However, with the radial deformation of air-pocket in the early pressure rise stage of the water filling process, the difference of air-pocket height is gradually eliminated, but it has an impact on the maximum pressure in the subsequent pressure fluctuation stage. In fact, it can be seen from Figure 4(a) that when air-pocket height is smaller, the pressure rises faster in the early stage of water filling process. In this paper, it is assumed that the axial deformation inertia force caused by different air-pocket height causing this phenomenon, but there is no corresponding evidence. In addition, the air-pocket separates into two parts in the early stage of water filling process, the inevitable condition triggering this phenomenon is difficult to measure, and the influence of multi air-pockets on subsequent pressure surges is also unknown.

The DSST system, a new type of urban waterlogging drainage facility, is prone to geysers or pressure surges under the conditions of heavy rainfall and rapid inflows. This work aimed to provide a general view on the relationship between the maximum pressures of water filling process and the entrapped air-pocket parameters. Although the two shafts and one tunnel model studied cannot represent the actual geometry of the DSST system, it can fully reflect the characteristics of air-pocket deformation and pressure surges. The relevant conclusions are expected to be extended to the study of multiple shafts and multiple air-pockets.

The CFD model has been widely used in the study of pipeline transient pressures with specific entrapped air-pocket position, but it is much more difficult to control the state and position of air in the experiment. Hence, this paper verifies the accuracy of the CFD model through a water filling experiment of the two shafts and one tunnel model. The pressure surges affected by initial length/diameter ratio, air volume fraction, and the center position of air-pocket are listed as follows:

• The maximum pressures with fixed air volume fraction increase with the decrease of initial length/diameter ratios, and the extremum one is achieved with the initial air-pocket blocking water flow. The existence of air-pocket with a certain scale (φ_a ≥ 0.7%) will produce pressures that are 60% higher than the static pressure.

• The influence of air volume fraction on the amplitude and frequency of pressure surges is obvious, especially the frequency, which is basically linear with the air volume fraction. In the range of φ_a = 0.2–5.0%, the maximum pressure increased by 29.0% and the cycle time T_s increased by 282.0% in total.

• When the initial position of the air-pocket is not in the center of tunnel, the transient pressure at the corresponding position will fluctuate violently in the initial stage and tend to be stable before 50 s. The maximum pressure of air-pocket at different initial positions shows a downward trend from the middle to two sides, while the pressure of air-pocket near higher shaft will be slightly higher.

• In this paper, only the influence of initial air-pocket parameters is considered, and there is a lack of detailed and verifiable explanation the air-pocket deformation and pressure transient change in the process of water filling. In the follow-up study, the effects of multi air-pocket parameters and inflow difference between two shafts on the pressure surges of entrapped air-pocket should be studied.

The authors acknowledge the financial support given by the National Natural Science Foundation of China (Nos 52179062 and 51879087) and the Postgraduate Research & Practice Innovation Program of Jiangsu Province (Grant No KYCX22_0652).

All relevant data are included in the paper or its Supplementary Information.

The authors declare there is no conflict.