**E**ff**ect of Overburden Height on Hydraulic Fracturing of Concrete-Lined Pressure Tunnels Excavated in Intact Rock: A Numerical Study**

#### **Moses Karakouzian 1, Mohammad Nazari-Sharabian 1,\* and Mehrdad Karami <sup>2</sup>**


Received: 2 June 2019; Accepted: 17 June 2019; Published: 19 June 2019

**Abstract:** This study investigated the impact of overburden height on the hydraulic fracturing of a concrete-lined pressure tunnel, excavated in intact rock, under steady-state and transient-state conditions. Moreover, the Norwegian design criterion that only suggests increasing the overburden height as a countermeasure against hydraulic fracturing was evaluated. The Mohr–Coulomb failure criterion was implemented to investigate failure in the rock elements adjacent to the lining. A pressure tunnel with an inner diameter of 3.6 m was modeled in Abaqus Finite Element Analysis (FEA), using the finite element method (FEM). It was assumed that transient pressures occur inside the tunnel due to control gate closure in a hydroelectric power plant, downstream of the tunnel, in three different closure modes: fast (14 s), normal (18 s), and slow (26 s). For steady-state conditions, the results indicated that resistance to the fracturing of the rock increased with increasing the rock friction angle, as well as the overburden height. However, the influence of the friction angle on the resistance to rock fracture was much larger than that of the overburden height. For transient-state conditions, the results showed that, in fast, normal, and slow control gate closure modes, the required overburden heights to failure were respectively 1.07, 0.8, and 0.67 times the static head of water in the tunnel under a steady-state condition. It was concluded that increasing the height of overburden should not be the absolute solution to prevent hydraulic fracturing in pressure tunnels.

**Keywords:** pressure tunnel; hydraulic fracturing; transient flow; finite element method (FEM); Abaqus Finite Element Analysis (FEA)

#### **1. Introduction**

High-pressure water tunnels, usually excavated in rock, convey water from the upstream reservoir toward the turbines in hydropower plants, and turn the energy of the flowing water into electricity. These tunnels can be either steel-lined or concrete-lined; steel-lined pressure tunnels are safer, but more costly to construct. Concrete-lined pressure tunnels are faced with serious challenges in design and construction. In these tunnels, due to the lining permeability in areas where the rock is jointed, water seeps into the joints and expands the cracks. In cases where the rock is intact and not jointed, high internal pressure can fracture the rock. This phenomenon is called hydraulic fracturing, which makes the rock mass on the tunnel unstable.

Researchers studied various aspects of hydraulic fracturing in pressure tunnels, some of which are mentioned below. In 1986, Schleiss investigated leakage from pressure tunnels in the presence of groundwater, and presented the leakage relationships under two conditions: cracked and non-cracked concrete linings. The author concluded that dynamic pressure caused by various factors, such as a water hammer or earthquake, had a significant impact on the rate of leakage and expansion of rock joints [1]. Fernández and Alvarez (1994) used the finite element method (FEM) to study water leakage from a tunnel and effective stresses in the surrounding rock. The authors considered the ratio of total stress to effective stress in the rock as a safety factor against hydraulic fracturing [2]. Furthermore, using the FEM and a closed-form analytical solution, Bobet and Nam (2007) studied effective stresses around pressure tunnels under steady-state operation. The authors assumed there was no groundwater presence around the tunnel, and presented a new analytical solution for stresses and displacements in the surrounding rock and the lining. Their solution considered the interaction between the lining, surrounding rock, and pore pressure in the rock environment, as a result of leakage through the lining [3]. Moreover, Hachem and Schleiss (2010) evaluated pressure wave velocities in steel-lined pressure tunnels, with and without considering the fluid–structure interactions. The authors showed changes in pressure wave velocity, due to the characteristics of the surrounding rock, the thickness of the lining, and wave frequency [4].

Using the FEM, and in order to optimize the thickness of the concrete lining, Olumide and Marence (2012) simulated a pressure tunnel in a two-dimensional elasto-plastic space in order to study the effects of seepage on the tunnel's concrete lining. The authors investigated cracks in the lining by coupling in situ stresses in the rock with the stresses caused by leakage, and validated their numerical model by the analytical method proposed by Schleiss (1986) [5]. Furthermore, Simanjuntak et al. (2014) studied displacements and stresses in the concrete lining and rock mass using the FEM, considering non-uniform in situ stresses around a pressure tunnel [6]. The authors showed that, in non-uniform stress conditions around the tunnel, when the coefficient of lateral rock pressure (*k*0) is ≤1, cracks form in the crown of the tunnel due to high internal pressure. In 2015, Zhou et al. considered the interactions between fluid and the concrete lining in pressure tunnels in order to study stresses in a tunnel lining. The authors considered the hydro-mechanical behavior of the surrounding rock and the tunnel lining by simulating the fluid inside rock joints and cracks in the concrete lining [7]. More recently, Pachoud and Schleiss (2016) investigated the impacts of the transient pressure wave in pressure tunnels due to water hammer. The authors presented analytical formulas for stresses and displacements in the environment around the tunnel [8]. In 2018, Zareifard presented an analytical solution for the design of pressure tunnels, considering seepage from a tunnel [9]. It should be noted that, in addition to the negative impacts of hydraulic fracturing in water resources and related structures, it is widely used in the oil and gas industry as a novel method for the exploitation of underground reservoirs [10].

Based on the literature review, many researchers investigated cracking in the rock or tunnel concrete lining in steady-state conditions with constant overburden, yet they did not consider the interactive effects of overburden height and transient pressures on hydraulic fracturing in the rock environment surrounding a tunnel. In the present study, in order to fill this gap in knowledge, the impact of overburden height on the bearing capacity of a tunnel in a steady-state condition was studied. For the next step, the effect of transient pressure caused by rapid gate operations on hydraulic fracturing in a tunnel in critical conditions (e.g., load rejection of the powerhouse) was investigated. The Mohr–Coulomb failure criterion was used to investigate the rock failure both in steady-state and transient-state conditions.

#### **2. Materials and Methods**

Using the FEM and considering the environment around the tunnel as an intact rock with elasto-plastic behavior, stresses around a tunnel were studied using the Abaqus Finite Element Analysis (FEA) software [11]. Because cracks were not explicitly modeled in the numerical model, principal stresses were used as indicators of potential fracturing. Since this solution involves the hydraulics of water flow in the tunnel and the impacts of internal pressures on the surrounding rock, the Hammer software and Abaqus FEA were linked for modeling and analyzing the system. The Hammer software [12], which works based on the method of characteristics (MOC), was employed to analyze the changes in the internal pressure inside the tunnel, as a function of time. The hydraulic

analysis results, from the Hammer software in transient-state conditions, were transferred to the Abaqus FEA in order to analyze stresses in the surrounding rock. It was assumed that, at first, the internal pressure was in the steady-state condition. Afterward, because of sudden gate closure, transient pressures built up in the tunnel. The flowchart below summarizes the main steps in this study (Figure 1).

**Figure 1.** The flowchart of this study.

#### *2.1. Governing Equations*

Using the general Navier–Stokes equation (Equation (1)) and simplifying assumptions—(a) fluid inside the tunnel is non-viscous, (b) fluid velocity fluctuations are ignored, and (c) water has linear compressibility—the hydrodynamic pressure inside the tunnel can be obtained using Equations (1) and (2) [13].

$$
\rho \left( \frac{\partial v}{\partial t} + v.\nabla v \right) = -\nabla p + \mu \nabla^2 v + B\_\prime \tag{1}
$$

$$\frac{\rho}{K} \frac{\partial^2 p}{\partial t^2} = \nabla^2 p\_\prime \tag{2}$$

where *v* is the flow velocity (m/s), μ is the fluid viscosity (Pa·s), *B* is the vector of body forces (N), *p* is the fluid hydrodynamic pressure in the tunnel (Pa), *K* is the bulk modulus (Pa), *t* is time (s), and ρ is the water density (kg/m3). The only boundary condition governing the fluid inside the tunnel is the boundary between the water and the lining. This boundary condition is defined as in Equation (3) [13].

$$\frac{\partial p}{\partial n} = -\rho \ddot{u}\_{s\_{n'}} \tag{3}$$

where *<sup>n</sup>* is a perpendicular vector to the shared surface between the water and the lining, and .. *us* is the second derivative of the lining elements' displacements in contact with the fluid (m/s2). The above equation is based on the assumption that the radial displacements of the fluid and lining are similar at the contact surfaces. The connection between the hydrodynamic pressures inside the tunnel {*p*} and forces applied to the model caused by the hydrodynamic pressure {*f*} are linked by matrix [*Q*] (Equation (4)) [14].

$$
\langle [Q] | p \rangle = \langle f \rangle. \tag{4}
$$

The coupling matrix [*Q*] links the hydrodynamic pressures due to gate closure, with forces generated in the concrete lining. This equation is solved by the Abaqus FEA, by applying the boundary and initial conditions. Finally, the displacements and stresses in the rock due to hydrodynamic pressures are calculated. The forces in the concrete lining due to water pressure are transferred to the surrounding rock. As a result, the equilibrium equation in the surrounding rock is as follows [14]:

$$[M]\{\ddot{u}\_r\} + [\mathbb{C}]\{\dot{u}\_r\} + [K]\{u\_r\} = \{R\_r\},\tag{5}$$

$$<\langle R\_{\mathbb{F}}\rangle = \langle f\rangle - \langle f\_{\mathbb{F}}\rangle,\tag{6}$$

where *M*, *C,* and *K* are mass, damping, and hardness matrices, respectively; *Rr* is the vector of external forces on the rock; .. *ur* and . *ur* are the second and first derivatives of the rock elements, respectively; and *ur* is the displacement of the rock elements. Equation (6) shows that part of the hydrodynamic force is resisted by the concrete lining (*f* s). Therefore, the relationship between hydrodynamic pressure and displacements in the surrounding rock is obtained by linking Equations (4) and (5).

#### *2.2. FEM and E*ff*ective Parameters*

In the FEM modeling, the following factors were considered in the simulations:


Using a trial-and-error procedure, the model boundaries were selected in a way so as to not affect stress distribution around the tunnel. Therefore, the model boundaries were set at 6d × 25d (d = tunnel diameter) (Figure 2). More details on the grid and boundary sensitivity analyses are presented in the Appendix A. Table 1 shows the mechanical properties of the lining. It was assumed that the tunnel was located above the groundwater level.

**Figure 2.** The geometry of the surrounding rock and the excavation section.


**Table 1.** Mechanical properties of the lining [15,16].

The only degree of freedom inside a water pressure tunnel is pressure in the fluid nodes. In Abaqus FEA, acoustic elements possess this feature and are appropriate for simulating the fluid movement inside a tunnel. Since acoustic environments are elastic, shear stresses do not exist in these environments and pressure is proportional to volume strain. In addition, the effects of inertia and compressibility are considered in acoustic elements [17].

In this study, the interactions between the fluid, lining, and surrounding rock were considered, and the fluid density was set at 1000 kg/m3, with a bulk modulus of 2.07 GPa. Taking into account the interaction between the fluid and lining, the structure surfaces and fluid that are in contact were tied together in Abaqus FEA. Applying the same reasoning, the external nodes of the lining were tied with the nodes of the surrounding rock due to structure–rock interactions.

The environment around the tunnel has in situ stresses. By defining the height, specific weight, and the lateral stress coefficients of the rock, equilibrium conditions between the horizontal and vertical stresses were achieved, which resulted in zero ground settlement before tunnel excavation.

Since the transient pressure inside the tunnel depends on factors such as gate operations, this pressure was applied to the model as a quasi-static load (by gradually increasing the loading, so that the impacts of inertia became negligible). This means that, at different times, different pressures were applied to the lining. In order to apply the time-variable pressure to all of the nodes in the acoustic environment, the amplitude function in Abaqus FEA was used.

The plane strain elements with linear interpolation function were employed to mesh the surrounding rock and concrete lining. In Abaqus FEA, these types of elements are called CPE4R (continuum plain strain element 4-node reduced integration). Meshing was such that, by approaching the tunnel section, the mesh size decreased to achieve more precise results. Moreover, plane stress family elements were used for meshing the lining. In Abaqus FEA, these elements are called CPS4R (continuum plain stress element 4-node reduced integration). Furthermore, in the lining and rock elements, the reduced integration method was implemented, and AC2D4 elements (acoustic continuum 2-dimension 4-node) were used to simulate water transient pressure inside the tunnel [11,18].

The first stepin this research study was simulatingin situ stresses. In order to have zero ground-surface settlement before tunnel excavation, these stresses should satisfy Equations (7) and (8) [19].

$$
\sigma\_{\mathcal{V}} = \gamma\_I h\_{\prime} \tag{7}
$$

$$
\sigma\_h = \gamma\_r h k\_{0\prime} \tag{8}
$$

where γ*<sup>r</sup>* is the specific weight of the rock, *h* is the height of the overburden rock, and *k*<sup>0</sup> is the coefficient of lateral rock pressure. In the next step, the tunnel excavation was simulated. It was assumed that the tunnel was excavated by a tunnel boring machine. At this stage, the release of stress and displacements occurred in the rock. After the tunnel cross-section, the concrete lining was simulated. In this step, a concrete lining with the specifications presented in Table 1 was modeled. Interaction between the rock and the lining was considered using tangential and normal stiffnesses. The final step was applying steady-state conditions and hydrodynamic water pressure to the lining. It was assumed

that the transient pressure inside the tunnel would occur due to the turbine gate closure in a typical hydropower plant (Figure 3).

**Figure 3.** Simplified hydropower plan.

Using the Hammer software, hydraulic analyses were performed in steady-state and transient-state conditions, by defining the characteristics of the duct (water tunnel), reservoir, and gate closure schedule. The parameters used for these analyses and their values are presented in Tables 2 and 3. Water velocity was considered to be 6 m/s, according to the recommended value by the United States Bureau of Reclamation [20]. Moreover, in order to create critical conditions in the tunnel to determine its bearing capacity, it was assumed that the surge tank was disabled. The maximum flow rate was considered to be 60 m3/s, and, according to the gate closure time in typical water tunnels, the gate closure times were considered to be 14 s (fast), 18 s (normal), and 26 s (slow). Furthermore, based on the duct characteristics, the pressure wave velocity according to the characteristics of the duct was determined by Equation (9), where *D* is the inner diameter of the tunnel, *E* is the elastic concrete modulus, and *tc* is the lining thickness [21].

$$a = \sqrt{\frac{1}{\rho\_{\text{lv}} \left( \left( \frac{1}{K\_w} \right) + \left( \frac{D\lambda}{Et\_c} \right) \right)}}. \tag{9}$$



It was assumed that the lining could expand in the longitudinal direction, which results in a λ (confinement coefficient) value of 1 [21].

The location for investigating stresses in the rock was considered to be near the penstock (Node C) in the tunnel path. Finally, the pressure values at this point, obtained by the Hammer software, were input to the Abaqus FEA, under hydrodynamic pressure loadings.


**Table 3.** Node elevations and properties of node E (turbine).

#### *2.3. Parametric Study in Steady-State Conditions*

In this study, the bearing capacity of the pressure tunnel under a steady-state condition was investigated. The bearing capacity of a concrete-lined pressure tunnel is dependent either on the tensile strength of the surrounding rock only, or both the concrete and the rock. Generally, in pressure tunnels that are excavated in intact rock, the concrete lining does not contribute to the stability of the tunnel during excavation and only provides an appropriate bed for water to be conveyed toward the turbine.

Crack formation in the surrounding rock (hydraulic fracturing), due to high internal pressure, shows that the pressure tunnel reached its ultimate capacity; therefore, a parametric study was performed to investigate the impact of overburden height on the ultimate bearing capacity of the pressure tunnel in a steady-state condition. Figure 4 shows the effective parameters in the bearing capacity of a tunnel. *Pe*(*pi*) is the water pressure on the outer surface of the concrete lining due to seepage. According to Equation (10), this parameter is dependent on internal pressure, and it was calculated using the FEM in the present study [22].

$$P\_{\mathbf{c}}(\boldsymbol{p}\_{i}) = \frac{\mathcal{P}\_{i}}{1 + \left(\frac{k\_{\mathrm{r}} \ln\left(\frac{\boldsymbol{p}\_{i}}{r\_{i}}\right)}{k\_{\mathrm{c}} \ln\left(\frac{\boldsymbol{p}\_{i}}{r\_{\mathrm{r}}}\right)}\right)'} \tag{10}$$

where *t* (thickness of the lining), *k*0, *kr* (hydraulic conductivity of the rock), and γ*<sup>r</sup>* (specific weight of the rock) were considered constant for all analyses, and their values are presented in Table 4. Moreover, the values for *Er*, *C*, ϕ, and *h* are presented in Table 5.

**Figure 4.** Effective parameters in the parametric study.


**Table 4.** Constant values in present study.

**Table 5.** Parameters and their values in parametric analyses [23–25].


In each analysis, the bearing capacity of the tunnel was determined by increasing the internal pressure (*pi*) to the point at which the first crack was formed in the surrounding rock. In this situation, the corresponding *pi* was considered as the failure threshold, which indicates the ultimate bearing capacity of the tunnel.

#### *2.4. Verification of the FEM*

Tunsakul et al. (2014) investigated fracture propagation in a rock mass surrounding a tunnel, under high internal pressure, in an experimental model [26]. In order to verify the numerical model developed in the present study, an experiment with *k*<sup>0</sup> = 0.5 was considered. Other assumptions based on the experimental model were as follows: no concrete lining exists (*t* = 0); the surrounding rock is impermeable (*Pe*(*pi*) = 0). According to Tunsakul et al. (2014), the cracking location depends on the coefficient of lateral rock pressure (*k*0). When *k*<sup>0</sup> < 1, cracking occurs in the crown of the tunnel, and, for *k*<sup>0</sup> ≥ 1, cracks form at the sides of the tunnel. The FEM results were in agreement with the findings of Tunsakul et al. (2014), since cracks formed in the crown of the tunnel (Figure 5).

**Figure 5.** Comparing the formation and location of cracks in the numerical model (**a**) vs. the experimental model by Tunsakul et al. (2014) (**b**).

#### **3. Results and Discussion**

#### *3.1. Changes in Pore Pressure*

Since concrete is not completely impermeable, water passing through the pores penetrates into the surrounding rock. The permeability of rock varies between 1 (completely permeable) and 10−<sup>12</sup> m/s (almost impermeable). Figure 6 shows a comparison between the pore pressures obtained from the numerical model and the value obtained from Bouvard's and Pinto's (1969) analytical solutions (Equation (10)) [22]. According to Figure 6, the pore pressure obtained using the analytical solution is greater than the value obtained from the numerical solution.

**Figure 6.** Numerical vs. analytical results for pore pressures in the rock.

The hydro-mechanical interaction is caused by jointing in rocks due to water pressure. As a result, water penetration toward the rock increases, which consequently leads to more leakage from the tunnel, and increases the pore pressure inside the rock joints. In analyses performed by Bouvard and Pinto, it was assumed that the rock was a porous medium, without effective stresses due to overburden rock. However, in the present numerical modeling, effective stresses due to overburden rock, which reduce the pore pressure in the rock, were considered as important factors in the simulation.

#### *3.2. Bearing Capacity of the Tunnel in Normal Operating Conditions (Steady-State Conditions)*

Based on the mechanical characteristics of the rock, and using the Mohr–Coulomb failure criterion, the failure of elements in the rock environment surrounding the tunnel was investigated. In this regard, maximum and minimum principal stresses were calculated and are presented in Figure 7. Since the problem was symmetric, half of the rock elements surrounding the tunnel were investigated.

In the Mohr–Coulomb failure criterion, the critical line divides the two-dimensional (2D) principal stresses into the safe zone (below the critical line), and the area where stresses cause failure in the element (above the critical line). Due to water pressure on the lining, and the transmission of part of this pressure to the rock, tensile stresses appear around the tunnel. Figure 7 shows that the combined stresses in elements surrounding the tunnel are located below the critical line, indicating that these elements do not fail. In Figure 7, the elements located near the critical line are considered as critical elements, which means higher water pressures will cause these elements to fail. Therefore, this internal pressure is the ultimate value that the rock can endure.

**Figure 7.** The principal stresses in elements surrounding the tunnel, and their locations relative to the critical line.

Considering the ultimate bearing capacity (*Pi <sup>U</sup>*) obtained from the FEM, the Norwegian design criterion (NC) was evaluated to determine the overburden required to prevent hydraulic fracturing, due to the internal static head in the pressure tunnel (Equation (11)).

$$h > \frac{\gamma\_w \cdot H}{\gamma\_{r^\*} \cos \beta'} \tag{11}$$

where *h* is the minimum overburden height to prevent hydraulic fracturing, *H* is the static water height in the tunnel, β is the slope of the valley that is equal to zero in this study, and γ*<sup>w</sup>* and γ*<sup>r</sup>* are the specific gravities of water and rock, respectively. The NC estimated the ultimate bearing capacity of the tunnel as 28 m water (for γ*<sup>w</sup>* = 10 kN, γ*<sup>r</sup>* = 28 kN, and *h* = 10 m).

Figure 8 shows (a) the ultimate bearing capacity of the tunnel under a steady-state condition based on the FEM results and the NC, (b) stress distribution in the rock, and (c) vertical displacements in the rock. According to the FEM results, the internal pressure equal to 8.7 bar (87 m water) is the maximum bearing capacity of the tunnel. According to Figure 8d, at first, the concrete lining cracks, and then cracking occurs in the surrounding rock. Afterward, a parametric study was performed and the effect of overburden height on the bearing capacity of the tunnel was investigated (Figure 9 and Table 6).

**Figure 8.** *Cont*.

**Figure 8.** (**a**) The bearing capacity of the tunnel; (**b**,**c**) maximum stress and vertical displacement contours in the surrounding rock in a steady-state condition; (**d**,**e**) crack formation in the lining and rock.

**Figure 9.** Effect of overburden height on the bearing capacity of the tunnel. (**a**) *h* = 10 m; (**b**) *h* = 20 m; (**c**) *h* = 30 m; (**d**) *h* = 40 m.


**Table 6.** The ultimate bearing capacity of the tunnel (*PU <sup>i</sup>* ), under different scenarios. FEM—finite elementmodel.

For validation of the results, an actual pressure tunnel project (Herlansfoss, Norway) was considered. In this project, the rock surrounding the tunnel is made of non-cracked schist with a friction angle equal to 30◦. At chainages 150–200 m of the tunnel path, the tunnel is horizontal, with overburden height and an internal pressure of 45 m and 136 m water, respectively. In this range, no failure of the rock due to the internal pressure is observed [23]. The *Pi <sup>U</sup>* values based on the FEM results (Figure 9), the NC, and in the actual project, are equal to 138, 126, and 136 m water, respectively, indicating that the FEM results are reliable. According to Table 6 and Figure 9, changes in *Pi <sup>U</sup>* with increasing the overburden height are more dependent on ϕ than *E* and *C*.

In this study, the NC estimated the ultimate bearing capacity of the concrete-lined tunnel with 20 cm of the lining thickness and 3.6 m of diameter in intact rock, less than the value obtained from the numerical model in most cases (Figure 10). According to Figure 10, for all friction angles, up to 30 m of overburden height, the NC estimated *Pi <sup>U</sup>* less than the FEM. By increasing the overburden height to values greater than 30 m, for all quantities of friction angles, *Pi <sup>U</sup>* obtained from the FEM results approached the value obtained from the NC. Therefore, the NC overestimates the overburden required to prevent hydraulic fracturing in pressure tunnels located near the ground surface (*h* ≤ 30 m).

**Figure 10.** Comparing *Pi <sup>U</sup>* based on the finite element model (FEM) results, the Norwegian design criterion (NC), and an actual pressure tunnel project, under steady-state conditions.

#### *3.3. Applying Pressure Fluctuations in Transient-State Conditions*

While the gate is being closed, the velocity upstream of the gate decreases, which results in a pressure build-up in the tunnel. The maximum pressure occurs when the gate is fully closed. Figure 11 shows the pressure fluctuations upstream of the gate during the total time of hydraulic analysis. In this study, it was assumed that the gate closure was done at a constant rate. Under sudden load rejection, which happens when the gate is closed quickly at a constant rate, critical conditions occur in the hydropower tunnel. This often occurs in special situations, such as a turbine generator shut down due to an earthquake. If the gate is closed in several stages, with different closure rates, the maximum transient pressure and the consequent damage to the lining will decrease.

**Figure 11.** Pressure fluctuations upstream of the gate during the hydraulic analysis.

In this study, the maximum pressure that was reached after complete gate closure was considered as the hydrodynamic load in the simulations; therefore, the rock resisted maximum water pressure, and the probability of hydraulic fracturing increased. The pressure–time curves at three different gate closure times are presented in Figure 12.

**Figure 12.** Changes in pressure upstream of the gate, in different scenarios: 14 s, 18 s, and 26 s.

Figure 12 shows that the curves have a rising trend. According to Figure 12, by increasing the gate closure time, the maximum hydrodynamic pressure was reduced; during the first few seconds, kinetic energy converts into pressure energy at a faster rate; during the last seconds of gate closure, as water upstream of the gate reaches a steady-state condition, the pressure energy decreases. Transferring these loads to the Abaqus FEA and employing the Mohr–Coulomb failure criterion, the stresses in the rock elements surrounding the tunnel were calculated for different gate closure times, and are presented in Figure 13. According to Figure 13, the stress in almost all of the elements is

above the critical line, which indicates that these elements failed. When the gate closes in 14 seconds, hydrodynamic pressures impose greater forces on the surrounding rock. Increasing the gate closure time decreases the pressures resulting from the transient-state condition. Figure 13c shows that some elements are in the safe zone, which demonstrates that increasing the gate closure time decreases the possibility of hydraulic fracturing in the rock.

**Figure 13.** Stress in the elements by closing the gate in (**a**) 14 s, (**b**) 18 s, and (**c**) 26 s.

The Effect of Increasing the Overburden Height on Preventing Hydraulic Fracturing

In order to investigate the impact of the overburden rock on hydraulic fracturing, the critical gate closure time of 14 s was considered. In this regard, the failure of the elements in the rock around the lining was investigated by considering the overburden height as 10, 20, 40, and 70 m. Figure 14 shows the stresses in these elements, with respect to the critical line, under different overburdens.

According to Figure 14, by increasing the overburden from 0.13H to 0.26H, stresses in some elements are higher than the critical value; therefore, these elements are located above the critical line. On the other hand, in some elements, stresses are less than the critical value, and they are located below the critical line. By increasing the overburden to 0.26H, as expected, the risk of hydraulic fracturing decreased as more elements appeared below the critical line. By increasing the overburden to 0.53H, more elements appeared below the critical line, but a number of elements were still above the critical line. For the case of an overburden equal to 0.93H, all elements appeared in the safe zone below the critical line, since the elements on the sides of the tunnel section are under the influence of stress concentration. As a result, they appear in the safe zone earlier than the other elements, due to an increase in the overburden.

**Figure 14.** The main stresses in surrounding rock by increasing overburden: (**a**) 0.13H, (**b**) 0.26H, (**c**) 0.53H, and (**d**) 0.93H (gate closure time in all cases was 14 s).

Furthermore, Figure 15 shows the relationship between the overburden and the percentage of fractured elements around the tunnel. According to Figure 15, an increase in the overburden is not a safe solution to prevent the hydraulic fracturing phenomenon. In other words, due to rock overburden, the environment around the tunnel has in situ stresses. If the combination of stresses caused by the transient pressure inside the tunnel and the rock in situ stresses reaches a critical state, the elements will fail. The occurrence of hydraulic fracturing in one element will lead to a critical state in other elements. Therefore, it should be noted that increasing the overburden does not necessarily put all elements in the safe zone.

**Figure 15.** Percentage of cracked elements due to increased overburden height for gate closure times of (**a**) 14 s, (**b**) 18 s, and (**c**) 26 s.

#### **4. Conclusions**

Different aspects of hydraulic fracturing in pressure tunnels were investigated in previous studies, and various strategies were proposed to prevent this serious issue. One strategy is to increase the height of overburden rock. In the present study, using the FEM, and considering the fluid–structure interactions, a pressure tunnel with an inner diameter of 3.6 m was modeled in intact elasto-plastic rock, and the Mohr–Coulomb failure criterion was employed to investigate failure in the rock environment surrounding the tunnel. A parametric study was performed to evaluate the bearing capacity of the tunnel, and the effect of transient pressure on the failure of the rock elements around the tunnel was investigated. Moreover, the applicability of the Norwegian design criterion was evaluated. Hydraulic analyses in transient hydraulic states were performed using Hammer software. The output of this software, including pressure–time values, were used in Abaqus FEA in order to analyze stresses in the rock.

The results showed the following:


In order to prevent hydraulic fracturing in the rock, construction of a grouted zone around the pressure tunnel, along with increasing the overburden height, is recommended. The grout integrates the rock, in case any cracks exist, and improves its mechanical properties. In addition, the grouted zone results in decreased pore pressure in the rock, and reduces the possibility of hydraulic fracturing.

**Author Contributions:** Conceptualization, M.K. (Karakouzian), M.K. (Karami) and M.N.-S.; methodology, M.K. (Karakouzian), M.K. (Karami) and M.N.-S.; software, M.K. (Karami) and M.N.-S.; validation, M.K. (Karami) and M.N.-S.; writing—original draft preparation, M.N.-S. and M.K. (Karami); writing—review and editing, M.N.-S.; supervision, M.K. (Karakouzian).

**Funding:** This research received no external funding.

**Conflicts of Interest:** The authors declare no conflicts of interest.

#### **Appendix A**

#### *Grid and Boundary Sensitivity Analyses*

Based on a trial-and-error procedure, a grid sensitivity analysis was performed in the model, in order to ensure convergence of the finite element equations, and the optimum pattern of meshing was determined. The mesh sizes used near the tunnel cross-section were smaller than those near the model boundaries (Figure 2).

Moreover, a sensitivity analysis was performed on the side boundaries of the model. According to Figure A1, for *x* ≥ 6D (D is the tunnel diameter) the side boundaries do not affect stresses in the rock environment surrounding the tunnel. In order to ensure preventing errors in the model calculations, this distance was considered twice the value obtained from grid sensitivity analysis (Figure 2).

**Figure A1.** Sensitivity analysis on the model boundaries: (**a**) changes in stress in the rock environment, vertically; (**b**) changes in stress in the rock environment, horizontally; *x* is the distance from the tunnel lining, σ*<sup>H</sup>* is the horizontal stress in the rock, σ*<sup>V</sup>* is the vertical stress in the rock, *PH* is the initial horizontal stress in the rock, and *PV* is the initial vertical stress in the rock.

#### **References**


© 2019 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).
