Flows in manholes are complex and may include swirling and recirculation flow with significant turbulence and vorticity. However, how these complex 3D flow patterns could generate different energy losses and so affect flow quantity in the wider sewer network is unknown. In this work, 2D3C stereo Particle Image Velocimetry measurements are made in a surcharged scaled circular manhole. A computational fluid dynamics (CFD) model in OpenFOAM® with four different Reynolds Averaged Navier Stokes (RANS) turbulence model is constructed using a volume of fluid model, to represent flows in this manhole. Velocity profiles and pressure distributions from the models are compared with the experimental data in view of finding the best modelling approach. It was found among four different RANS models that the re-normalization group (RNG) k-ɛ and k-ω shear stress transport (SST) gave a better approximation for velocity and pressure.
Manholes are one of the most common features in urban drainage networks. They are located at changes in slope and orientation of the sewer pipes, as well as at regular intervals along the pipes to enable maintenance. The flow pattern in a manhole is complex, especially during high flows, and involves several hydraulic phenomena such as local flow contraction, expansion, rotation, recirculation as well as possible air entrainment and sediment mixing. These flow phenomena can control the overall energy loss, transport and dispersion of solute and particulate materials in the manhole structure. PIV measurement can provide a good representation of the complex velocity field of a manhole. Previously Lau (2008) studied two-dimensional particle image velocimetry (PIV) in a surcharged scaled manhole. Attempts to measure stereo PIV data in a scaled manhole are, however, new and, to the authors' knowledge, have not been done before. Several researchers studied flow patterns in surcharged manholes using computational fluid dynamics (CFD) models. Use of different Reynolds Averaged Navier Stokes (RANS) modelling approach like re-normalization group (RNG) k-ɛ model (Lau et al. 2007), realizable k-ɛ (Stovin et al. 2013), k-ω model (Djordjević et al. 2013) have been reported. A little research study has been conducted on how these flow patterns could affect flow movement and flow quality in the wider piped network. In the current work, the flow phenomena of a scaled manhole are measured by stereo PIV and modelled numerically using OpenFOAM® CFD tools. Four different RANS models, i.e. RNG k-ɛ, realizable k-ɛ, k-ω shear stress transport (SST) and Launder–Reece–Rodi (LRR) were used, and the differences in flow structures among them were compared and discussed.
METHODS AND MATERIALS
The experimental facility was installed in the Hydraulic Laboratory of the University of Sheffield. It consists of a transparent acrylic circular scaled manhole, linked to a model catchment surface (Rubinato et al. 2017). The manhole has an inner diameter (Φm) of 240 mm and connected with 75 mm diameter (Φp) inlet–outlet pipes. Both pipes are co-axial, and the pipe axis passes through the centre of the manhole vertical axis (Figure 1). Two butterfly valves (one at >48Φp upstream of the manhole, the other one at >87Φp downstream of the manhole) are used to control the inflow and the water depth of the manhole, respectively. The inflow was monitored using an electromagnetic MAG flow meter fitted within the inlet sewer pipe (10Φp from the butterfly inlet valve). The ratio of the manhole diameters to inlet pipe diameters (Φm/Φp) is 3.20. Two Gems series pressure sensors (product code 5000BGM7000G3000A, serial number 551362, range 0–70 mb) were installed vertically at the inlet and outlet pipes, the first one at 350 mm upstream from the centre line of the manhole, and the second one, 520 mm downstream from the centre line of the manhole by making a hole in the pipe of Φ = 5 mm. They can measure piezometric pressures for both free surface and pressure flow conditions within the inlet–outlet pipes. These transducers were calibrated such that transducer output signal (4–20 mA) can be directly related to gauge pressure. For the calibration of the pressure transmitter, 10 different water levels were measured (range 50 to 500 mm, no discharge) and checked with a point gauge. Rubinato (2015) found an overall accuracy can be defined as ±0.72 mm. A transparent measuring tape was attached vertically to the manhole side in order to check the manhole water levels during the experiments. The tape position was at the other side of the camera, keeping an equal distance from both inlet and outlet.
A stereo PIV measurement setup was installed. Two Dantec FlowSense EO 2M cameras and an Nd:YAG pulsed laser were placed at opposite sides of the scaled manhole. Each camera resolution was 1,600 × 1,200 pixels (Figure 1) and was set at the same distance from the manhole, making more than 45° angle at the vertical centre of the measuring plane. To reduce error due to refraction through the curved manhole wall, a transparent acrylic tank was constructed around it and filled with water, keeping flat surfaces parallel to both camera lens axes. The laser was directed from the bottom of the acrylic manhole as a laser sheet with the help of a flat mirror set at 45° to the horizontal direction. The laser sheet thickness was around 4 mm.
Conventional 2D PIV can give the velocity vectors perpendicular to the camera direction only, which is typically parallel to the laser sheet (known as the in-plane velocity (Leandro et al. 2014)). However, the use of two cameras can give the reading of the third component of the velocity vector (referred as out of plane velocity) with proper regeneration of the 2D velocities from each camera's images, provided that appropriate calibration is done beforehand. As both of the cameras see the same stationary image from two different angles, the calculated velocity vectors from the two successive snapshots of a moving particle will also orient differently on both cameras. This change in orientation is the result of the different camera position and out of the plane velocity vectors. If the camera position is known and distortion of a still image at each camera is known, the out of plane velocity component can be calculated. Standard calibration plates were used in this study to calibrate the cameras.
When the manhole surcharge level is below a certain limit (typically around 0.2Φm), then the manhole inlet flow reacts vigorously with the surface, creating irregular flow pattern and a very high head loss. This surcharge limit is known as threshold surcharge (Stovin et al. 2013). At higher surcharge, the flow pattern inside the manhole is regular. This work focuses to analyse above threshold surcharge flow. The PIV measurements' condition was chosen as 4 dm3/s of inflow through the inlet pipe and a water level of 310 mm, which resulted in a surcharge level of 235 mm (s) at the manhole centre, making surcharge to manhole diameter ratio (s/Φm) 0.98. This was well over the estimated threshold surcharge level. Initial inspection showed that the measurement zone has two different distinctive velocity characteristics. One part of the measuring plane is approximately in line with the inlet–outlet pipes and is characterized by a fast, slowly expanding jet flow. The second part is outside the jet flow zone and is characterized by a recirculation in which the velocity magnitude is around 10% of the jet flow. For these two distinctive velocity zones, data were taken at different image time separation intervals ranging from 250 ms to 4,000 ms, so that velocities of both zones can be estimated accurately using the PIV cross-correlation algorithm.
For seeding, 100 μm polyamide 12 particles were chosen (density = 1,010 kg/m3). The particles were mixed with water and kept in a seeding tank with continuous circulation so that they remain in suspension. The particles were pumped from the seeding tank at a constant rate into the inlet pipe approximately 40Φp upstream on the manhole so that they were well mixed before entering the manhole.
The seeding rate was adjusted by checking the PIV images to facilitate at least five particles in any selected interrogation area. Data were recorded at three vertical planes; one passing through the central axis of the inlet-outlet pipes and the other two at 50 mm offset from it (see Figure 1). Each data set was measured for 5 min, at a rate of 8 image pairs per second, totalling 2,400 pairs of images.
The data were analysed using Dantec Dynamics' DynamicStudio v3.31 software. The collected data were pre-processed after masking the area of interest. The fluid velocity was calculated using an adaptive cross-correlation technique keeping an interrogation area of 128 × 128 pixels with consideration of 50% overlap between two adjacent areas. Median correction post-processing was applied to remove erroneous vectors and which removed approximately 2 to 8% vectors from each measurement set.
Due to the resolution and positioning of the laboratory setup, neither of the cameras was able to cover the whole manhole height. Emphasis was given to the incoming jet to see how velocity is distributed over the length of the manhole. Hence, the data were recorded covering the lower zone of the manhole, from the manhole bottom until the height of 150 mm of the manhole, which is two times of the inlet-outlet pipe diameter.
The open source CFD model tool OpenFOAM®v4.1 was used in this work. The solver interFoam is selected as it can predict the velocity patterns and the free-surface for sharp interfaces. This solver uses a single set of Navier-Stokes/Reynolds-Average equations where the velocity is shared by both phases and a volume of fluid (VOF) method (Hirt & Nichols 1981) captures the free-surface position. The length of the inlet pipe was chosen as 1,000 mm (more than 13Φp) based on some other previous works (Lau 2008; Stovin et al. 2013), and the outlet pipe was kept as 400 mm, which is until the position of the pressure sensor at the downstream of the manhole (Figure 2(a)). The computational mesh for the simulation was prepared with hexahedral Cartesian mesh using cfMesh (Juretić 2015).
Some pre-analysis of CFD simulations were performed in order to test the mesh independence. For this case, three computational meshes were constructed having dx = 2 mm (Mesh 1), 4 mm (Mesh 2) and 6 mm (Mesh 3), respectively, keeping the global refinement ratio as 3. Number of cells at these meshes are: 2.4 million, 861,500 and 380,000, respectively. The inlet flow condition was prescribed as constant discharge of Q = 4 dm3/s. The meshes were simulated using k-ɛ turbulence model and the velocity profiles at the manhole centre and at the outlet pipe were extracted from the results (Figure 2(c)). The mesh analysis was performed applying Richardson extrapolation (Celik et al. 2008). The meshes gave similar results at the manhole jet zone and at the pipe. However, the velocity profiles showed different results closed to the manhole water surface (at around z = 0.29 m to 0.31 m) and close to the bottom (around z = 0 m). Mesh 3 predicted slightly slower velocity at the near surface zone. The apparent order (p) and the Grid Convergence Index (GCI) was calculated at each grid point of the meshes. The average value of p at the manhole centre and pipe were found 2.76 and 2.32, respectively. The GCI values were found higher close to the surface and the walls. Analysis showed that 50% cells in the manhole has GCI value below 10% when comparing Mesh 2 and Mesh 3; while 65% of cells showed below 10% GCI, in the case of comparing Mesh 1 and Mesh 2. However, in case of results at the pipe, 70% and 76% cells showed GCI value below 10% in case of comparing Mesh 3-Mesh 2 and Mesh 2-Mesh 1 respectively. It was apparent that the results go towards mesh independence and Mesh 1 and Mesh 2 show almost similar results. Average approximate relative error between Mesh 2 and Mesh 1 was found to be 2.7%, whereas the simulation time requirement for Mesh 1 was more than three times to that of Mesh 2. Considering the accuracy level and computational time required, Mesh 2 with dx = 4 mm was found to be best suited for this work (y+ is around 5).
Standard wall function was considered in k-ω SST model only as according to the model description, it requires wall function for y + >1. The rest of the models require wall function when 30 < y + <300.
During the simulations, the adjustableRunTime code was used, keeping the maximum Courant–Friedrichs–Lewy (CFL) number at 0.95. Cluster computing system at the University of Coimbra was used to run the simulations using MPI mode. Each simulation was run for 65 s. The steady condition was ensured by checking the residuals of p, alpha.water, k, epsilon and Rxx (where applicable). The first 60 s were required to reach steady state condition and the results of 101 velocity profiles in the last 5 s were averaged and all the numerical analysis were made using averaged data of these mentioned time step results.
RESULTS AND DISCUSSION
The processed PIV velocity data were compared with velocity data of the manhole model. The PIV measurement was taken at the central vertical plane (CVP) along with the left vertical plane (LVP) and right vertical plane (RVP) (Figure 1). However, the position of LVP was much further from the camera, and light ray coming from the two edges of the LVP have to travel through the manhole inlet–outlet pipe. Although the scaled manhole model is made of transparent acrylic, the joints between the model components were semi-transparent to non-transparent. Due to this limitation, the LVP image edges had bad data in some cases, as the camera could not see the seeding particles due to the obstruction made by the model joints. So in this work, PIV data comparison is only done to the CVP and RVP of the manhole. Figure 3 shows the axial (Vx) and vertical (Vz) velocity components comparison between the PIV data and numerical model results. All the comparison is shown as dimensionless velocities as a ratio of the average inlet velocity (Vavg), where Vavg is the ratio between inlet discharge (Q) and pipe cross section (Ap). Both PIV data and CFD data are showing temporal mean velocities from the measurements and simulation results. PIV data at the CVP near the manhole inlet and outlet pipe were not collected as it was not visible clearly by the cameras due to the joints between pipe and manhole.
It can be seen from the velocity comparison at the CVP of PIV (1st and 3rd column of row 1 at Figure 3), that the jet flow in the experimental results starts expanding slowly as it proceeds from the inlet towards the outlet. The CVP axial velocity (Vx) (1st column of Figure 3) reaches up to 110% of Vavg at the experimental results. It dampens down at the manhole centre. The maximum velocity near the outlet pipe is the same as Vavg. However, in the CFD data, this damping effect is not seen. The high-velocity core stays as 110% of Vavg until it reaches the outlet. While the jet flow trying to escape through the outlet pipe, it hits the manhole wall at the top of the outlet (near x = −100 mm and y = 75 mm) and creates a vertically upward velocity component which is up to 20% of Vavg (3rd column of Figure 3). However, in all the CFD, the vertical velocity component at this zone reaches up to 30% of Vavg. Comparing the velocity contours of Vx and Vz at the CVP, it can be seen that LRR model could not predict the velocity profile properly. This model shows different axial velocity contour at CVP which is dissimilar to PIV measurements.
Comparing the velocity towards the pipe axis at the RVP (2nd column of Figure 3), RNG k-ɛ and realizable k-ɛ models slightly underestimate the axial velocity component compared to the PIV experiment. The PIV data showed axial velocity up to 0.3 Vavg, whereas, in these two numerical models, the highest velocity is found 0.25 Vavg. The axial velocity at this plane is properly estimated by k-ω SST model. The shape of the velocity contour is almost similar to that of the PIV data. The vertical velocity component (Vz) (4th column of Figure 3) at PIV is observed between −0.2 Vavg and 0.2 Vavg. The numerical models show similar results of Vz near the outlet of the manhole. However, at the upper part of the measuring plane near the inlet, the Vz was measured in PIV as around −0.2 Vavg, which was not predicted by the numerical models. Only k-ω SST model in this plane shows a negative (downward) velocity near the inlet, however, still underestimated than the PIV.
The velocity comparison was not done for the plane component (Vy) as at CVP, this component is very close to zero and considered not significant (Figure 4).
The standard deviation of the data was also compared at both CVP and RVP for axial and vertical velocity components (Figure 5). The comparison was also made on dimensionless velocities as a ratio of Vavg.
It can be seen that for axial velocity component (Vx) at CVP and RVP (1st and 2nd column at Figure 5), the PIV data show high standard deviation near inlet, outlet and jet expansion zones. High standard deviation can also be seen near the manhole floor of CVP. Vertical velocity component (Vz) showed lower deviation compared to Vx. The standard deviation was found to be significantly low at both planes in CFD results compared to PIV data. All CFD results show marginally higher standard deviation values close to the jet expansion zone. The realizable k-ɛ model shows the minimum fluctuation in the velocity fields, resulting almost zero standard deviation. The k-ω SST model shows the highest standard deviation among all the four RANS models, still, the value is lower than that of the PIV data. As these numerical models are formulated from RANS equations, a big part of the turbulence variabilities is averaged out from the results already. Perhaps Large Eddy Simulation (LES) model could show better standard deviation match with the PIV data. Moreover, for this work, only 5 s of numerical simulation data were considered. In case the longer period of data was taken into consideration, probably the standard deviation in the CFD results would come higher. Turbulence coming from pipe rather than the standard turbulence inlet conditions could also be a reason for the differences in numeric.
Examining the velocity contours at Figures 3 and 5, it can be stated that k-ω SST model creates the closest approximation of the manhole velocity field followed by RNG k-ɛ model. The LRR model overestimates the axial velocity component in this plane and the velocity contour is significantly different than the PIV measurement.
To understand the manhole flow more clearly, streamline and water level fluctuation was analysed from the numerical model results (Figure 6). It can be seen that the inflow jet originating from the inlet, marginally spread in the manhole (termed as the diffusive region). When exiting the manhole through the outlet pipe, the diffusive region impinges the manhole wall and, after that, part of the jet flow moves vertically upward. This flow region reaches the manhole water surface, starts moving opposite to the core jet direction and makes a clockwise circulation. Due to this constant circulation, the free surface of the manhole fluctuates slightly. The fluctuation was found more towards the inlet of the manhole. In this case, the water level was observed varying between 0.295 m to 0.320 m. The fluctuation level was found different in the four numerical models. In RNG k-ɛ model, the water level varies from 0.300 m to 0.310 m. While the water level fluctuations in the other three manhole models were found as realizable k-ɛ: varying between 0.300 m and 0.305 m; k-ω SST: varying between 0.300 m and 0.310 m and LRR: varying between 0.310 m and 0.320 m. In the experimental work, the average water level at the manhole was observed as 0.310 m. The fluctuation of the manhole level was not possible to measure in the experimental work as it would require installing a pressure sensor at the manhole bottom, which would create an obstacle in the laser ray path line and hence was not used.
The pressure distributions in different CFD models were also compared with experimental data and can be seen in Figure 7. All distances showed at the horizontal axis are measured from the manhole centre. The inlet and outlet pipes are connected at distance of 0.12 m and −0.12 m, respectively. The bottom pressure in between these two distances also represents free surface water level inside the manhole. Two box plots represent pressure data recorded during the experimental measurement using pressure sensors. The left box plot shows the pressure at the outlet pipe, whose average value was used to generate boundary condition of the numerical models. As the downstream pressure is the same for all the four models, all the model results pass through the average value of this box. Pressures at the upstream side were calculated from the models. All the line plots are showing maximum, minimum and average bottom pressure from each of the four CFD models. The circular marker at the manhole centre (at x = 0 m) is showing the product of recorded average water height at the manhole during the experimental works (h = 0.310 m), water density (ρ) and gravitational acceleration (g).
From Figure 7, it can be seen that when the upstream flow through the inlet pipe enters the manhole, the flow experiences a pressure drop, which is due to the expansion of the flow. The second pressure drop can be observed when the flow exits the manhole and enters the outlet pipe. This drop is due to the flow contraction and much bigger in magnitude. The average pressure line of RNG k-ɛ, realizable k-ɛ and k-ω SST models produce a similar pressure pattern throughout the computational domain. The LRR model overestimates the bottom pressure of the manhole, although the difference found is in the range of a few millimetres of water column head. The maximum, minimum and the average bottom pressure from RNG k-ɛ model shows almost the same line which presents that the RNG k-ɛ model shows almost zero pressure fluctuation. Apparently, the LRR model overestimated the pipe loss in compared to the other three models. This results in an overestimation of upstream pressure at the inlet direction.
As the numerical model reached steady state before extracting any results, and as both inlet and outlet pipes were full, the temporal averaged velocity at each pipe can be considered equal. In this case, a difference in bottom pressure would give the same value as head loss. To compute the pressure drop at the manhole centre for a certain CFD result, each line showing the average bottom pressure (Figure 7) was projected to the manhole centre from both inlet and outlet pipe. The vertical difference of pressure value between these two lines at the manhole centre gives the value of pressure drop for the manhole, which is later divided by and considered as .
The value of manhole head loss coefficient has been reported in different literature. This is directly related to the structural mould types of the manhole, manhole to pipe diameter ratios, as well as manhole surcharge ratio. A comparable analysis of the coefficients of head loss (K) with those of the values reported in different literature are shown in Table 1. It should be noted that only the research works reporting the same manhole mould type are considered here. In the work of both Arao & Kusuda (1999) and Lau et al. (2008), authors reported high head loss coefficient at below threshold surcharge conditions and comparably lower coefficient at above threshold surcharge condition. As in this research, the manhole surcharge condition is comparable to above threshold surcharge, only the coefficient range covering this condition are shown.
|Works done .||Head loss coefficient .||Experimental condition (Φm/Φp) .||Surcharge ratio range|
|s/Φm .||s/Φp .|
|Arao & Kusuda (1999)||0.18–0.58||3.60||0.55–1.67||2–6|
|Lau et al. (2008)||0.28–0.69||9.08||0.65–0.82||5.90–7.45|
|Works done .||Head loss coefficient .||Experimental condition (Φm/Φp) .||Surcharge ratio range|
|s/Φm .||s/Φp .|
|Arao & Kusuda (1999)||0.18–0.58||3.60||0.55–1.67||2–6|
|Lau et al. (2008)||0.28–0.69||9.08||0.65–0.82||5.90–7.45|
Table 1 shows that the four models calculate the manhole flow differently and, hence, give different values of manhole head loss coefficient. It is reported by different authors that the head loss coefficient becomes higher when manhole to pipe diameter ratio (Φm/Φp) is high and vice versa (Bo Pedersen & Mark 1990; Stovin et al. 2013). The manhole reported at Arao & Kusuda (1999) has a similar Φm/Φp. Comparing the findings from the literature, it is apparent that the realizable k-ɛ model gives rather low head loss coefficient for this case. From the remaining three models, the coefficient given by k-ω SST model is almost 50% more than that of the value given by RNG k-ɛ, however, both values lie within the range specified by other researchers.
In this work, two-dimensional three component (2D3C) stereo PIV measurement was done on a scaled inline manhole with manhole to pipe diameter ratio of 3.20, in order to evaluate CFD model constructed in OpenFOAM® and four different RANS models with VOF method. From the analysis, it can be apparent that each model calculates the velocity inside manhole differently. Comparison with PIV measurement at the CVP showed similar velocity as compared to all the numerical models. However, comparison of the velocity at another vertical plane 50 mm offset to the centre, showed that all the CFD models slightly underpredicts the axial velocity. The velocity and locations of vortex structures centres were found marginally different among the models. The k-ω SST model showed the closest approximation of velocity contour followed by the RNG k-ɛ model. All the models showed very good approximations of the average water surface level at the manhole. However, the LRR model could not quite capture the velocity profile in compared to PIV data. This model predicts considerably higher axial velocities for CVP as compared to PIV. Nonetheless, the temporal standard deviations of the axial and vertical velocity components were found significantly low when compared to those of experimental measurement through PIV. As a RANS model is formulated based on time-averaged turbulence data, which could be the reason of having a lower standard deviation in the CFD model.
Bottom pressure analysis through the computational domain shows that the average pressure line is almost similar at RNG k-ɛ, Realizable k-ɛ and k-ω SST models. However, the comparison could be made with data from only two pressure sensors installed at the inlet and outlet pipe, respectively. RNG k-ɛ model showed almost no pressure fluctuation while maximum pressure line predicted by the realizable k-ɛ model was found much higher than the measurement. The calculated head loss coefficients were compared with the values reported in the literature. It was seen that the realizable k-ɛ model shows much lower value compared to the values reported.
Considering all the aspects of the four models analysed here, it can be said that both RNG k-ɛ model and k-ω SST models give a very good approximation of manhole hydraulics. However, it should be noted while using k-ω SST model, the wall boundary cell size must be made considerably small for the proper formulation of the model.
The work presented is part of the Quantifying Uncertainty in Integrated Catchment Studies (QUICS) project. This project has received funding from the European Union's Seventh Framework Programme for research, technological development and demonstration under grant agreement No. 607000. The laboratory facility was made possible through EPSRC, project EP/K040405/1. All the numerical results here showed were performed on the Centaurus Cluster of the Laboratory for Advanced Computing of University of Coimbra, Portugal.