Both canadian and international experts agree that an underground repository is the single best solution for ultra-long-term storage of nuclear waste, both from a longevity and safety perspective. Numerous countries have carried out extensive research and preliminary design analyses associated with site assessment for high-level, in-rock fuel storage facilities known collectively as deep geological repositories (DGRs).
For the safe disposal of Canada’s high level nuclear waste, an “Adaptive Phase Management” plan has also been undertaken by Canada’s Nuclear Waste Management Organization (NWMO). A vertical shaft provides access to the DGR, which consists of an extensive horizontal layout encompassing centrally located access drifts, as well as perimeter drifts. Perpendicularly oriented tunnels will serve as the entry points to the nuclear fuel container storage rooms. These placements rooms will form a vast array of panels (Figure 1).
Most recently, a design known as the “Mark II” has been developed for preliminary design of high level repositories. As shown in Figure 1, used fuel canisters would be placed across the width of a rectangular placement room in a staggered fashion, surrounded by compacted bentonite (Guo 2016).
CATEGORIZING EXCAVATION – INDUCED DAMAGE
In order to accurately classify the degree and extent of damage induced around the underground nuclear waste shafts, tunnels and storage voids, it is vital to understand the terminology associated with the excavation-induced damage.
The basic premise of the excavation damage zone (EDZ) categorization lies in the fact that the density and connectivity of fractures decrease in a radial fashion outwards from the surface of the excavation (Perras and Diederichs 2016). It is important to note that immediately surrounding the excavation wall is the construction damage zone (CDZ), which is independent of the EDZ. Although sometimes impossible, the construction induced damage around an underground excavation can be mitigated by altering the excavation method to better suit the in-situ conditions (Siren et al 2015).
The generalized EDZ, which can develop independent of excavation methodology, is composed of concentric annuli (Siren et al 2015; Lanyon 2011; Tsang et al. 2005; Emsley et al. 1997), which are classed as the innermost “highly damaged zone” (HDZ), the inner and outer EDZi and EDZo, and the outermost excavation influence zone (EIZ). These regions are schematically represented in Figure 2.
The outermost region of the EDZ defines the mechanical limit of influence of the excavation. The excavation influence zone or EIZ consists of rock that has been elastically strained (Perras and Diederichs 2016). The outer “damaged” zone (EDZo) undergoes small scale crack damage with isolated unconnected fractures, subsequently increasing the permeability of the surrounding rockmass by up to tenfold. The inner “fractured” zone (EDZi) contains partially continuous fractures that increase the bulk permeability by tenfold to a hundredfold. It should be noted that measurable dilation (inelastic volumetric strain) is only observed within the EDZi (Perras and Diederichs 2016).
Immediately adjacent to the excavation at higher stress levels, there may exist a zone known as the highly damaged zone (HDZ). Due to stress-induced damage and compounding structural and geometric effects within the HDZ, the rockmass in this zone typically presents a high degree of interconnected fractures. Spalling is also commonly observed (Perras and Diederichs 2016). These fractures pose a serious issue as they result in irreversible damage and further promote fluid conductivity through the “fractured” zone (Massart and Selvadurai 2014; Bossart et al. 2002; Souley et al. 2001), with equivalent permeability increases of greater than four orders of magnitude.
In terms of the mechanics of stress-induced damage (in the various EDZs in Figure 2), it is helpful to consider the stress space plotted on the right side of Figure 2. If the excavation does not move the stress state into either tension or above the compressive “crack initiation threshold” defined as CI (Diederichs and Martin 2010), no damage can occur, as only an EIZ is present. If tensile rupture is generated however, or if shear failure takes place at higher confinements, an HDZ will form. Due to the effects of extension crack propagation at lower confining stresses, HDZ in the form of spalling can occur near the excavation boundary at stresses below the yield envelope for lab samples (critical damage limit or CD). This phenomenon is shown by the non-linear “S-shaped” curve bounding the red zone in Figure 2. The physics behind this S- shaped curve, damage initiation, propagation and the spalling limit are discussed extensively by Diederichs (2003, 2007) and Diederichs et al. (2004). Between damage initiation at CI and yield at CD, a zone of progressive distributed damage exists (EDZo to EDZi), with fracture intensity and connectivity increasing at higher differential stress levels.
NUMERICAL MODELLING, INPUT PARAMETERS AND MODEL SET-UP
Strength models and material Properties
To reproduce the non-linear “S-shaped” spalling model shown in Figure 2, Diederichs (2007) presents a simplified constitutive model (Damage Initiation and Spalling Limit, or DISL model) that can be used with analysis software to accommodate a strain weakening constitutive model, using either the generalized Hoek-Brown model (Hoek et al. 2002) or an equivalent Mohr-Coulomb approach. Equations and parameters for the DISL Hoek-Brown approach are shown in Figure 2. When used with a non-linear finite element or finite difference code, the DISL model is capable of predicting brittle rock spalling as a function of confinement, within underground excavations. Damage begins when the stresses exceed the “crack initiation threshold” (“CI” for unconfined conditions) represented by “peak” or initial strength parameters in DISL. The limit coinciding with “CD” in Figure 2 is the upper bound yield strength obtained in the lab. The “spalling confinement limit” (defined by post yield parameters in DISL) controls the strain weakening behaviour, representing spalling damage at low confinements. Equivalent, Mohr-Coulomb peak and residual envelopes were also determined by matching the non-linear DISL envelopes formed from the generalized Hoek-Brown parameters. A similar but alternative approach based on Mohr-Coulomb parameters is described by Hajiabdolmajid et al (2002).
Table 1 depicts the material properties of the sedimentary rock unit used during the modelling exercises. These values were based on the results of extensive testing carried out by the NWMO pertaining to the Cobourg Limestone sedimentary rock unit (Al et al. 2011). An intact Young’s Modulus of 42 GPa and a Poisson’s’ Ratio of 0.2 were used.
Numerical modelling software
As members of the RocScience suite, RS2 (Phase2) version 9.0 and RS3 version 1.0, finite element method (FEM) programs with continuum modelling were used for 2D and 3D applications respectively (see www.rocscience.com for more info).
Model geometry
Various tunnel shapes were assessed to gather an understanding of their impact on the EDZ propagation within a brittle rockmass at depth. The first of the three tunnel geometries consist of a rectangular tunnel 3.2m in width and 2.2m in height, replicating the Mark II nuclear waste storage rooms (Guo 2016). To explore the effects of rounded corners on the predicted decrease in EDZ propagation, the corners of the aforementioned rectangular profile were rounded by 0.5m on all edges. Finally, a horseshoe shaped tunnel was also assessed, with the tunnel dimensions kept similar to the Mark II rectangular storage tunnel design (Guo 2016), to draw comparative conclusions in terms of EDZ propagation.
Stresses and Depth
Each of the various models were assessed at a depth of 700m below ground surface. Major and minor principal stresses were calculated based on two in-plane stress ratios (KHh) of 1.5 and 2.0, to depict a suitable range within the Canadian context (Al et al. 2011).
NUMERICAL MODELLING RESULTS
2D shape effects and impact on EDZ
As anticipated, a larger KHh ratio resulted in a greater EDZ development, regardless of the tunnel profile. To assess the impact of the shape effects (corner geometries) on tunnel yield, the extent of yielded elements was measured and compared. In general, the rounding of the right angle corners act to reduce the overall yield by 1 – 20 per cent. The extent of damage in the floor of the horseshoe shaped tunnel numerically corresponds to that of the rectangular tunnel. Typical results are shown in Figure 3.
One striking aspect of Figure 3 is the well-defined block created above the two rectangular profiles and created in the floor of all the tunnels, compared to the continuous yielded zone above the roof in the arched case. This creation of a curved yield zone encompassing a minimally damaged block of rock is referred to in brittle rock engineering as the “baggage zone” (Kaiser et al.1996). It is the portion of un-yielded rock between the tunnel roof/floor and the concentrically yielded rockmass. This phenomenon is a result of the fact that rectangular tunnels rupture after CI is reached and move towards a more circular geometry. Per Kaiser et al. (2006), this “baggage zone”, while requiring support, also acts as an energy absorbing barrier to dynamically yielding rockmass, and can serve to increase tunnel safety if kept in place. The authors have observed such “baggage zones” in mining tunnels.
Hoek-Brown DISL versus
Mohr-Coulomb DISL
The impact of using generalized Hoek- Brown or equivalent Mohr-Coulomb failure criterion at different (increasing) stress regimes was also explored. For staged excavations simulated in RS2, aforementioned results prove that a lower stress regime is far more sensitive to the use of different failure criteria. It was discovered that for the case of a KHh of 1.5, the tunnels modelled with generalized Hoek-Brown failure criterion present approximately double to triple the amount of yield around the tunnel corners, as compared to tunnels assessed with equivalent Mohr-Coulomb parameters.
When modelling the various geometries with a KHh of 1.5 and Mohr-Coulomb failure criterion, EDZ development did not occur either above the roof or within the floor for both the rectangular tunnel and the tunnel with rounded corners; damage was limited to the corners. Conversely, these geometries experienced concentric yield about the roof and floor when modelled with generalized Hoek-Brown failure criterion (while still using a KHh of 1.5), generating a “baggage zone” (see Figure 4).
When using an increased KHh of 2.0, the tunnel geometries presented nearly identical results when modelled with either generalized Hoek-Brown or equivalent Mohr-Coulomb failure criteria. Unlike results obtained when using a KHh of 1.5, a “baggage zone” occurred when modelling with either of the two failure criteria at a KHh of 2.0. For both cases, the yield zone around the roof of the horseshoe tunnel was relatively insensitive to the yield criterion selected.
Impacts of 2D shape on EDZ
Categorization
Perras and Diederichs (2016) have developed a system for delineating the extent of the EDZ regions. The EDZ development around a brittle circular tunnel was categorized based on results obtained from numerical modelling of a limestone. It has been proven that principal stresses (?1 and ?3), yielded elements, volumetric strain and maximum shear strain aid in delineating the EDZ regions (Perras and Diederichs 2016). These parameters were all measured radially from the tunnel surface, aligned with ?3 (direction of most damage). The HDZ/EDZi boundary is defined by the onset of positive ?3, while the EDZi/EDZo boundary is attributed to the reversal in the volumetric strain (lowest value). Finally, the extent of plastic yield (yielded elements) serves as the EDZo extent. These boundaries have been highlighted in Figure 5.
Modelling results using the Hoek-Brown failure criterion with a KHh of 2.0 in Figure 5 suggest that the Diederichs and Perras (2016) guidelines are not entirely appropriate for non-curved geometries. Modelling results of the rectangular tunnel with rounded corners present similar results to the rectangular tunnel. The existence of the “baggage zone” raises the question as to its inclusion within the HDZ, when categorizing EDZ for nuclear waste disposal purposes. Secondly, there are measurement issues when delineating the HDZ/EDZi and EDZi/EDZo boundaries, since the onset of positive ?3 has been noted to occur after the volumetric strain reversal.
Furthermore, multiple reversals in positive-to-negative ?3 are noted to occur within the yielded rockmass due to this “baggage.” Modelling of horseshoe tunnel geometry supports the results presented by Diederichs and Perras (2016). Note the lack of “baggage zone” in Figure 5b, due to the arched roof geometry, which mimics that of a circular excavation.
Overall, the generalized guidelines for model interpretation need to be revised to encompass corner and “baggage” issues.
2D vs 3D modelling and EDZ development
In using both RS2 and RS3, it is possible to compare the results of EDZ development for staged excavations when modelled respectively in 2D and 3D views. The RS3 model is shown in Figure 6.
It is evident that a 2D simulation cannot capture the complex stress rotations and variations that occur as the tunnel advances. It has also been proven that the degree of stress rotation can dramatically affect fracture propagation and subsequent damage development (Eberhardt 2001 and Diederichs et al. 2004). To understand modelling discrepancies with respect to EDZ development in brittle rock, a comparison between the results obtained using RS2 (2D) and RS3 (3D) has been conducted. A fixed point in the tunnel roof will serve to provide the data for the following discussion.
In RS3, the true 3D excavation sequence was modelled with 1.0m blast rounds. The 2D model staging (reduction in internal pressure to simulate tunnel advance) was completed according to the longitudinal displacement profile (LDP) calibration approach of Vlachopoulos and Diederichs (2009, 2014). Two options for calibrating the LDP for a 2D analysis are compared to a 3D analysis in Figure 7. The use of a circular tunnel with a radius of 1.6m, along with isotropic stress conditions (LDP1) was used to provide a calibration by relating 2D model stages to representative distances from the face in the rectangular test model. In addition, a rectangular tunnel with anisotropic stress conditions (KHh of 2.0) was used to calibrate the 2D rectangular test model (LDP2).
For these analyses, the evolution of boundary stress versus tunnel advance is shown in Figure 7.
Note however, the difference in the stress increase experienced in the 3D model as the face passes the reference point (shown as “0” m from the face). This difference is a real phenomenon that is captured only in the 3D model. This spike in stress can result in damage that would not appear in the 2D model. This phenomenon is evident in the final damage (yield and shear strain) plots from the 2D and 3D analysis in Figure 8. Upon quantifying the EDZ development about each of the three tunnel geometries, it was determined that the RS2 models present 1.4 to 2.4 times the extent of damage (measured radially from the excavation surface) compared to the damage measured in RS3. Consequently, there are clear discrepancies in brittle damage development using 2D versus 3D models in high stress cases.
CONCLUSIONS
This study highlights several issues and challenges encountered when modelling brittle EDZ evolution with a continuum model. The authors have emphasized the differences between 2D and 3D tunnel simulation and have demonstrated some of the impacts of excavation geometry on brittle damage development.
Corners and flat boundaries create potentially undamaged “baggage zones” that have been observed in mining tunnels, where stresses increase after development due to mining. Three dimensional modelling of single advancing tunnels demonstrates however that in some cases, stress changes around the advancing face creates a continuous brittle damage zone, while the 2D simulation predicts localized fracture development bounding an undamaged “baggage zone” in the roof or walls. The true nature of EDZ in this case depends on the geometry and the blast round length (excavation sequencing) and is the subject of future investigation.
The choice of approximation for the DISL brittle damage model (generalized Hoek-Brown or equivalent Mohr- Coulomb parameters) impacts the nature and extent of EDZ for moderate stress states.
Finally, it is also important to be aware of the process of stress path recreation techniques and limitations of the LDP (Longitudinal Displacement Profile) when calibrating a non-circular 2D model.
