## Abstract

An analytical solution was proposed for the groundwater flow through a defective pipe, which can be used to estimate the water flow rate into the pipe and predict the pore water pressure distribution in surrounding soils. This analytical solution was verified by comparing with experimental results, and the predicted pressure distribution around the defective pipe is proved to be consistent with numerical simulations using the finite element method. From the parametric analysis, the infiltration rate increases as the defect position changes from top to bottom on the pipe, and the effect of defect position is not significant if the water head above the defect is 10 times greater than pipe radius. An approximated solution for estimating the groundwater flow infiltration rate through a circular orifice on the pipe is proposed as well. From the verification and parametric studies, this proposed analytical solution is proved to be an efficient approach for the estimation of groundwater infiltration through a defective pipe.

## INTRODUCTION

Defects can develop in sewer pipes because of the pipe deterioration. Groundwater infiltration through the defect into sewer pipe will increase the cost of wastewater treatment, and may also lead to soil erosion around the defective sewer pipe. Sinkhole or ground collapse frequently occurred as a result of the gradual soil loss with serious consequences (Davies *et al.* 2001; Guo *et al.* 2013a). Furthermore, the untreated groundwater infiltration may result in the contaminant intrusion if the defective pipe is drinking water pipe with negative pressure (Yang *et al.* 2016). Therefore, an effective method to estimate the groundwater infiltration rate into a defective pipe is necessary.

Traditional methods of infiltration rate estimation through the defect on a sewer pipe are developed using the simplified orifice flow model without considering the effect of soil around the defective pipe (Karpf *et al.* 2011). Theoretically, groundwater flows in soils under the hydrostatic water pressure, but water will be infiltrated into the defective pipe due to the change of pressure boundary condition at the defect. The water flow in porous media is driven by the pore water pressure gradient, while the relationship can be linear (Darcy's law) or nonlinear (Ergun's law or other models). The parameters in these models are determined by surrounding soil properties. Guo *et al.* (2013b) proposed a two-dimensional analytical solution for the infiltration rate estimation based on Darcy's law. From the parametric studies, the water infiltration rate through the pipe defect increases as the defect size or the pipe size increases. Yang *et al.* (2016) developed an analytical model by combining the nonlinear relationship, Ergun's law, with an assumed ‘negative pressure zone’, which was verified by comparing with the experimental results.

The governing equation is too complicated to be solved with the combination of the nonlinear model with mass conservation. The governing equation becomes the Laplace equation by combining Darcy's law with the mass conservation, which has an analytical solution for some specific boundary conditions. It has been shown that it is reasonable to simulate groundwater flow using Darcy's law (Bear 1972), and the complicated nonlinear model will make it difficult to determine relevant parameters. Verruijt (1997) first proposed a conformal mapping method using Mobius transformation to change the semi-infinite domain around a circular pipe to the circular complex plane. After that, this conformal mapping method was used by many researchers to solve the seepage issue around a drained circular tunnel (El Tani 2003; Kolymbas & Wagner 2007; Park *et al.* 2008; Huangfu *et al.* 2010). The analytical solution of the Laplace equation can be obtained for the regular shape domain, and the solution for the original domain can be obtained by the reverse transformation.

Current exact analytical solutions were derived for either zero-pressure or constant total water head boundary condition at the circular shape boundary. It has been found that the leakage on the tunnel commonly occurred, and different leaking positions may lead to different water pressure distributions in the surrounding soils (Zhang *et al.* 2015). Because of the complexity of the boundary conditions on the leaking tunnel or pipe, the exact solution is difficult to be obtained. Guo *et al.* (2013b) simplified the problem by assuming that an equivalent circumference was located at the defect, and the perimeter of this circumference was equal to the size of the defect. This approximated solution can predict the water infiltration rate considering the effects of leaking positions, defect size and soil properties. The water flow domain above the pipe defect can be assumed as a sphere from the numerical simulations by Collins & Boxall (2012). Yang *et al.* (2014) conducted experiments to investigate the water infiltration rate through a circular orifice, and the effects of soils, orifice size and flow Reynolds number were analyzed.

Previous studies provided an approach for infiltration rate estimation, while pore water pressure distribution around the defective pipe is not clear. It is necessary to predict the water pressure distribution accurately, as the mechanical behavior of the defective pipe is affected by the surrounding pore water pressure. In this paper, an analytical model is proposed based on the pressure distribution prediction since the water flow velocity in porous media should be proportional to the pressure gradient from Darcy's law. This analytical approach will be verified in comparison with experimental results. To verify the prediction of pore water pressure, numerical simulations were conducted. The effects of defect size and defect position on the infiltration rate are investigated from the analytical and numerical results. Furthermore, this analytical model is extended to the infiltration rate estimation of groundwater flow through a circular orifice on a pipe, and the analytical results are compared with the experimental results by Yang *et al.* (2016).

## ANALYTICAL METHOD

*ϕ*is the total water head, which is given by the sum of the pressure head

*p*/

*γ*and the elevation head

*y*, and

*p*is the pore water pressure, and

*γ*is the unit weight of water.

*ϕ*

_{out}=

*h*where

_{w}*h*is the water layer height, while the boundary condition at the defect is

_{w}*ϕ*

_{in}=

*y*where the pressure head is zero. Because the boundary conditions at the defective pipe are characterized by high complexity, the exact solution is hard to be obtained. Since the water flow is driven by water pressure gradient, it will be possible to estimate the infiltration rate through the defect if the pore water pressure distribution can be accurately predicted. As shown in Figure 1, an assumed zero-pressure boundary was located at the center of the defect, and this assumed boundary is slightly larger than the real zero-pressure boundary in a circular shape. The assumed zero-pressure boundary is a circle with a radius

*r*′ =

*D*/2, where

*D*is the defect size. After introducing the assumed zero-pressure boundary, the burial depth of zero-pressure circle

*h*′ and the radius of zero-pressure circle

*r*′ can be obtained from the geometrical relationship: where

*h*is the burial depth of the defective pipe;

*r*is the radius of the defective pipe;

*α*is the defect position angle;

*D*is the defect size;

*β*is the angle indicating defect size and is equal to

*D*/

*r*.

*et al.*2010). If neglecting the effect of the impermeable zone inside the pipe on the water head distribution, the analytical solution by Huangfu

*et al.*(2010) can be used: Based on the total water head determined from Equation (4), the pressure distribution can be simply calculated:

*p*(

*x*,

*y*) =

*γϕ*(

*x*,

*y*)-

*γy*.

*q*, through the defect under two-dimensional conditions can be calculated by: where

*k*is the permeability of soil in the groundwater flow domain;

*θ*

_{1}and

*θ*

_{2}are determined as

*θ*

_{1}=

*π*/2-

*α*,

*θ*

_{2}=

*π*/2 +

*α*as shown in Figure 1. From Equations (4) and (5), the total water head

*ϕ*can be expressed in terms of

*ρ*and

*θ*, while the flow rate can be determined using Equation (6). After integration and simplification, the flow rate

*q*is expressed as:

In the derivation of this analytical method, the soil domain is assumed to be saturated. For the scenario that water table is below the ground surface, this proposed analytical method can be applied to the saturated zone neglecting the effect of the unsaturated zone, and the soil surface is assumed to be located at the water table.

## RESULTS AND DISCUSSION

### Verification of the analytical method

Experiments under two-dimensional conditions using the same setup as Tang *et al.* (2017) are conducted to verify this analytical solution, and the schematic of experiment setup is shown in Figure 2(a). The model box is constructed with a dimension of 500 × 80 × 500 mm (length × width × height), and the defective pipe is 50 mm in diameter with a 3-mm slot at the top of the pipe. Coarse sand with a mean particle size *d _{p}* = 1.52 mm is used, and the uniformity coefficient

*C*=

_{u}*D*

_{60}/

*D*

_{10}is 1.31 indicating a relatively uniform distribution of sand particles.

*D*

_{60}is the particle size such that 60% of the distribution is finer than this size, and

*D*

_{10}is the particle diameter such that 10% of the distribution is finer than this diameter. The specific gravity of sand particles is 2.6, and the sand layer is maintained fully saturated during the experiment with the porosity of 0.4. The water flow rate is measured using time/weight method by collecting outflow. In this experiment, the sand layer height

*h*is kept at 0.25 m while the water layer height

*h*is increased from 0.01 m to 0.22 m.

_{w}*C*is the shape factor of the granular material;

_{S}*T*is the tortuosity factor. In porous media containing approximately uniform pore sizes,

*C*has a value of about 2.5 and

_{S}*T*is about .

*S*is the surface area per unit volume of soil solids, which is 6/(

_{S}*φ*) and

_{s}d_{p}*φ*is the sphericity of sand particle;

_{s}*μ*is the dynamic viscosity of water and can be taken equal to 0.001 Pa·s;

*e*is the void ratio of sand, which is calculated from the porosity

*n*,

*e*=

*n*/(1-

*n*). Although Equation (8) provides an approach to estimate the sand permeability, coefficient

*φ*is difficulty to obtain.

_{s}To accurately determine the sand permeability in this test, the numerical simulation was firstly conducted using the finite element method (ABAQUS 6.14), and the numerical setup to simulate the experiment is shown in Figure 2(b). In this numerical setup, the boundary condition at the model surface sets the pressure equal to the constant total water head, which corresponds to *h _{w}* if the model surface is selected as the datum plane. At the pipe defect, the pressure head is zero, which indicates that the total water head is equal to the elevation head -

*y*. The permeability is calibrated to be 0.011 m/s in comparison with the measured water flow rate. From back-calculation using Equation (8), sphericity

*φ*of this sand was found to be 0.71, which is a reasonable value from the estimation of sphericity (Clayton

_{s}*et al.*2009) in comparison with the sand particle image in Figure 2(a).

Using this analytical approach, the infiltration rates under various water layer height *h _{w}* are shown in Figure 2(c) in comparison with the experimental and numerical results. The results using the analytical method by Guo

*et al.*(2013b) are also plotted in Figure 2(c) as a dotted line for comparison purposes. The error bar in Figure 2(c) indicates the discrepancy between the experimental results and analytical calculations using this proposed method. It is shown that the result using this proposed analytical method agrees well with the experimental and numerical results, only being approximately 15% lower than the experimental result. Estimations using the analytical model by Guo

*et al.*(2013b) are approximately 50% larger than experimental measurements. The difference between the analytical and numerical results is likely due to the prediction of water head distribution. The method by Guo

*et al.*(2013b) assumed an equivalent circumference at the defect center, and the perimeter of this circumference is equal to the defect size. This equivalent circumference is too small compared with the defect size, and the pressure gradient around the defect is overestimated leading to the difference between analytical predictions and experimental measurements. Although there is a slight difference between the experimental measurements and this analytical estimation, the proposed analytical method still shows acceptable accuracy (15%) in Figure 2(c).

The discrepancy between water infiltration analytical and numerical results can be attributable to the assumption of the semi-circular zero-pressure boundary in this analytical model. This assumption will decrease the pressure gradient at the pipe defect because of the increase of zero-pressure boundary size, and the water drainage at both sides of the defect also decreases the pressure gradient. Because the water flow is proportional to the pressure gradient in Darcy's law, this assumption will lead to the underestimation of infiltration rate comparing with the numerical and experimental results in Figure 2(c). Additionally, the analytical model is developed based on the assumption of semi-infinite soil domain. The dimensions of the numerical setup are the same as those of the experimental setup, restricted by the left and right boundary of the studied area.

The water infiltration rate through the pipe defect is mainly affected by the soil permeability, which is dependent on the soil particle size and the soil porosity. As either the particle size or the soil porosity increases, soil permeability increases, resulting in the increase of infiltration rate. As expected, the particle size of 1.2 mm and 1.8 mm, and the porosity of 0.35 and 0.45, will lead to the higher differences of permeability between the analytical and experimental results, as shown in Figure 3.

The pore water pressure distribution around the defective pipe can be calculated using Equation (4), and the contours with different defect positions are shown in Figure 4 in comparison with numerical simulations. If the defect is at the top of the pipe, angle *α* has a value of *π*/2 as shown in Figure 1, while *α* is 0 for the horizontal defect and −*π*/2 for the bottom defect. In Figure 4, the left half of the contour for each plot is from numerical simulation while the right half is from the analytical prediction, and the analytical results are essentially consistent with numerical simulations.

As shown in Figure 4(a), the pore water pressure obtained from the numerical simulation changes significantly close to the defect, while the pressure change in the analytical prediction is slightly small. From the numerical simulations in Figure 4(a), the water pressure is zero at the pipe defect, and the pressure is about 1 kPa at the zone above the defect (−0.2 m < *y* < 0.15 m). On the other hand, the water pressure from analytical results is less than 1 kPa in this area (−0.2 m < *y* < 0.15 m), which indicates that the pressure slowly changes close to the defect, while the numerical simulation shows a rapid change. This discrepancy is because of the assumed zero-pressure semi-circle in the analytical model, which expands the zero-pressure boundary and results in the decrease of pressure gradient close to the pipe defect. This discrepancy is also shown in Figure 4(b) with a horizontal defect on the pipe, while the difference is very small in Figure 4(c) with a bottom defect. The difference in pressure distribution indirectly indicates the difference of flow rate estimation in Figure 2(c). In this proposed analytical model, the assumed zero-pressure boundary is a circle at the defect, while pore water pressure can actually only be dissipated through half of the circle. Therefore, the predicted pressure gradient close to the defect is less than that in the actual condition.

Although there exists some difference between this analytical prediction and experimental measurements or numerical results, the difference is still acceptable from the groundwater infiltration rate estimation in Figure 2(c) and the pressure distribution prediction in Figure 4. In Figure 2(c), the estimated flow rate is about 15% lower than the experimental measurements. The predicted pressure distribution around the defective pipe in Figure 4 is similar to the numerical simulation except for the position adjacent to the defect. Guo *et al.* (2013b) stated that the defect angle *β* was normally less than *π*/18, otherwise, the surrounding soil would be eroded into the pipe and eventually lead to the pipe failure. If the defect is small, this analytical method will cause less error since the difference as shown in Figure 4 is only at the position close to the defect. In comparison with the experimental and numerical results, this analytical method can effectively predict the infiltration rate through the defect and water pressure distribution around the defective pipe.

### Parametric analysis

Davies *et al.* (2001) reported that the defect on a pipe could occur at various positions. The effect of defect position on groundwater infiltration rate is investigated using this analytical method and numerical simulation. As shown in Figure 5(a), the infiltration rate is slightly increased from the top (*α* = *π*/2) to bottom position (*α* = −*π*/2), and the analytical estimation agrees well with the numerical results. Studies by Zhang *et al.* (2015) on the leaking tunnel showed the leakage at the bottom of a tunnel (tunnel invert) would cause the most serious consequences. The water head Δ*h* above the defect increases as the defect changes position from top to bottom, and the effect is not significant due to the small size and the large burial depth of the pipe. In Figure 5(b), the infiltration rate with various defect positions is compared for different burial depths of pipes where *q*_{T} and *q*_{B} are the infiltration rate through a top and a bottom defect, respectively. The ratio *q*_{T}/*q*_{B} is close to one when the burial depth is more than 10 times the pipe radius, which indicates the small effect of defect position. On the other hand, for shallow buried pipes, the effect of defect position can be significant.

Since this analytical model is developed assuming a small defect, it is necessary to investigate the effect of defect size on the accuracy of this model. As shown in Figure 5(c), the estimated infiltration rate agrees well with the numerical results for the defect size angle *β* from 0 to *π*/9, and the error between the analytical and numerical results is less than 15%. Therefore, the proposed analytical model is proved to be accurate enough for various defect sizes.

### Application in three-dimensional conditions

*et al.*2001). The proposed analytical model is derived based on two-dimensional conditions and will be extended to three-dimensional conditions with the simplification shown in Figure 6. If the flow domain is still assumed to be semi-infinite, an assumed zero-pressure boundary is located at the center of the defect, and the assumed boundary is a small sphere rather than the circle in two-dimensional conditions. If neglecting the impermeable effect of the pipe in the longitudinal direction, the total water head at any cross section through the center of this sphere can be simplified to be satisfied by Equation (4). Therefore, the infiltration rate

*Q*in three-dimensions can be calculated using Equation (9) considering the integration along the half spherical surface: where

*θ*

_{1}and

*θ*

_{2}correspond to the angle indicating the defect size and position, which can be determined using the same equations in two-dimensional conditions (Figure 6).

Yang *et al.* (2016) conducted experiments to study the infiltration rate through a circular orifice at the top of a pipe. Sand with the mean particle size of 0.47 mm was used, and the porosity of sand layer was 0.36. The sand layer height above the orifice is 0.155 m, and water infiltration rate with various water layer height *h _{w}* was measured. Through numerical simulation using ABAQUS 6.14, the permeability of sand is calibrated to be 7.8 × 10

^{−4}m/s. The numerical results agree well with experimental measurements as shown in Figure 7. From back-calculation using Equation (8), the sphericity of sand particle is 0.79. Using this calibrated permeability, the groundwater infiltration rate can be estimated from the proposed analytical method in Equation (9) for the circular defect, and the results are consistent with experimental measurements as shown in Figure 7. In this analytical method for three-dimensional conditions, the assumed zero-pressure boundary is in the shape of semi-sphere. The pipe affects the water pressure distribution in the longitudinal direction, which is neglected in this assumption and may lead to the discrepancy between the analytical and numerical results in Figure 7. Also, this proposed analytical method for three-dimensional conditions is based on the assumption of semi-infinite soil domain, while the dimensions of experimental and numerical setups are restricted in the studied area, which can cause the difference between the results in Figure 7.

## CONCLUSIONS

An analytical method was proposed to estimate the pore water pressure around a defective pipe and the groundwater infiltration rate through the defect on the pipe. This analytical method is developed based on the exact solution of seepage into a circular drained tunnel, and an assumed zero-pressure boundary was proposed to account for the partly drained zone at the defective pipe. This proposed analytical solution was verified by comparing the infiltration rate with experimental results, and numerical simulations were carried out to examine the pore water pressure distribution. From these comparisons, the analytical predictions of infiltration rate and pressure distribution agree very well with the numerical and experimental results. From the parametric analysis using this analytical method, the infiltration rate is found to increase with the defect position from top to bottom position on the pipe, and the effect of defect position is not significant if the water head above the defect is about 10 times greater than the pipe radius. In comparison with the numerical simulations for various defect sizes, this analytical method is valid for defect size up to *π*/9. A modified approach based on this analytical solution is proposed to estimate the infiltration rate through a circular orifice on a pipe, which is verified by comparing with the experimental results by Yang *et al.* (2016). This proposed analytical method provides an effective and efficient approach to estimate the infiltration rate through the defective pipe and water pressure distribution in surrounding soil.

## ACKNOWLEDGEMENTS

This research is supported by Alberta Innovation Technology Futures (AITF) Graduate Student Scholarship, Mitacs Globalink Research Award, and China Scholarship Council (CSC). The authors also acknowledge the support of the Discovery Grant program of the Natural Sciences and Engineering Research Council of Canada (NSERC).

## REFERENCES

*.*