Static and Dynamic Inﬂuence of the Shear Zone on Rockburst Occurrence in the Headrace Tunnel of the Neelum Jhelum Hydropower Project, Pakistan

: Rockburst is an unstable failure of a rock mass which is inﬂuenced by many factors. During deep excavations, the presence of nearby geological structures such as minor faults, joints, and shear zones increases the likelihood of rockburst occurrence. A shear zone has been observed in the headrace tunnel in the Neelum Jhelum Hydropower Project, Pakistan, which has played an important role in major rockburst events in the project’s history. A rockburst is a seismic event that involves the release of a great amount of energy as the dynamic wave radiated from the seismic source reaches the excavation boundary. In this paper, the FLAC 2D explicit numerical code has been used to simulate the dynamic phenomenon of rockburst near the shear zone in a headrace tunnel. The behavior of the rock mass around the tunnel has been studied under both static and dynamic loading. According to modeling results, rockburst signiﬁcantly a ﬀ ected the upper left quadrant of the tunnel similar to the actual failure proﬁle with a depth of approximately 5 m. The dynamic impact of rockburst has also a ﬀ ected the loading conditions of the support system in the adjacent tunnel. This study elucidates one of the most important rockburst controlling factors through numerical analysis and recommends yielding support measures that can withstand the dynamic impacts of rockburst in deep, hard rock tunnels.


Introduction
The rockburst is a dynamic form of damaging phenomenon that is linked with the violent ejection of rock blocks. It is hazardous during the excavation of deep tunnels [1]. Managing this risk in deep and highly stressed environments is challenging. The severity of the damage caused by rockburst varies depending on different parameters, such as excavation depth, stress level, the quality of the rock mass near the excavation, excavation shape, excavation method, geological structures, and dynamic disturbance. The impact of geological structures on rockburst occurrence in deep hydropower tunnels is discussed in this paper. There are different geological structures such as faults, shear zones, joints, dykes, and discontinuities that encounter during deep tunneling [2].
In the literature, there is sufficient evidence that geological structures have a significant influence on rockburst occurrence, and deep excavations in mining are more prone to rockburst when they approach any structural plane. For example, Durrheim et al. [3] studied 21 rockburst events in deep South African mines and showed that regional structures such as faults and dykes are the principal controlling factors of rockburst. Haile [4] described rockburst events in several deep South African mines influenced by fault slip, and he performed mechanistic evaluations and a design of tunnel support systems. Andrieux et al. [5] described fault slip mechanisms in North American mines that triggered a large rockburst. Hedley [6] discussed the influence of faults on rockburst in hard rock mines in Ontario, Canada. Snelling et al. [7] studied the relationship between shear zones, stress, and seismicity in the Creighton nickel mine in Canada and explained how the major shear zones influence stress flow around the excavation. Morissette et al. [8] showed that excavation-induced stresses along shear zones are the main contributing factor in Creighton and Copper Cliff mines in Canada. From these different case studies, it can be concluded that stresses are normally concentrated near the geological structures. When the excavation approaches such a zone, these stresses overcome the normal stress, causing an unclamping of the structural plane, resulting in an enormous amount of energy release due to shear slip along the structure plane. There has been little focus in the literature related to the influence of the shear zone on rockburst, and nobody has studied the dynamic influence of the shear zone on rockburst at great depth. The localized stress concentrations near the shear zone contribute to seismic activity, which further increases the damage zone around the excavation as the seismic wave reaches the excavation boundary.
It has also been recognized that, in recently completed deep civil and hydropower tunnels, rockburst is influenced by these geological structures, which has been discussed earlier. For example, in the Jinping II hydropower project's tunnels, in China, more than 1000 rockburst events have occurred [9], and the most severe rockbursts occurred where geological structures were present nearby. Zhou et al. [10] studied rockburst in deep tunnels of the Jinping-II hydropower station and determined that the structural plane is the controlling factor for rockburst there. Jeon et al. [11] performed tests to study the effect of faults and weak planes on the stability of a tunnel; they found that deformation increases as the distance from the fault to the tunnel decreases. Loew et al. [12] explained that a steeply angled, brittle fault caused a major rockburst during the construction of the Gotthard Base Tunnel through the Alps in Switzerland. During the construction of the Neelum Jhelum Hydroelectric Project (NJHEP) in Pakistan, 879 rockburst events were recorded [13]; only one intense rockburst occurred, near a shear zone, and Naji et al. [14] described the concentration of stresses near that shear zone, which resulted in the devastating rockburst event of 31 May 2015.
The occurrence of intense rockburst in deep tunnels is inevitable when geological structures are present in deep massive rock mass, and it normally has dynamic characteristics. In such conditions, the stability of the underground excavation is critical. In previous work, static numerical modeling using FLAC was done to find only the influence of the shear zone on rockburst [14]. The dynamic effect of the shear zone on rockburst occurrence both near the boundary of a tunnel and on its support system is still unclear. In this paper, a FLAC 2D dynamic numerical modeling has been done to study the mechanism of a rockburst at great depth. The actual field measured parameters have been used during simulation. It is believed that the most intense rockburst event of 31 May 2015 in NJHEP was due to a slip along a shear zone, which resulted in severe damage to the excavation boundary due to its dynamic impact. We have numerically investigated the mechanism of rockburst in the headrace tunnels of the NJHEP, which have been subjected to dynamic loading, and we have also evaluated the dynamic impact of rockburst on the installed support in the adjacent tunnel.

Neelum Jhelum Hydroelectric Project (NJHEP) Overview
NJHEP is located in the Muzaffarabad district of Azad Jammu and Kashmir, Pakistan, as shown in Figure 1. The project is constructed in a geologically young mountain range of the Himalayas [15,16]. This area has complex geological regime, and having major tectonic features such as faults, shear zones, synclines, and anticlines can be seen in Figure 2. Twin headrace tunnels were excavated with two TBMs (tunnel boring machines), under a rock cover ranging from 1000 to 1900 m. These TBM tunnels passed through the Murree Formation, which comprises alternate beds of sandstone, siltstone, mudstone, and shale [17,18].
Energies 2019, 12, x 3 of 16 [15,16]. This area has complex geological regime, and having major tectonic features such as faults, shear zones, synclines, and anticlines can be seen in Figure 2. Twin headrace tunnels were excavated with two TBMs (tunnel boring machines), under a rock cover ranging from 1000 to 1900 m. These TBM tunnels passed through the Murree Formation, which comprises alternate beds of sandstone, siltstone, mudstone, and shale [17,18]. Among these rock units, the sandstone unit is massive and competent under a highly stressed environment. This is the strongest unit of the Murree Formation, and the bedding thickness varies from a few meters to 50 m. This rock unit has a mean uniaxial compressive strength (UCS) of 90 MPa during laboratory testing [19]; subsequent testing after the extreme rockburst produced very high UCS values, in the range of 130-170 MPa [13]. This massive and strong rock stores more energy due to its high stiffness as compared to other rock units in the area. Due to this high-strength property, this rock unit can have high stress concentrations. Plumb [20] found that the stresses are 20-40% higher in hard rocks compared to weak rocks. Yang et al. [21] recently studied the concentration of stresses in strong and weak rocks and concluded that in-situ stresses are very high in the strong layer of the NJHEP TBM tunnels.
The over-coring in-situ stress measurements were conducted in the sandstone beds because of frequent strain burst events in this rock unit during the construction of the headrace tunnels. The maximum principal stress (SH) has a horizontal direction, while the minimum principal stress (Sv) is vertical. In addition, the hydro-fracturing program of the nearby planned Kohala hydropower project has also indicated the same stress state trend of the area; SH > Sh > Sv [22]. These high horizontal stresses are due to the active tectonics of the area [23]. The major principal stress has a sub-horizontal orientation and is nearly perpendicular to the tunnel direction. The ratio (K0) between horizontal and vertical stress is up to 2.9, as shown in Figure 3. The yellow dots show high overburden stresses within the dotted lines limits, and the red dots show high abnormal stress influenced by high horizontal tectonic stress. Therefore, under great depth, stresses became concentrated in hard sedimentary rock under high stress concentrations due to tectonics, which has made conditions favorable for rockburst occurrence.   Among these rock units, the sandstone unit is massive and competent under a highly stressed environment. This is the strongest unit of the Murree Formation, and the bedding thickness varies from a few meters to 50 m. This rock unit has a mean uniaxial compressive strength (UCS) of 90 MPa during laboratory testing [19]; subsequent testing after the extreme rockburst produced very high UCS values, in the range of 130-170 MPa [13]. This massive and strong rock stores more energy due to its high stiffness as compared to other rock units in the area. Due to this high-strength property, this rock unit can have high stress concentrations. Plumb [20] found that the stresses are 20-40% higher in hard rocks compared to weak rocks. Yang et al. [21] recently studied the concentration of stresses in strong and weak rocks and concluded that in-situ stresses are very high in the strong layer of the NJHEP TBM tunnels.
The over-coring in-situ stress measurements were conducted in the sandstone beds because of frequent strain burst events in this rock unit during the construction of the headrace tunnels. The maximum principal stress (S H ) has a horizontal direction, while the minimum principal stress (S v ) is vertical. In addition, the hydro-fracturing program of the nearby planned Kohala hydropower project has also indicated the same stress state trend of the area; S H > S h > S v [22]. These high horizontal stresses are due to the active tectonics of the area [23]. The major principal stress has a sub-horizontal orientation and is nearly perpendicular to the tunnel direction. The ratio (K 0 ) between horizontal and vertical stress is up to 2.9, as shown in Figure 3. The yellow dots show high overburden stresses within the dotted lines limits, and the red dots show high abnormal stress influenced by high horizontal tectonic stress. Therefore, under great depth, stresses became concentrated in hard sedimentary rock under high stress concentrations due to tectonics, which has made conditions favorable for rockburst occurrence.

The Presence of Shear Zones in the Project Area
The NJHEP is located in a tectonically active range of the Himalayas, which is pierced by a series of regional thrust faults: the main front thrust, the main boundary thrust, the main central thrust, and the main mantle thrust [24][25][26]. These structures have been originated by the collision between two plate boundaries from the Mesozoic to late Cenozoic period. The headrace tunnels pass through the Murree Formation, which has a highly deformed geology due to tectonic stresses. As a consequence, different types of geological structures such as folds, faults, and shear zones are present in the area [27,28]. There are three local faults in the area: the Tanda fault, the Murree fault, and the Muzaffarabad fault. Details about these local faults are beyond the scope of this study, but the initial geological mapping program [29] of the NJHEP project reported 15 local, small-scale, sharp fault planes (denoted with "F" in series) and 26 shear zones present along the headrace tunnels. Among these faults and shear zones, few structures are expected to pass near a rockburst-susceptible environment, such as F2 and F3 as shown in Figure 2. F2 has an expected invert level chainage of 10 + 360 (m) with an orientation of 62/058 (dip/dip direction), while F3 has an expected invert level chainage of 9 + 350 (m) with an orientation of 25/080. The NJHEP headrace tunnels also pass near many shear zones due to these complex geological conditions. Some of these are present at the contact

The Presence of Shear Zones in the Project Area
The NJHEP is located in a tectonically active range of the Himalayas, which is pierced by a series of regional thrust faults: the main front thrust, the main boundary thrust, the main central thrust, and the main mantle thrust [24][25][26]. These structures have been originated by the collision between two plate boundaries from the Mesozoic to late Cenozoic period. The headrace tunnels pass through the Murree Formation, which has a highly deformed geology due to tectonic stresses. As a consequence, different types of geological structures such as folds, faults, and shear zones are present in the area [27,28]. There are three local faults in the area: the Tanda fault, the Murree fault, and the Muzaffarabad fault. Details about these local faults are beyond the scope of this study, but the initial geological mapping program [29] of the NJHEP project reported 15 local, small-scale, sharp fault planes (denoted with "F" in series) and 26 shear zones present along the headrace tunnels. Among these faults and shear zones, few structures are expected to pass near a rockburst-susceptible environment, such as F2 and F3 as shown in Figure 2. F2 has an expected invert level chainage of 10 + 360 (m) with an orientation of 62/058 (dip/dip direction), while F3 has an expected invert level chainage of 9 + 350 (m) with an orientation of 25/080. The NJHEP headrace tunnels also pass near many shear zones due to these Energies 2019, 12, 2124 5 of 17 complex geological conditions. Some of these are present at the contact points of different rock units, while others are within specific rock units. A shear zone (S5) with an orientation of 30/048 has an expected invert level chainage of 9 + 430 (m), which is near the intense rockburst zone. The thickness of this shear zone is 0.5-1 m.
In deep, hard rock tunnels, these structures are critical and are closely associated with rockburst occurrence. They can be at various locations, such as at the working face, on the crown, at the invert, at shoulders, and on the side walls of a tunnel, resulting in a concentration of abnormal stresses near these areas. Bawden [30] found one of these shear zones in a massive sandstone bed in the wall of headrace tunnel exposed after the extreme event of rockburst on 31 May 2015. Such rockburst influenced by the small-scale fault slip also occurred during the construction of the Tianshengqiao II hydropower station on the Nanpan River in China [31].

Extreme Event of Rockburst in NJHEP and Its Effects
According to Kaiser and Cai [32], during the TBM excavation, constant stress and the mechanical properties of the rock are responsible for the strain burst potential, which can result in constant bursting or no bursting. They believe that geological factors such as geological structures or stiff dykes control the local system stiffness and are the main cause of severity during rockburst, relating to their frequency and damage. According to this concept, the NJHEP headrace tunnels were built in vulnerable conditions, with strong rocks under an abnormal stress regime, which resulted in frequent small strain burst events at different locations in hard sedimentary sandstone, but no major rockburst event. The frequent strain bursts were due to the brittle failure of strong rocks. The extreme rockburst event in the NJHEP was influenced by one of the many shear zones extant on the path of the tunnels, and Bawden [30] has reported a shear zone exposed in a wall of a tunnel after the intense rockburst event, as shown in Figure 4a. Just prior to this most intense rockburst of the project, the excavating TBM in Tunnel 696 approached the shear zone, causing an unclamping of the existing shear plane resulting in an intense rockburst when the cutter head of the TBM was at chainage 09 + 706 at 11:35 p.m. on 31 May 2015 [33]. The shear zone near the tunnel boundary lessened the local system stiffness, resulting in this intense rockburst event. This structural plane was perpendicular to the in-situ principal stress, and, according to Feng [9], a large rockburst with a deep damage pit can be triggered and similar occurred in the NJHEP tunnel, resulting in failure of rock mass in the upper left shoulder of the tunnel, as shown in Figure 4b. This event was equal to a 2.4 magnitude earthquake and severely damaged the TBM 696 and the tunnel support system. Three workers died, and several others were injured after this event, and a large area of the tunnel was collapsed.
After this deadly event, the area was studied thoroughly to understand the rockburst mechanism. There were some abnormal conditions prevailing in the area. Different geological structures were present in the subsurface as shown in Figure 2. During excavation, the geological conditions before the major event of rockburst, were normal; for example, the bedding planes were perpendicular to the axis of the tunnel. When Tunnel 696 crossed the chainage 09 + 790, suddenly there was no bedding along the tunnel axis in the next 40 m. This was the area where abnormal geological conditions started during excavation, having changed orientation of the bedding plane with abnormal stress conditions. Earlier the bedding plan was perpendicular to the tunnel axis, and this has changed and is now parallel to the tunnel axis [30,34]. Additionally, tunneling has also disturbed the equilibrium conditions of the rock mass after excavation. Therefore, these complex geological settings led to more stress concentrations, which resulted in the development of a shear zone in the rock mass that finally triggered the most severe rockburst.

Numerical Simulation of the Extreme Rockburst
Rockburst is a dynamic phenomenon that involves the unstable failure of rock. The risk of severe rockburst is high when a geological structural plane is present near any tunnel along with change in equilibrium status of the area due to tunnel excavations. Different numerical studies have been done to explain the effect of these structures on rockburst occurrence. Zhang et al. [35] numerically evaluated the failure of a rock mass in the Jinping-II hydropower station. The blocking effect of the fault caused intense stress concentration, which resulted in increased shear strain energy near the fault which, in turn, caused severe seismic activity and energy release. Zhang et al. [36] have used the failure approach index (FAI) during numerical simulation and determined that structural planes led to local stress concentration that caused rock mass failure in the hanging wall. In the current study, the FLAC 2D numerical simulation has been used to evaluated the influence of the shear zone on tunnel stability, and its possible effect on the rockburst failure mechanism around Tunnel 696 and on the support system in the adjacent Tunnel 697 subjected to static and dynamic loading.

Modeling Scheme
The numerical model was constructed according to on-site conditions. The twin headrace tunnels, each 8.53 m in diameter, were built to actual size and were named 696 and 697 (based on TBM); looking upstream, they are also referred to as the right and left tunnels, respectively, as shown in Figure 5. The distance between these two tunnels was set to 33 m, which is same as it was in the field before the major rockburst. According to the rockburst investigation report [30], a shear zone was present in the Tunnel 696 wall after the extreme rockburst event, as shown in Figure 4a. It can be noticed that this and other, similar shear zones were overlooked during the TBM tunnel excavation. We tried to simulate such shear zones during numerical modeling. An interface element (a collection of triangular elements) representing the shear zone-type structures was used to simulate a shear zone in FLAC 2D. A nearly vertical interface element was placed near Headrace Tunnel 696 to reflect the actual field conditions, as shown in Figure 5, while there was no interface element near Tunnel 697. The purple color lines show the fixed boundary condition of the static model.

Numerical Simulation of the Extreme Rockburst
Rockburst is a dynamic phenomenon that involves the unstable failure of rock. The risk of severe rockburst is high when a geological structural plane is present near any tunnel along with change in equilibrium status of the area due to tunnel excavations. Different numerical studies have been done to explain the effect of these structures on rockburst occurrence. Zhang et al. [35] numerically evaluated the failure of a rock mass in the Jinping-II hydropower station. The blocking effect of the fault caused intense stress concentration, which resulted in increased shear strain energy near the fault which, in turn, caused severe seismic activity and energy release. Zhang et al. [36] have used the failure approach index (FAI) during numerical simulation and determined that structural planes led to local stress concentration that caused rock mass failure in the hanging wall. In the current study, the FLAC 2D numerical simulation has been used to evaluated the influence of the shear zone on tunnel stability, and its possible effect on the rockburst failure mechanism around Tunnel 696 and on the support system in the adjacent Tunnel 697 subjected to static and dynamic loading.

Modeling Scheme
The numerical model was constructed according to on-site conditions. The twin headrace tunnels, each 8.53 m in diameter, were built to actual size and were named 696 and 697 (based on TBM); looking upstream, they are also referred to as the right and left tunnels, respectively, as shown in Figure 5. The distance between these two tunnels was set to 33 m, which is same as it was in the field before the major rockburst. According to the rockburst investigation report [30], a shear zone was present in the Tunnel 696 wall after the extreme rockburst event, as shown in Figure 4a. It can be noticed that this and other, similar shear zones were overlooked during the TBM tunnel excavation. We tried to simulate such shear zones during numerical modeling. An interface element (a collection of triangular elements) representing the shear zone-type structures was used to simulate a shear zone in FLAC 2D. A nearly vertical interface element was placed near Headrace Tunnel 696 to reflect the actual field conditions, as shown in Figure 5, while there was no interface element near Tunnel 697. The purple color lines show the fixed boundary condition of the static model. The interface element is characterized by four different parameters, including shear stiffness (ks), normal stiffness (kn), cohesion (c), and friction angle (ϕ). The FLAC manual [37] clearly explains how to calculate the stiffness parameters. The following relationships were used to calculate the normal stiffness (kn) and shear stiffness (ks) [38].
The normal stiffness can be calculated as where kn = joint normal stiffness; E = rock mass Young's modulus; Er = intact rock Young's modulus; s = joint spacing. The shear stiffness can be calculated as where ks = joint shear stiffness; G = rock mass shear modulus; Gr = intact rock shear modulus; s = joint spacing. Different rock masses have different joint spacing [39]. The massive rock mass has very high chances of rockburst and has very few joints. On the other hand, according to [40,41], fractured rock mass has much less chance of rockburst occurrence because a high amount of stored strain energy is dissipated in fracturing the massive rock mass into the fractured one, having less strain energy available for rockburst occurrence. The stiffness values of the interface element were calculated by the aforementioned formulas and Palmstr m's recommended joint spacing. In this analysis, the values for the shear and normal stiffness of the interface element were calculated as ks = 3.77 × 10 2 MPa/m and kn = 9.43 × 10 2 MPa/m, respectively.
For this analysis, an elasto-plastic, nonlinear constitutive material model was used based on hard sandstone rock unit parameters. After simulating the geometric model, all sides of the model were fixed. For studying the stress characteristics of rock materials during loading and unloading conditions, the average stress values of in-situ stresses (σxx, σyy, σzz) were used. The maximum principal stress (σxx) was established as 60 MPa, the intermediate principal stress (σyy) was 37 MPa, and the minimum principal stress (σzz) was 35 MPa [21]. The mechanical properties of the rock mass applied during the numerical simulation are shown in Table 1, and the model was run for initial conditions. The stresses were distributed normally according to input parameters. The interface element is characterized by four different parameters, including shear stiffness (k s ), normal stiffness (k n ), cohesion (c), and friction angle (φ). The FLAC manual [37] clearly explains how to calculate the stiffness parameters. The following relationships were used to calculate the normal stiffness (k n ) and shear stiffness (k s ) [38].
The normal stiffness can be calculated as where k n = joint normal stiffness; E = rock mass Young's modulus; E r = intact rock Young's modulus; s = joint spacing. The shear stiffness can be calculated as where k s = joint shear stiffness; G = rock mass shear modulus; G r = intact rock shear modulus; s = joint spacing. Different rock masses have different joint spacing [39]. The massive rock mass has very high chances of rockburst and has very few joints. On the other hand, according to [40,41], fractured rock mass has much less chance of rockburst occurrence because a high amount of stored strain energy is dissipated in fracturing the massive rock mass into the fractured one, having less strain energy available for rockburst occurrence. The stiffness values of the interface element were calculated by the aforementioned formulas and Palmström's recommended joint spacing. In this analysis, the values for the shear and normal stiffness of the interface element were calculated as ks = 3.77 × 10 2 MPa/m and kn = 9.43 × 10 2 MPa/m, respectively.
For this analysis, an elasto-plastic, nonlinear constitutive material model was used based on hard sandstone rock unit parameters. After simulating the geometric model, all sides of the model were fixed. For studying the stress characteristics of rock materials during loading and unloading conditions, the average stress values of in-situ stresses (σ xx , σ yy , σ zz ) were used. The maximum principal stress (σ xx ) was established as 60 MPa, the intermediate principal stress (σ yy ) was 37 MPa, and the minimum principal stress (σ zz ) was 35 MPa [21]. The mechanical properties of the rock mass applied during the numerical simulation are shown in Table 1, and the model was run for initial conditions. The stresses were distributed normally according to input parameters.

Static Analysis Results
The twin tunnels were excavated via FLAC 2D when the model was in an equilibrium after initial conditions. Tunnel 697 was excavated first, and Tunnel 696 was excavated later as it was a trailing tunnel in the field. The numerical simulation results included contours of maximum and minimum principal stresses, a failed zone, displacement vectors near the failed zone, and shear displacement in the shear zone. The contours of the principal stresses (maximum and minimum) are shown in Figure 6a,b. These results showed that the maximum principal stress was in the horizontal direction, so the stress remained concentrated in the crown and the invert in Tunnel 697. Due to this stress concentration, both the crown and the invert were yielding, and the crown was heavily guarded with a rock bolt support system to ensure its stability and to decrease the influence of the potential failure zone.
According to the TBM tunnels rock support classification, 14 rock bolts of length 3.85 m were installed where the TBM excavation encountered a massive sandstone, which had a very high probability of rockburst. During the Tunnel 697 excavation, the construction was smooth and without any unstable failure. On the other hand, during the Tunnel 696 excavation, the maximum principal stress was concentrated near the extreme points of the shear zone because it acted as a barrier to equal stress distribution around the tunnel boundary. As the tunnel reached near the shear zone, a shear slip rockburst occurred due to the unclamping effect of Tunnel 696, resulting in shear displacement of the interface element. The yielded profiles of the tunnels are shown in Figure 6c; there was more yielded area near the extreme points of the shear zone due to the high stress concentration there. This high stress concentration near these points was released in the form of seismic energy, which damaged the tunnel badly due to its dynamic effect, which is explained subsequently.

Dynamic Analysis Results
Many researchers have discussed that mining-induced stress concentration near the shear zone contributes to seismic activity, which is a dynamic phenomenon. To capture the actual failure zone as it is in the field, additional dynamic loading was carried out, which further increased the stresses around the excavation. For dynamic simulation, quiet boundaries were assigned on the vertical, top, and bottom boundaries of the model to avoid the reflection of the outgoing wave back into the model. During this simulation, a seismic wave reached the boundary of the tunnel after travelling through the rock mass from the source point of high stress concentration near the extreme point of the shear zone as shown in Figure 6a, causing a dynamic stress rise on the already fractured zone after static stress analysis. The dynamic wave generated by a seismic event in rocks usually has a low frequency. In this case, a sinusoidal synthetic compressive wave was used as an input to pass through the model at a point where the stresses were concentrated near the upper boundary of the shear zone. The dynamic input was applied in the horizontal direction because, in the field, the maximum principal stress was horizontal. This resulted in the sudden release of stored strain energy in the rock mass, causing unstable rock failure with the violent expulsion of the failed rock in the tunnel.  During construction of the headrace tunnel, the rockburst event of May 31, 2015, badly damaged the upper left portion of Tunnel 696; therefore, hereafter only the upper section of this tunnel is discussed. The yielded area around Tunnel 696 after dynamic simulation is shown in Figure 7, which illustrates the shear and tensile failures around the boundary of the tunnel. The tunnel had a very large yielded area in the upper left quadrant. The dynamic simulation showed a 5 m depth of the failure zone, which is similar to that reported in the field, and this also confirms that a similar shear zone as reported by [30] was present near the wall of Tunnel 696, which caused this extreme rockburst event. This failure resulted in a large release of energy and caused enormous damage to installed supports 28 to 50 m behind the shield of the TBM. In contrast, the adjacent Tunnel 697 experienced no such major damage-only some minor damage to the installed supports. illustrates the shear and tensile failures around the boundary of the tunnel. The tunnel had a very large yielded area in the upper left quadrant. The dynamic simulation showed a 5 m depth of the failure zone, which is similar to that reported in the field, and this also confirms that a similar shear zone as reported by [30] was present near the wall of Tunnel 696, which caused this extreme rockburst event. This failure resulted in a large release of energy and caused enormous damage to installed supports 28 to 50 m behind the shield of the TBM. In contrast, the adjacent Tunnel 697 experienced no such major damage-only some minor damage to the installed supports. In underground excavation, the displacement of the surrounding rock mass is important for its stability. In that regard, it is essential to study the vector field of displacement around a tunnel. The presence of a geological structure near a tunnel increases the displacement field in a rock mass. After this dynamic simulation, the NJHEP Headrace Tunnel 696 had shown an increasing trend of displacement on the side where the shear zone was present, in contrast to the adjacent tunnel, which had no shear zone nearby. The displacement vectors corresponding to Figure 8a   In underground excavation, the displacement of the surrounding rock mass is important for its stability. In that regard, it is essential to study the vector field of displacement around a tunnel. The presence of a geological structure near a tunnel increases the displacement field in a rock mass. After this dynamic simulation, the NJHEP Headrace Tunnel 696 had shown an increasing trend of displacement on the side where the shear zone was present, in contrast to the adjacent tunnel, which had no shear zone nearby. The displacement vectors corresponding to Figure 8a showed larger deformations on the side of the shear zone. The high displacement vectors in the upper left quadrant of the tunnel confirmed the failure area where the extreme rockburst event of May 31, 2015, occurred, as shown in Figure 4a. The displacement pattern after a dynamic impact of rockburst mimics the actual rockburst damage. The shear displacement in the interface element (shear zone) can also be seen in Figure 8b. displacement on the side where the shear zone was present, in contrast to the adjacent tunnel, which had no shear zone nearby. The displacement vectors corresponding to Figure 8a showed larger deformations on the side of the shear zone. The high displacement vectors in the upper left quadrant of the tunnel confirmed the failure area where the extreme rockburst event of May 31, 2015, occurred, as shown in Figure 4a. The displacement pattern after a dynamic impact of rockburst mimics the actual rockburst damage. The shear displacement in the interface element (shear zone) can also be seen in Figure 8b.

The Effect of the Rockburst on Installed Supports
The designed support system in the NJHEP tunnels consists of expansion shell rock bolts, TH beams, shotcrete, and wire mesh containment. The seismicity due to the rockburst caused the most severe loading conditions for this system. After the major rockburst event of May 31, 2015, the damage was visible almost 70 m behind the TBM in Tunnel 696, and there was severe damage to the TBM, the tunnel profile, and permanent rock supports 28 to 50 m behind the TBM shield [33]. The most severe over break due to the rockburst occurred in the upper left quadrant of the tunnel. The FLAC 2D numerical analysis validated this observation. It is evident from the damage that the installed supports did not bear well such adverse rockburst conditions in a deep headrace tunnel. All the support systems were severely damaged, and most of the support elements were completely sheared off, as shown in Figure 9. The dynamic impact of rockburst has caused a complete failure of the support system, resulting in the total collapse of the excavation due to seismic energy. In this area of the tunnel, a large volume of rock failed, and the rock bolts lost their anchorage, indicative of high ground velocities. The yielded zone increased rapidly due to the seismic action of the concentrated dynamic stress, which had been confirmed by dynamic numerical analysis. In such an environment, where the rock mass is usually massive, the greater interaction between the rock bolts and the rock mass resulted in tensile failure of the rock bolts, particularly under the dynamic loading conditions. Significant debonding along the length of the rock bolts occurred, but the rock bolts remained anchored within the stable rock mass, and the same is observed in the field as shown in Figure 9a. The shearing of the TH beam ring is shown in Figure 9b. dynamic stress, which had been confirmed by dynamic numerical analysis. In such an environment, where the rock mass is usually massive, the greater interaction between the rock bolts and the rock mass resulted in tensile failure of the rock bolts, particularly under the dynamic loading conditions. Significant debonding along the length of the rock bolts occurred, but the rock bolts remained anchored within the stable rock mass, and the same is observed in the field as shown in Figure 9a. The shearing of the TH beam ring is shown in Figure 9b.  Damage was inflicted by the dynamic impact of the rockburst resulting in deformed and destroyed support system. The depth of the failure zone (5 m) was greater than the length of the installed rock bolts (3.85 m), which means that this support design methodology was completely flawed for such extreme conditions. The design of the support systems was based on empirical design guidelines, which are often not applicable in such a dynamic environment. This event was so severe, after this extreme rockburst, the adjacent Tunnel 697 incurred minor shotcrete spall due to the dynamic impact of the event, which was initially thought to be wide enough that the excavation of one tunnel would not affect the adjacent tunnel (33 m away). Nevertheless, according to Gary Peach [33], the event of May 31, 2015, was unusual in that the adjacent tunnel was affected due to its high intensity and dynamic impact. To check the dynamic impact of rockburst on the installed support of nearby tunnel, a subsequent support analysis was run for rock bolt support. The performance of the support was analyzed in FLAC 2D for Tunnel 697 only, because it was not possible to check the supports in Tunnel 696 behind the face in 2D.
The numerical analysis results showed that the axial forces in the rock bolts increased from the static to the dynamic stage, shown in Figure 10. Initially, there was a significant increase in the axial forces in those rock bolts, which were present on the side of the trailing Tunnel 696, which was excavated and passed by the existing Tunnel 697 that was already excavated and supported. The rock bolts did not reach their yielding point, but the bolts experienced higher axial loads than the previous stage. The dynamic wave due to the intense rockburst struck the new Tunnel 696, and it also affected the existing support system in the adjacent tunnel too, which is dynamically evaluated in FLAC 2D.
During this support analysis, the properties of the rock mass and the rock bolts remained the same as those used during static loading conditions. During dynamic analysis, the seismic wave approached the twin tunnels, which resulted in a sudden increase in rock mass yield around Tunnel 696 and an increase in axial forces in the adjacent tunnel's installed rock bolt support.
The numerical analysis results showed that the axial forces in the rock bolts increased from the static to the dynamic stage, shown in Figure 10. Initially, there was a significant increase in the axial forces in those rock bolts, which were present on the side of the trailing Tunnel 696, which was excavated and passed by the existing Tunnel 697 that was already excavated and supported. The rock bolts did not reach their yielding point, but the bolts experienced higher axial loads than the previous stage. The dynamic wave due to the intense rockburst struck the new Tunnel 696, and it also affected the existing support system in the adjacent tunnel too, which is dynamically evaluated in FLAC 2D. During this support analysis, the properties of the rock mass and the rock bolts remained the same as those used during static loading conditions. During dynamic analysis, the seismic wave approached the twin tunnels, which resulted in a sudden increase in rock mass yield around Tunnel 696 and an increase in axial forces in the adjacent tunnel's installed rock bolt support.

The Static and Dynamic Impact of the Shear Zone on Rock Mass
Geological structures act as a barrier, and stresses are normally accumulated near their boundaries. Usually, intense rockburst occurs in deep hard rock tunnels influenced by these structures. Manouchehrian and Cai [2] have discussed the influence of geological weak planes on rockburst occurrence. These planes are present near the tunnel boundary and cause high energy release, indicating rock failure near the excavation boundary. Feng et al. [42] analyzed the locations, intensities, and shapes of rockburst, which are controlled by such structural planes in the Jinping-II hydropower tunnels; they found that different angles of structural planes, with the tunnel axis and tangential stress, define the intensity of rockburst. In our previous work [14], it was shown that stresses are concentrated near the shear zone, which triggers rockburst statically when excavations

The Static and Dynamic Impact of the Shear Zone on Rock Mass
Geological structures act as a barrier, and stresses are normally accumulated near their boundaries. Usually, intense rockburst occurs in deep hard rock tunnels influenced by these structures. Manouchehrian and Cai [2] have discussed the influence of geological weak planes on rockburst occurrence. These planes are present near the tunnel boundary and cause high energy release, indicating rock failure near the excavation boundary. Feng et al. [42] analyzed the locations, intensities, and shapes of rockburst, which are controlled by such structural planes in the Jinping-II hydropower tunnels; they found that different angles of structural planes, with the tunnel axis and tangential stress, define the intensity of rockburst. In our previous work [14], it was shown that stresses are concentrated near the shear zone, which triggers rockburst statically when excavations approach such a highly stressed environment. The stress concentration near such a zone also contributes to seismic activity [8,43]. The rockburst is a seismic phenomenon that involves dynamic disturbance due to shear slip along the shear zone that can severely damage the excavation boundary and installed support. This is solely discussed in this paper. Therefore, the presence of the shear zone is the main controlling factor for rockburst occurrence, which has been further studied with the help of dynamic numerical simulation to find its dynamic effect on excavation boundary and on installed support.
Due to active tectonics, the area of the NJHEP project has a very high stress regime, and the maximum principal stress is in the horizontal direction (i.e., perpendicular to the tunnel axis). This is the first unfavorable condition for tunneling in the area. Secondly, the presence of strong, massive rock under the complex geology in the NJHEP area also led to large strain energy storage, which was finally released in the form of the extreme rockburst event of 31 May 2015. The project's twin tunnels passed through a series of synclinal and anticlinal structures; the stresses usually concentrate in synclines in already abnormal stress zones, and this made the conditions more vulnerable to rockburst. Finally, the presence of shear zone structures, which was confirmed by a post-event investigative report [30], acted as a barrier resulting in a higher stress concentration near the tunnel's boundaries. Such conditions were modeled during numerical simulation, and the model has a shear zone nearby Tunnel 696 only. During initial conditions, the in-situ stresses remained in a virgin state. When this tunnel passed near this shear zone structure, the in-situ principal stress increased up to 120 MPa near the shear zone, as shown in Figure 6a; which was initially 60 MPa as measured through an over-coring in-situ measurement program. On the other side, there was low stress concentration (80 MPa) near the boundary of the adjacent Tunnel 697. This stress variation reflected the significance of the shear zone's presence near the tunnel boundary, which was indicated in [8] in the case of deep mines. On the other hand, the minimum principal stress increased from 35 MPa up to 50 MPa as shown in Figure 6b near the shear zone, while there was no such stress concentration near the adjacent tunnel.
This localization of principal stresses caused the slip along the shear zone, which had a dynamic effect due to the release of seismic energy, which resulted in a violent ejection of rock. Tunnel 696 has tensile and shear zones near its boundaries, which can be seen in Figure 7, while the adjacent tunnel has no such failures (not shown). During the initial static phase, the shear displacement in the shear zone was 100.5 mm, as shown in Figure 8b. During dynamic loading, the maximum principal stress increased more due to the seismic impact of rockburst as the dynamic wave approached the tunnel boundary. The yielding of rock mass around Tunnel 696 also increased due to the dynamic impact of rockburst, and its depth increased up to 5 m, as shown in Figure 7, in the upper left quadrant, which was the same as the actual field condition after the intense rockburst event (compare with the small picture shown in Figure 7). The displacement vectors also confirmed that the maximum displacement concentration occurred in the upper left quadrant of Tunnel 696, as shown in Figure 8a. We have compared our findings with the results of famous case study of rockburst in the Jinping-II drainage tunnel carried out by [44].

The Static and Dynamic Impact of the Shear Zone on Installed Support
During static analysis, the initial supports were modeled in Tunnel 697. The rock bolts installed in the crown exhibited maximum axial force due to the high horizontal principal stress of the area. However, when Tunnel 696 was excavated nearby, the axial forces in the installed rock bolts increased on the right wall of Tunnel 697, which was alongside the trailing Tunnel 696, and forces in the crown rock bolts also increased due to the relaxing effect of Tunnel 696. During dynamic simulation, some rock bolts installed in Tunnel 697 showed an increase in axial forces due to the dynamic influence of rockburst, as shown in Figure 10, while in Tunnel 696 the rock support behind the face was badly damaged and affected the rock bolts, shotcrete, and wire mesh, as shown in Figure 9a,b. Most of these supports were very stiff and did not absorb the dynamic energy and thus failed in tension. Therefore, in such conditions there is a dire need to install yielding supports to absorb the extra energy produced by the dynamic impact of rockburst.

Conclusions
Rockburst is the sudden failure of rock mass associated with stress concentration along a shear zone type geological structure. This paper has discussed the extreme rockburst event of May 31, 2015, at the NJHEP project in Pakistan, which caused enormous loss of life and property and has also delayed the excavation of Headrace Tunnel 696. This project has the deepest hydroelectric tunnels in the Himalayan region. The devastating rockburst event of 2015 was influenced by the shear zone near the boundary of the tunnel. Such a zone could be easily investigated if it is present in the drill and blast section of a tunnel, but during TBM tunnel excavation such planes can be easily overlooked. Therefore, the presence of such structures must be checked for a good understanding of the likelihood of rockburst occurrence and its influence on a tunnel's support system.
In this study, FLAC 2D explicit coding was used to study the influence of a shear zone on stress concentrations and damage occurrence around twin tunnels, under static loading initially and then dynamically. During the static analysis, the principal stresses were concentrated near the boundaries of the shear zone, and the failure zone was also influenced by this structural plane. During seismic loading, more damage occurred to the rock mass around the boundary of Tunnel 696 due to its dynamic impact. The most damage occurred in the upper left quadrant of that tunnel, and numerical modeling results captured the failure pattern of the extreme rockburst event well. The failure zone of 5 m depth captured by the dynamic numerical analysis was the same as in the field. The displacement vectors were also concentrated at the same location where the maximum deformation occurred around the boundary of the tunnel. This intense rockburst badly damaged the support system in this tunnel and affected the adjacent tunnel due to its dynamic impact.
Therefore, the presence of a shear zone near the headrace tunnel has a very damaging impact on its boundary and support system because when a high-energy seismic wave hits the deep tunnel, it magnifies, and the damage caused by its dynamic impact increases the failure zone by many times. The numerical analysis results in this paper successfully modeled the dynamic response of the rock mass, and the approach presented in this paper can be useful for the dynamic analysis of rockburst rock influenced by the shear zone. Therefore, it is recommended that energy-absorbing support systems should be installed in future projects at great depth in the Himalayas, Pakistan.