Analysis of Sub-Atmospheric Pressures during Emptying of an Irregular Pipeline without an Air Valve Using a 2D CFD Model

: Studying sub-atmospheric pressure patterns in emptying pipeline systems is crucial because these processes could cause collapses depending on the installation conditions (the under-ground pipe covering height, type, ﬁll, and pipeline stiffness class). Pipeline studies have focused more on ﬁlling than on emptying processes. This study presents an analysis of the following variables: air pocket pressure, water velocity, and water column length during the emptying of an irregular pipeline without an air valve by two-dimensional computational ﬂuid dynamics (2D CFD) model simulation using the software OpenFOAM. The mathematical model predicts the experimental values of the study variables. Water velocity vectors are also analysed within the experimental facility, assessing the sensitivity of the drain valve to different openings and changes in water column length during the hydraulic phenomenon.


Introduction
Air pockets in pressurized pipelines must be studied because they can cause different problems in the system, such as sub-atmospheric pressure situations, flow regime interruptions, reduced pump and turbine efficiency, fatigue or rupture of the pipeline, seals and/or valves of the hydraulic installation, changes in the properties of the fluids, and even erroneous readings in water meters [1][2][3].
The emptying process is a frequent operation, which should be examined in detail because sub-atmospheric pressure patterns can cause problems to pipelines depending on its characteristics, such as cover depth, pipeline stiffness and backfill type, especially when air is not admitted during the emptying operation in pressurized systems [4]. Not only sub-atmospheric pressure conditions can be reached during emptying operations but also hydraulic events with a single phase (water) during pump stoppages, regulating valve closures, and rupture of hydraulic installations. In these situations, negative pressure can be reached in installations, such as water supply networks, sewer systems, hydropower plant, among others [4,5].
These sub-atmospheric pressure values can also affect different hydraulic elements (accessories, pumps, valves, and turbines) when not taking these pressures into account. Sub-atmospheric pressures are generated as air pockets which expand during this transient phenomenon. To avoid this problem, air valves should be installed, especially at high points of the pipeline profile, but if these valves are not properly designed, the risk remains in the pipelines. Therefore, the compressible fluid dynamics of entrapped air in pipeline systems must be known to avoid these undesirable scenarios.
Recent studies have analysed pipeline emptying and filling operations to determine the behaviour of these processes, thereby enabling planners and designers to make decisions to ensure enhanced safety by providing the required system reliability [6][7][8]. The current literature includes few studies on the emptying process in water pipelines. Laanearu et al. [9] developed a semi-empirical model for studying the behaviour of the drainage or emptying process using pressurized air in the installation to estimate local energy losses. While Karadžić [10] studied similar developments on drainage manoeuvres in a pipeline and Coronado-Hernández et al. [11] developed a mathematical model to analyse this operation. Fuertes-Miquel et al. [12] developed a one-dimensional mathematical model for simulating the emptying process of a simple pipeline in two scenarios: with and without this vacuum valve. Coronado-Hernández et al. [4] validated the one-dimensional mathematical model for a pipeline with an irregular profile and with an air valve.
In turn, Besharat et al. [13] proposed 2D computational fluid dynamics (CFD) models as a resource to understand emptying processes. Martins et al. [14] studied the filling process in a pipeline by 3D CFD simulation. Zhou et al. [15] studied the dynamics of an entrapped air pocket in a pipeline during the filling process by using a 3D CFD model.
The present study is supported on previous measurements of emptying a pipeline with an entrapped air pocket in the upstream end of the installation without admitting air [13]. In practice, the scenario under study corresponds to the situation where the air valve can present failures in its operation due to lack of maintenance, when the manhole where the air valve was installed is flooded, or when this device was not installed. During the research, the effect of the air pocket size and the percentage and time of drain valve opening was also assessed [16]. Few studies report simulations of the emptying process of pipelines with irregular profiles and without the vacuum valve using computational fluid dynamics. For this reason, this line of research is further pursued in this study.
Emptying processes were simulated using one-dimensional mathematical model [11,12]. In this study, the drainage process of a pipeline with an irregular profile and without an air valve is simulated using two-dimensional computational fluid dynamics (2D CFD) with the OpenFOAM software reducing the computational simulation time in comparison with previous models [13]. This research contains the following sections: (i) considerations for the application of the one-dimensional mathematical model [11], which was used in previous researches; (ii) a description of the used experimental set-up used for validating the proposed model in the current research; (iii) the formulation and assumptions of the 2D CFD model for studying the emptying process, which has not been before reported in the literature; and (iv) results and discussion about the implementation of the 2D CFD are presented.

One-Dimensional Mathematical Model
The one-dimensional mathematical model developed by Coronado-Hernández et al. [11] analyses the emptying process of a pipeline with an irregular profile and without air valves. This condition could generate low sub-atmospheric pressure values, which could cause the collapse of the system. These mathematical models present the following assumptions: (i) the rigid column model is used to simulate the water phase; (ii) the friction factor is constant throughout the transient phenomenon; (iii) the polytropic model is used to determine the air pocket thermodynamics; (iv) piston flow is assumed for the mobile air-water interface; and (v) a symmetry plane is used to reduce by half the computational effort.
The drain valves located in the downstream ends of the installation are considered as boundary conditions. With respect to the initial conditions, the water velocity is considered static, and the pressure of entrapped air pockets in the system is at atmospheric pressure. When opening the drain valves, the air pockets depressurize rapidly, reaching low values of atmospheric pressure.
The equations that describe the hydraulics and thermodynamics of the system are presented on Table 1. Table 1. Equations for an irregular pipeline system.

No. Equations
1. Water phase equation 2. Air phase equation 3. Air-water interface equation where P * i is the absolute pressure of the air pocket i, i can be 1, 2, 3 . . . ; ρ w is the density of water, P * atm is the atmospheric pressure (101,325 Pa), ∆Z j is the elevation difference, A is the cross-sectional area, Q t is the total discharge in the drain valve, K s is the flow factor of the drain valve s, L j is the total length of the pipeline, where j can be 1, 2, 3 . . . ; L w is the length of the water column and V is the velocity.

Experimental Model
To perform the experiments, an installation was prepared by building a polyvinyl chloride (PVC) pipeline with an irregular profile and with a total length of 7.3 m and with a nominal diameter of 63 mm. The condition under analysis presents two water columns emptied by two drain valves located at the downstream ends. The valves have a nominal diameter of 25 mm.
The experimental installation was built at the Civil Engineering Research and Innovation for Sustainability (CERIS), School of Engineering (Instituto Superior Técnico-IST), University of Lisbon (Universidad de Lisboa-ULisboa), in Portugal [11]. Table 2 outlines the characteristics of the experimental runs, which are reflected in Figure 1. * The valve position is completely closed at t = 0. After that, a valve opening manoeuvre is performed until reaching a partial or full opening with time varying from 0.15 to 0.75 s, where these conditions remain constant until the end of the hydraulic event. These manoeuvres were simulated in the 2D CFD model. ** Where X 1 corresponds to the total air pocket size (X 1 = X 1,1 + X 1,2 + 0.04).

D-CFD Model with OpenFOAM
Computational fluid dynamics (CFD) is a technique used to study systems involving gas and liquid flow for biphasic flow modelling [17]. Finite volume theory was used for this purpose. This theory results from the continuity and conservation of momentum equations by Navier-Stokes [18]. Table 3 outlines the equations used for the simulation.

No.
Formulation Equation

D-CFD Model with OpenFOAM
Computational fluid dynamics (CFD) is a technique used to study systems involving gas and liquid flow for biphasic flow modelling [17]. Finite volume theory was used for this purpose. This theory results from the continuity and conservation of momentum equations by Navier-Stokes [18]. Table 3 outlines the equations used for the simulation. Table 3. Equations proposed by Navier-Stokes.

No.
Formulation Equation

Moment conservation equation
2. Air-water interface equation where ρ m is the mixed density of the fluid, V is the velocity, P is the pressure, → g is the gravitational acceleration, T is the temperature, µ m is the mixed viscosity of the fluid, t is the time, and ∇ is the vector gradient.
In the CFD simulation, the software OpenFOAM, 2006 version (2020), was used to model the interaction between two fluids, water and air. A multiphase flow solver was used to represent the experimental conditions. In the simulation, water was considered a compressible, nonisothermal, immiscible fluid. An interface capture approach based on VOF (volume of fluid) phase fraction was used, with optional mesh movement and mesh topology changes including adaptive remeshing. The simulation of the ball valves in the installation required developing dynamic meshing, which was supported by the selected solver.
The simulation runs for 3 s where the valve manoeuvre time for all runs is 0.7 s. It is necessary to clarify that the full opening of 100% corresponds to a movement of 90 • , based on the design of the valve; in this sense, the partial opening of 6% corresponds to 5.4 • .
The geometry of the model consists of 14 blocks defined by point coordinates in an XYZ plane. It is necessary to clarify that the model is in 2D, however, it covers three dimensions because they are necessary by OpenFOAM, so the z plane is imposed with a thickness of 0.005 m. The model has 10 boundaries defining the water-air inlets and outlets. The hydraulic system is 7.3 m long, with a nominal diameter of 63 mm and with a 25 mm reduction of the nominal diameter at the ends (23.66 mm internal diameter), that is also used as the lateral dimension. Figure 2 shows the experimental configuration. A symmetry plane was used to simplify the hydraulic system to a pipe with a single branch and a single valve to assess the results of the complete system when applying the mirror function (see Figure 2b).
The meshing of the model consists of blocks cubically defined by structured cells, which were gradually refined in a descending direction towards the areas near the pipe walls. A difference range of 20% is maintained between the area of each cell of a block and the area of the cells of its adjacent block. The final mesh has 52,984 cells, with a block cell size in the valve of 0.00147875 m × 0.00124526 m (see Figure 2a,b). The drain ball valve was represented by defining rotation values with a moving screen to determine the opening of the valve. It should be noted that the main point of comparison by which it was determined that this was the most accurate mesh for the project was the observed evolution of the pressure curve in the air pockets of the pipe.
The following initial conditions were used: for the air pockets, an absolute pressure of 10.33 m (101,325 Pa); an ambient temperature of 20 • C (293.15 • K), which was taken when performing the experiments; and an initial water velocity of 0 m/s, for starting the fluid at rest. Based on this, it was determined that the boundary conditions for each of the variables would be as shown on Table 4. Table 4. Boundary conditions for the variables of the turbulent model.

Variables
Outlet Walls inletOutlet KqRWallFunction and the area of the cells of its adjacent block. The final mesh has 52,984 cells, with a block cell size in the valve of 0.00147875 m × 0.00124526 m (see Figure 2a,b). The drain ball valve was represented by defining rotation values with a moving screen to determine the opening of the valve. It should be noted that the main point of comparison by which it was determined that this was the most accurate mesh for the project was the observed evolution of the pressure curve in the air pockets of the pipe. The following initial conditions were used: for the air pockets, an absolute pressure of 10.33 m (101,325 Pa); an ambient temperature of 20 °C (293.15 °K), which was taken when performing the experiments; and an initial water velocity of 0 m/s, for starting the fluid at rest. Based on this, it was determined that the boundary conditions for each of the variables would be as shown on Table 4. Table 4. Boundary conditions for the variables of the turbulent model.

Variables
Outlet The motion of a turbulent fluid is an irregular flow condition in which the various quantities show a random variation with the coordinates of time and space, so that statistically different average values can be discerned [19]. The most used turbulent models in CFD are the models made up of two differential equations that take into account variables that handle the energy generated and dissipated by the phenomenon. Among these are the k-ε model [20], the k-ω model [21], and the k-ω SST-Shear Stress Transport model [22], which combines robustness of the k-ԑ with the accuracy of the k-ω model, The motion of a turbulent fluid is an irregular flow condition in which the various quantities show a random variation with the coordinates of time and space, so that statistically different average values can be discerned [19]. The most used turbulent models in CFD are the models made up of two differential equations that take into account variables that handle the energy generated and dissipated by the phenomenon. Among these are the k-ε model [20], the k-ω model [21], and the k-ω SST-Shear Stress Transport model [22], which combines robustness of the k-ε with the accuracy of the k-ω model, which makes it applicable to a wide range of flows. Additionally, it is much better at predicting flows with strong pressure gradients.
This model was chosen because it shows low sensitivity in the results of the variables studied (pressure and velocity) and improves the precision in predicting flows with strong adverse pressure gradients [22]; ensuring better transitions within the flow.
During the calculation process, the different equations that govern the flow over all the control volumes (finite) of the domain are solved, as can be seen in Table 5, through discretization processes that allows obtaining the solutions by means of iterative methods. Table 5. Turbulence k-ω SST model equations.

Formulation Equation
where k is the turbulent kinetic energy, ω is the dissipation frequency, F 1 is a blending function, ρ and u i correspond to the density and velocity of the specific flow, respectively; µ and µ t correspond to the laminar and turbulent dynamic viscosity, respectively; ν t is the turbulent kinematic viscosity; β and β* are constants considered by the Model, where β* = 0.09 and β depends on β 1 = 0.075 and β 2 = 0.0828. Table 6 shows the variables added in the initial and boundary conditions due to the turbulent model of the system (α t , k, ω y ν t ). The geometry at its lower left end presents a ball valve, which is represented in a dynamic mesh block as seen in Figure 2a. The dynamicMeshDict file, created to control the transformation and deformation of the mesh during the simulation, uses a rigid body motion function that allows the valve to rotate on its own axis through the z-axis. The valve opening is controlled by a rotation vs. time data tab in the datafile.dat file. A synthetic manoeuvring was considered during simulations.

Results
Considering the cases outlined in Table 2, the simulations were performed with OpenFOAM using a 2D CFD model.

Analysis of the Variables
The variables absolute pressure of the air pocket and water velocity were analysed to determine the fit of the 2D CFD model to the experimental results.
The variables were analysed for the two types of patterns identified in the experiments (see Table 2).

Air Pocket Pressure Type A Pattern: Partial Opening of the Ball Valve
Partly opening the ball valve at 6% produces a minimum variation in the water-air interface. For simulations A1, A2, and A3, there is no backflow air intrusion through the ball valves, and the absolute pressure of the air pocket decreases without oscillations. The 2D CFD OpenFOAM model correctly predicts the behaviour of the air pocket pressure, as shown in Figure 3. Two experimental tests were performed for all cases. The conducted comparison only shows the average experimental results. Figure 4 presents repetitions of case A2. In the presentation of results, the data are presented as average experimental.
For pattern A1, a minimum air pocket pressure of 9.71 mwc was reached in the experiments (at 1.35 s), while for the 2D CFD model, a value of 9.68 mwc was reached (at 1.14 s). After reaching the minimum values of air pocket pressure, small oscillations persist until the water column flow comes to a complete stop (see Figure 3a).
For pattern A2, the experimental measurement reaches a minimum pressure value of 9.88 mwc in 2.1 s, while for the 2D CFD model, a minimum pressure of 9.86 mwc is reached in 1.63 s, showing the accurate prediction of this tool. Figure 3b presents the sub-atmospheric pressure patterns, which tend to a value of 9.88 mwc.
The experimental results of pattern A3 reach a minimum air pocket pressure of 10.07 mwc (at 2.41 s), while for the 2D CFD model, a minimum value of 10.05 mwc is reached in 1.82 s. In this case, the pressure values are 0.28 mwc below the atmospheric pressure, with no fluctuations in this pattern (see Figure 3c). ball valves, and the absolute pressure of the air pocket decreases without oscillations. The 2D CFD OpenFOAM model correctly predicts the behaviour of the air pocket pressure, as shown in Figure 3. Two experimental tests were performed for all cases. The conducted comparison only shows the average experimental results. Figure 4 presents repetitions of case A2. In the presentation of results, the data are presented as average experimental.

Corrections
1. Where k is the turbulent kinetic energy, ω is the dissipation frequency, F_1 is a blending function, ρ and u_i correspond to the density and velocity of the specific flow, respectively; μ and μ_t corresponds to the laminar and turbulent dynamic viscosity, respectively; ν_t is the turbulent kinematic viscosity; β and β * are constants considered by the Model, where β * = 0.09 and β depends on β1 = 0.075 and β2 = 0.0828. Table 6 shows the variables added in the initial and boundary conditions due to the turbulent model of the system (α_t, k, ω y ν_t). It's correct. 2. For figure 4, this is the new revised version: 3.
4. On table 3, please replace the old one for this one:

Type B Pattern: Full Opening of the Ball Valve
The comparison between experimental measurements and the 2D CFD OpenFOAM model shows very similar results, indicating robustness in the application of our model. During the analyses, the size of the air pocket affected the variation of air pocket pressures. Figure 5 shows that the most critical scenario (minimum pressures) is reached for small air pockets. When the ball valves are fully opened (type B pattern), the pressure values are lower than the pressures found in the type A pattern (partial opening of 6%). Figure 5 shows that the 2D CFD OpenFOAM model correctly predicts the behaviour of the air pocket pressure for patterns B4 to B6.
For case B4, the experimental data reached a minimum air pocket pressure of 9.40 mwc (at 0.

Water Velocity
Water velocity was measured in the horizontal section of the pipeline at a distance of 1.15 m from the location of the ball valve. The water velocity was measured using an Ultrasonic Doppler Velocimeter (UDV).
Type A Pattern: Partial Opening of the Ball Valve The velocity flow in the type A pattern could not be experimentally measured be-  Figure 6 presents the comparison between measured and computed water velocity. The 2D CFD model can predict extreme values of the water velocity pattern; however, it cannot capture oscillations. The UDV works considering the reflection principle, using suspended small seeding particles inside the pipeline installation. Water velocities are very low during the analysed experiments. In this sense, there are uncertainties related to the measured water velocity pulses.  Figure 7 shows the velocity vectors produced by the water-air interaction during the water emptying process without an air valve. In Figure 7, red arrows represent the air velocity vector, while blue arrows represent the water velocity vector. At 0.5 s (partly open valve), the internal air flow has not entered the pipeline with a nominal diameter of 25 mm, while at 1.0 s, air enters the hydraulic system (see Figure 7b).  Figure 7 shows the velocity vectors produced by the water-air interaction during the water emptying process without an air valve. In Figure 7, red arrows represent the air velocity vector, while blue arrows represent the water velocity vector. At 0.5 s (partly open valve), the internal air flow has not entered the pipeline with a nominal diameter of 25 mm, while at 1.0 s, air enters the hydraulic system (see Figure 7b).
(c) Figure 6. Comparison of water velocity patterns between the experimental model and the 2D CFD OpenFOAM symmetry model for cases B: (a) run B4; (b) run B5; (c) run B6. Figure 7 shows the velocity vectors produced by the water-air interaction during the water emptying process without an air valve. In Figure 7, red arrows represent the air velocity vector, while blue arrows represent the water velocity vector. At 0.5 s (partly open valve), the internal air flow has not entered the pipeline with a nominal diameter of 25 mm, while at 1.0 s, air enters the hydraulic system (see Figure 7b).

Water Column Analysis
As the water exits the pipeline, the air pocket lengthens; for this reason, the variation of the length of the water column (Lw) is analysed, as shown in Figure 7. The water column in case A2 has an initial length of 3.16 m and a length of 3.14 m at the end of the

Water Column Analysis
As the water exits the pipeline, the air pocket lengthens; for this reason, the variation of the length of the water column (L w ) is analysed, as shown in Figure 7. The water column in case A2 has an initial length of 3.16 m and a length of 3.14 m at the end of the hydraulic process, thus indicating that, by restricting the ball valve opening, only 0.016 m can be drained in 2 s during this process. In contrast, the water column of case B5 begins with a length of 3.16 m (at t = 0 s) and ends with 3.13 m, indicating that 0.21 m can be drained in 2 s during the process. And in the same way, the process is studied for all the following cases which observe how cases A1 and B4 have an air pocket of 0.28 m and a water column of 3.49 m at the initial time; however, by time 2 s it has a minimum water column of 3.48 m and 3.46 m.
For cases, A3 and B6, an initial air pocket of 0.98 and a water column of 2.79 were stipulated, and by time 2 s there will be a column of 2.78 and 2.77, showing that the cases of rapid opening tend to decrease much more water column during the emptying process. The results can be analyse on Figure 8.

Discussion
To assess the effect and importance of correctly opening the ball valve, this operation was subjected to sensitivity analysis in values ranging from 2 to 10%, as shown in Figure 9. Less critical results were found at larger openings, also generating pressure oscillations, which do not occur at smaller openings, whose pressure data remain more stable. In both cases, the variation of the model over time is also described as the system reaches the pressures assessed when opening the valve at 6%. Figure 10 shows a comparison between the results using the 2D CFD model (Open-FOAM symmetry) and the one-dimensional mathematical model. Figure 10a,b presents information regarding cases A1 and B5, respectively. Both models can predict the extreme drop of sub-atmospheric pressure conditions and follow the air pocket pressure pattern.
The 2D CFD OpenFOAM symmetry model has some limitations. When defining the limit of symmetry, which functions as a mirror surface, some cases show different lengths of air pockets on the (left and right) sides of the pipeline, and in other cases, the pipeline is not symmetrical. It should be noted that, as mentioned above, the experimental opening times may vary among themselves, allowing the pressures to vary in future projects. Furthermore, these can vary according to the type of valve designed and the opening mode used in the simulation.
To assess the effect and importance of correctly opening the ball valve, this operation was subjected to sensitivity analysis in values ranging from 2 to 10%, as shown in Figure 9. Less critical results were found at larger openings, also generating pressure oscillations, which do not occur at smaller openings, whose pressure data remain more stable. In both cases, the variation of the model over time is also described as the system reaches the pressures assessed when opening the valve at 6%.  Figure 10 shows a comparison between the results using the 2D CFD model (OpenFOAM symmetry) and the one-dimensional mathematical model. Figure 10a,b presents information regarding cases A1 and B5, respectively. Both models can predict the extreme drop of sub-atmospheric pressure conditions and follow the air pocket pressure pattern. The 2D CFD OpenFOAM symmetry model has some limitations. When defining the limit of symmetry, which functions as a mirror surface, some cases show different lengths of air pockets on the (left and right) sides of the pipeline, and in other cases, the pipeline is not symmetrical. It should be noted that, as mentioned above, the experimental opening times may vary among themselves, allowing the pressures to vary in future projects. Furthermore, these can vary according to the type of valve designed and the opening mode used in the simulation.

Conclusions
The 2D CFD OpenFOAM symmetry model developed in the present study for type A and type B patterns shows a good fit compared to previous experimental results and the one-dimensional model, demonstrating the accuracy of the analysis of 2D computational fluid dynamics in predicting sub-atmospheric pressure and flow velocity in the

Conclusions
The 2D CFD OpenFOAM symmetry model developed in the present study for type A and type B patterns shows a good fit compared to previous experimental results and the one-dimensional model, demonstrating the accuracy of the analysis of 2D computational fluid dynamics in predicting sub-atmospheric pressure and flow velocity in the system. The above is numerically evidenced by the following findings: - The type A pattern is the least critical scenario in the pipeline, because the subatmospheric pressure pattern reaches values close to atmospheric pressure. The hydraulic event begins with an atmospheric pressure (10.33 m) and then decreases rapidly until reaching sub-atmospheric pressure lows of 9.683 m (at 1.14 s), 9.857 m (at 1.63 s), and 10.045m (at 1.82 s) in runs no. 1 and no. 2 and o. 3, respectively. This pattern shows no oscillations in the first two seconds due to the partial opening of the valve. Subsequently, minimal fluctuations are observed in the pressure stabilisation that occurs in different cases as a function of the air pocket size. The more water left to exit the pipeline, the more its pressure will vary. - In the type B pattern, the valves are completely open. The sub-atmospheric pressure pattern in the hydraulic event begins under atmospheric conditions and then decreases rapidly until reaching minimum values of 9.387 m (at 0.38 s), 9.624 m (at 0.46 s) and 9.90 m (at 0.53 s) in runs no. 4, no. 5 and no. 6, respectively. After reaching minimum values, oscillations occur throughout the hydraulic event. As the air pocket increases, the pressure fluctuations within the pipeline become increasingly closer to atmospheric pressure, decreasing the amplitude of pressure changes and increasing the time between the different sub-atmospheric drops.
An adequate simulation of regulating valve opening is needed, since the difference of the prediction of air pocket pressure pattern were found in the current research; finally, it can be said that the simulation of the 2D model is a simplification with respect to a fully 3D transient model, but the benefits of reduced computational cost are very important in 2D modelling. The main advantage of this 2D approach with respect to the 1D approach is the ability to simulate the valve opening. The 1D model represents the experimental results in a better way than the 2D model. In this sense, 3D CFD models are required to improve results considering higher computational times.