Rock mechanics and DFN models in the Swedish Nuclear Waste Disposal Program

Discrete Fracture Network (DFN) is a modeling framework for fractured rocks. The core element is the description of geological medium as a network of discrete fractures that can be either generated from statistical distributions or imported as deterministic surfaces. It is an alternative to continuum methods with both advantages of easily integrating the statistical properties of fracture networks, and of not assuming any homogenization scale. DFN has been extensively used to describe fracture network flow properties supported by the fact that connectivity, which is a constitutive element of the network organization, is a key element of fluid percolation. Application of the DFN modeling framework to geomechanics is also promising and, conversely, DFN models will benefit from rock mechanics integration. Integration between DFN and rock mechanics modeling is in expansion in many fields and broad contexts. This includes prediction of mechanical effective properties, increased understanding of the fracture scales and indicators that control these properties, distribution of block sizes and shapes for block fall risk analysis, potential wave attenuation effect and fracture shear displacements caused by and within the fracture network induced by an earthquake, or hydromechanical effects for flow and transport predictions. These applications are relevant only if DFN models involve the right complexity and provide a reliable description of the fracture networks. DFN models also benefit from rock mechanics concepts to improve their realism as it is done with genetic models that mimic the growth and arrest of fractures according to stress conditions prevailing at the time of their formation.


Introduction
The Swedish Nuclear Fuel and Waste Management Co (SKB) is in charge of handling and safe disposal of nuclear waste in Sweden. More than 30 years of research have led to the KBS-3 method for the final disposal of spent nuclear fuel. This method is based on three protective barriers which perform two safety functions: isolation and retention. The spent nuclear fuel is encapsulated in copper canisters with a cast iron insert. Embedded in bentonite clay for isolation and mechanical protection, the canisters are deposited in vertical deposition holes in the floor of horizontal tunnels at a depth of around 500 m in crystalline bedrock (Figure 1). The deposition tunnels are also backfilled with bentonite and sealed with concrete plugs. Rock caverns, transport and main tunnels are backfilled and a top sealing is put in place after the disposal is completed. Identifying, understanding, conceptualizing and bounding the main rock mechanics and thermohydromechanical aspects that can influence the construction of the repository and the performance of the barrier system in the short-, mid-and long term is of critical importance for the safe geological disposal of nuclear waste [1].
Numerous processes in hard, brittle, fractured rock masses are controlled by the fracture network. This is why the Discrete Fracture Network (DFN) approach is increasingly being incorporated into rock mechanics analysis ( Figure 2). Additionally, stress state influences fracture aperture and consequently flow in fractures and it is therefore being considered in coupled hydromechanical flow and transport analysis. At the same time rock mechanics is lately being used in DFN related research as a way to further improve our understanding and representation of fracture networks. The fracturing process is being considered through the use of concepts such as fracture nucleation, propagation and interaction based on fracture mechanics. This is a way of possibly improving the representation of fracture interaction (i.e. fracture crossing, termination, etc) and its consequences on the connectivity of the fracture network.  This paper focuses specifically on the use of DFN within rock mechanics and coupled hydromechanical analysis and the aspects of rock mechanics which are being incorporated and driving forward the DFN approach within the Swedish nuclear waste disposal program. A more extensive conceptual presentation of the relation between rock mechanics and DFN is presented in [2] and [3].

Discrete Fracture Network (DFN) approach
The DFN modelling approach describes a fractured rock mass as a population of individual and idealized tabular objects, including their geometrical and physical properties. It is particularly suited for sparsely fractured crystalline rocks and multiscale fractured systems, as are often the host rocks selected for deep storage of nuclear fuel.
The generic term fracture refers to the idealized, fracture-like objects, whatever their size and geological nature, that are the discrete entities of the whole fracture population. The geometrical description of a DFN model by statistical distributions of geometrical properties is the foundation of the approach and a basis for the subsequent definition of the fracture hydraulic and mechanical properties of the population of fractures as a whole [3].
For each fracture of a DFN, the geometrical description encompasses at minimum a location, size, shape and orientation, as given by the spatial density distribution of the whole population, including the size, shape and orientation statistical distributions. Additional fracture characteristics, e.g., detailed in-plane variations of aperture and surface morphology and quality, potentially contribute to define the fracture mechanical, hydraulic or hydromechanical properties.
A typical DFN representation combines deterministic, semi-deterministic and stochastic fractures. If some rather large-scale features can be deterministically defined (from integration of geophysical, surface-based boreholes, outcrop, monitoring, and underground characterisation), for most fractures it is not possible to determine more than a small fraction of their properties and spatial distributions. Hence, they are described in terms of probability distributions, with statistical models and through correlation structures.
The stochastic aspect of a DFN description is necessary to predict the nature of fracture networks, in a statistical sense, between observations, that are below the resolution of geophysical methods. Similarly, it is necessary to use a statistical approach to describe fractures that are known to exist but for which there are no direct observations.
The DFN density distribution model [3,4] refers to the most general definition of a DFN recipe. It defines the number of fractures per unit volume of any given range of fracture properties. It embraces declinations of the model as an ensemble of fracture orientation sets and size scale dependency [5] to define the relative density evolution between "large" and "small" fractures at different modelling scales.
Simpler quantitative metrics can be derived from the density distribution model. For any given range of fracture parameters, one may derive the total fracture surface area per unit volume of rock mass (commonly noted following the notation introduced by [6]), a quantity often relative to borehole logging. One may derive too the total DFN equivalent volume per unit volume of rock mass (each fracture contributes as its size to the cube, as an equivalent volume) per unit volume of rock mass, a quantity, known as the DFN percolation parameter, which defines the connectivity level of a DFN [7] or the rock mass effective elastic properties under some conditions [8].
A DFN model is the recipe for creating DFN realizations and numerical representations of a rock mass from which a wide range of processes and applications can be simulated (e.g., example in Figure  3a). Depending on the process under study, the numerical model may rely only on the DFN representation (boundary element methods) or a synthetic rock mass specimen when both discrete fractures and embedding rock matrix mechanical properties are represented in the numerical modelling. In any case, the DFN representation is the medium model representation that is the closest to the reality. If the contribution of DFN modeling to complex geomechanical modeling is now well-established [2, 3, 9], the reverse is also true, in the sense that using principles of rock mechanics and fracture network growth truly improves the DFN geological realism. DFN models solely based on statistical distribution of geometrical parameters are often criticized for their inability to reflect the fractured system spatial organization (spatial correlation between fractures [10,11].
Process-based DFN models, fully geomechanical [12,13] and genetic ones [14], significantly improve the geological realism. They both built on physical and mechanical arguments of fracture and fracture network growth, yet with very different degrees of complexity. Fracture mechanics modelling and truly geomechanical models are still limited in their capacity to produce DFN systems sufficiently large to address applications relevant to geological systems. Genetic models define the same fracture mechanic process but in a much simpler way so that sufficiently large and complex DFN can be generated. The simplified fracturing process involves mainly three stages of a fracture life cycle, the nucleation, growth, coalescence and arrest processes. At each stage, interactions are considered between a growing fracture and the surrounding network.
The reliability of genetic models relies on the underlying physical rules, the ability to understand the geological processes that led to the observed fracturing patterns and the capacity to generate realistic density distributions and fracture pattern [3].
Finally, using fracture and rock mechanics arguments potentially brings up physical arguments to understand observations and therefore guide the DFN modelling process from data sampling to conceptual DFN model definition and calibration to specific site conditions.

Influence of background fractures on the response of co-seismic secondary fractures displacements for different earthquake scenarios
Large faults rupture and consecutive earthquake are one identified source of risk for the future repository, as it may induce significant shear displacements on off-fault pre-existing secondary fractures and damage the canisters. DFN modeling and rock mechanics intervene in this issue where one needs to identify which fractures of the fracture system (distance to the fault, size, orientation) will potentially display shear displacements sufficiently high to threaten canister integrity (Figure 4). In preliminary developments of scenarios and numerical modeling, the acceptable shear displacement limit is 50 mm and the fractures likely to be affected have a diameter of 300 m.  Figure 4. Left) schematic view of a canister damage due to fracture shear displacement, from [15], right) visualization of the highest shear displacements (from 10 cm to 20 cm) in the vicinity of the primary fault, from [16].
In this context, dynamic modeling for simulating the consequence of an earthquake on some fractures is challenging. Numerically the challenge is to adequately generate and propagate seismic waves in a rock mass medium including an explicit representation of the fractured system with at least one major large fracture (the Deformation Zones or the faults) and target secondary fractures, but also potentially a wider part of the fractured system (the "background fractures"). An additional challenge is to have a resolution high enough to adequately simulate the shear displacement on secondary offfault fractures. [17] initially presented a method to calculate the shear displacements with the Distinct Element Code 3DEC, a well-established code for modelling the mechanical behavior of jointed rock masses. The numerical implementation and the representation of the rock mass have continuously improved since then. Numerically, [17] focused on the use of 3DEC and [16] focused on the use of the coupled environment FLAC3D/PFC.
For initial modeling assumptions where the only explicit fractures in the model are the primary fault and secondary fractures with the worst possible orientation (to the primary fault) to maximize the shear displacements, both approaches are consistent and predict comparable shear displacements on the secondary fractures.
Thorough analyses on the response of the fracture system to an earthquake, with a focus on which fractures are ultimately significantly affected, is still in its infancy. More complexity is added step by step in the model to get closer to realistic conditions. [16] have tested the impact of including more fractures explicitly in the model, for fracture sizes in a range covering roughly one order of magnitude and maximum value of the order of the so-called here secondary fractures. Their results confirm the capacity of the numerical set-up to include more discrete fractures to adequately model the DFN response to the earthquake wave propagation. The preliminary simulations performed predict a rather negligible effect of the background fractures on the shear displacements of the secondary fractures and maximum shear displacements in the background fractures of the order of 10 to 20 cm (Figure 4). [18] have also tested the influence of background fractures over a different range of fracture sizes for the background (between 20 and 160 m in equivalent diameter). They reach similar conclusions for the target fractures effect and lower shear displacement on the background fractures. the effect seems to be related to the fracture size. All these results are preliminary.

Tunnel rock fall
Risk of rock fall rely on the geometrical structure of the fracture network and on the fracture mechanical properties. When the geometry of blocks can be identified, deterministically or statistically, their stability can be tested under dynamic conditions. However the number of tests is limited by the burden associated to numerical testing, and a limited number of tests may be insufficient to catch the variability of the true block distribution.
A first step forward is to establish a better knowledge of the block shapes statistical distributions and their relationship with the DFN statistical properties. Although DFN and block assembly can be seen as a dual representation of the same rock mass, the correspondence between them is not straightforward. In common cases blocks are not fully formed and can be seen as a collection of "almost" formed blocks ( Figure 5) connected to each other by rock bridges of various sizes. In this context a method based on image analyses techniques can be used to identify blocks assuming that rock bridges below a certain threshold value can be broken.
In this context, the input of the DFN modeling provides the basis of the block analyses and is adapted to evaluate blocks distributions defined from given rock bridges and range of fracture sizes, together with true tunnel shape.

Rock mass effective mechanical properties
Building a quantitative approach to derive effective rock mass properties from DFN models, or directly from DFN model(s) parameters, is a ground-breaking research recently developed by different groups (see below). This evolution is only briefly recalled below. An alternative approach to bridge the gap between an increasing understanding of the complexity of naturally fractured systems on the one hand (the DFN model(s)), and of the rock mass deformability and strength on the other, is developed since a few years. It consists in building the quantitative relationship between the DFN relevant geometrical and mechanical properties and the rock mass effective elastic and strength properties. Following the fundamental work in Kachanov and Laures [20], Kachanov [21], Guéguen and Kachanov [22] for populations of frictionless cracks, Davy, Darcel [23] are currently further developing the approach for multiscale DFN models considering the Coulomb slip model to describe the mechanical behavior of each fracture. In short, the method, at present, takes as input an entire DFN embedded in a matrix rock (defined by elastic properties), including for each fracture the mechanical normal and shear stiffness and friction angle (and optionally also cohesion, tensile strength and dilation angle). The rock mass equivalent elastic properties (stiffness tensor) are thus analytically calculated at a chosen scale. [23] show how the latter tend to vary with scale (size of the rock mass), until a threshold is reached, beyond which the scale effect vanishes. The transition scale, between the scale dependent and independent regimes, results from a combination between fracture sizes, stiffnesses and the rock matrix modulus. This analytical DFN based method has several advantages:  It provides a direct prediction of rock mass properties, based on individual fracture mechanical properties and on DFN structure.  It is adapted to quantify scale, anisotropy and spatial variability effects.  It is a tool for upscaling of properties in the context of multiscale DFN systems. Its development is still ongoing [24]. Combination between the alternative analytical-based [23] approach and the SRM numerical modelling is a key component for future modelling of stress distribution at different scales. This alternative approach for rock mechanics is very close to the logic behind DFN flow modelling applied to hydrogeological modelling. It is quantitative, based on the complete DFN description, and adapted to upscale properties, including their scaling and anisotropy effects.

Rock mass permeability dependency on stress
A commonly observed dependency of permeability with depth (Figure 3b), its sensitivity to injection pressure and earthquake are all evidence that permeability depends on stress. In-situ stress conditions at depth are heavily dependent on lithostatic pressure, tectonic stresses but also stress fluctuations related to large deformation zones and to the overall spatial distribution of the fractured system. Equivalent rock mass permeability is hence expected to be both depth dependent and anisotropic with potential consequences on flow and transport that are critical for the safety assessment of a repository at depth. More specifically for sparsely fractured rocks and almost impervious intact rock matrix, flow can occur only through the fractures. Hence, for such conditions the effective stress permeability dependency is likely heavily related to the properties of the DFN and to the impact of the stress on the DFN structure and the transmissivity of the fractures. The key condition for further understanding the role of the stress on the equivalent permeability is to combine hydro-, mechanical and DFN modelling [25]. In practice, rock mechanics is primarily used to predict the stress loading on each fracture of a DFN as an input for defining individual fracture transmissivity from the present stress state and its evolution.

Concluding remarks
The use of DFN in rock mechanics, and in particular in nuclear waste disposal, to explicitly capture and study the roll of the rock structure is a rising field of activity. Through a series of examples, we have highlighted the benefit of resorting to DFN modelling to address many applications in the field of rock mechanics, and also emphasized the significant contribution of rock and fracture mechanics to the building of geologically sound DFN models. Rock mass stress permeability dependence is a relevant example of the necessity to combine DFN models in both rock mechanics and hydromechanics applications.