Hydrogeologic Characterization of Fractured Rock Masses Intended for Disposal of Radioactive Waste

There are currently 441 nuclear power reactors in operation or under construction distributed over 30 countries (International Atomic Energy Agency, 2011). The global radioactive waste inventory reported as storage in 2008 was approximately 17.6 million cubic meters: 21% short-lived, lowand intermediate-level waste, 77% long-lived, lowand intermediate-level waste and 2% high-level waste (International Atomic Energy Agency, 2011). There is a consensus among most of the scientific community that geologic repositories offer the best solution for the long-term disposal of radioactive waste. In the United States, for example, geologic disposal is considered the only technically feasible, long-term strategy for isolating radioactive waste from the biosphere without active management (Long & Ewing, 2004; National Research Council, 2001; Nuclear Energy Agency, 1999).

velocity zones. These processes are present for both unsaturated and saturated hydrologic conditions. Thus, hydrogeologic characterization of a potential repository must provide estimates of advective transport velocity based on permeability and porosity of the host rock; delineate the surrounding groundwater flow system including regional flow directions, hydraulic gradients, and zones of recharge and discharge; and determine the potential significance of radionuclide retention mechanisms.
In this Chapter, we focus our attention on the hydrogeologic characterization of fractured rock masses intended for the disposal of radioactive waste. The emphasis on fractured rock masses is two-fold. First, many repository sites proposed for high-level radioactive waste disposal are comprised of low-permeability hard rock (e.g., volcanic, igneous and sedimentary), although softer geologic materials such as salt domes and thick clay sequences are also under consideration. Low-permeability rock masses have little or negligible matrix porosity and permeability, with connected networks of discontinuous fractures imparting secondary porosity and permeability. Second, predicting radionuclide transport in fractured media is a formidable challenge as the spatially-discontinuous nature of fracture networks, along with high degrees of heterogeneity within fracture properties, lead to highly anisotropic flow systems with complex patterns of fluid flow and subsequent radionuclide migration (de Dreuzy et al., 2001;Klimczak et al., 2010;Neuman, 2005;Reeves et al., 2008a;b;c;Schwartz et al., 1983;Smith & Schwartz, 1984). Detailed radionuclide transport predictions are typically reliant on numerical simulations that incorporate site-specific fracture data (e.g., Arnold et al., 2003;Bodvarsson et al., 2003;Cvetkovic et al., 2004;Pohll et al., 1999;Pohlmann et al., 2002;2004;Reeves et al., 2010;Robinson et al., 2003;Smith et al., 2001), though some analytical techniques for first-cut approximations have been developed (Reeves et al., 2008b;Zhang et al., 2010).
This Chapter is not designed to provide a full treatise on the characterization of rock masses for development of complex numerical models used to predict radionuclide transport. Rather, we present an alternative approach where site-specific fracture network properties can be used to infer flow and transport characteristics of fracture networks by expanding the framework proposed by Reeves et al. (2008c). This framework can be used to qualitatively evaluate the suitability of candidate rock masses intended for the disposal of high-level radioactive waste based on fracture statistics. In stark contrast to current evaluation approaches that rely on costly field investigations to supply data to numerical models for radionuclide transport predictions, this type of evaluation promotes both time and economic savings by screening candidate fractured rock masses according to relatively simple criteria obtained from fracture characterization efforts.

Fractured rock characterization
Fractures are spatially discontinuous features that exhibit strong variability in geometrical and hydraulic properties. This variability is a result of the complex interplay between current and past stress fields, rock mechanical properties (i.e, Young's modulus, Poisson's ratio), mechanical fracture interaction and distributions of flaws or weakness in a rock mass. Fractured rock masses are typically characterized during field campaigns that measure fracture attributes from a number of sources including boreholes, rock outcrops, road cuts and tunnel complexes. Seismic techniques can also be used to image fault structures in the subsurface. Hydraulic properties of fractured media can be either inferred from fracture aperture or hydraulic tests performed on boreholes. Full characterization of fractured rock masses is not possible since known fracture locations and their attributes consist of an extremely small sample of the overall fracture network, i.e., any fracture characterization effort grossly undersamples a field site due to limited accessibility to the fractures themselves. Fracture data, however, can be used to generate representative, site-specific fracture networks through the derivation of probabilistic descriptions of fracture location, orientation, spacing, length, aperture, hydraulic conductivity/transmissivity and values of network density (Figure 1). With the exception of adsorption, these are the statistical properties that form the basis of the geologic repository screening framework. Statistical analysis of fracture attributes will be extensively covered in this section.

Orientation
Fracture networks will typically have two or more fracture sets characterized by fracture orientation (e.g., Barton, 1995;Bonnet et al., 2001;Bour & Davy, 1999;Pohlmann et al., 2004;Reeves et al., 2010). The presence of at least two intersecting sets of fractures reflects the physics of rock fracture propagation where two sets of fractures can arise from a single stress field (Jaeger et al., 2007;Twiss & Moores, 2007). Unless fractures are very long such as regional-scale faults, it is important from a flow perspective to have at least two intersecting fracture sets to promote connectivity across a rock mass.
The orientation of fracture planes is denoted by strike and dip convention. Analysis of fracture orientation begins with projecting the poles to fracture planes onto a stereonet and using contours of pole density to identify fracture sets ( Figure 2). Upon identification of fracture sets, mean orientation and the variability of fracture poles for each fracture set can be determined.
The distribution of fracture orientation is usually modeled using a Fisher distribution (Fisher, 1953): where the divergence, x (degrees), from a mean orientation vector is symmetrically distributed (− π 2 ≤ x ≤ π 2 ) according to a constant dispersion parameter, κ. The Fisher distribution is a special case of the Von Mises distribution, and is similar to a normal distribution for spherical data (Mardia & Jupp, 2000). The extent to which individual fractures cluster around the mean orientation is described by κ where higher values of κ describe higher degrees of clustering. It is our experience that values of κ are commonly in the range of 10 ≤ κ ≤ 50 for natural fracture networks. Stochastic simulation of Fisher random deviates in the discrete fracture networks following this section is based on the method proposed by Wood (1994).
The Bingham distribution provides an alternative to the Fisher distribution for cases in which fracture strike and dip are asymmetrically clustered around mean fracture orientations (Bingham, 1964): where κ 1 , κ 2 and κ 3 are dispersion coefficients that satisfy the condition: κ 1 ≤ κ 2 ≤ κ 3 = 0, M 1 , M 2 and M 3 are the column vectors of matrix M, E is the eigenvector matrix, and F(1/2; 3/2; E) is a hypergeometric function of the matrix argument. The probability distribution function described by (2)   Illustration showing the correspondence between two-and three-dimensional fracture networks. The three-dimensional network (top) is generated according to two fracture sets with significant variability about mean fracture orientations, a power-law length distribution exponent of a = 2.0 and a relatively sparse density. The two-dimensional network at the bottom left is computed by projecting all fractures onto the yellow horizontal slice located in the center of the three-dimensional DFN. The two-dimensional network on the bottom right is the result of identifying the hydraulic backbone by eliminating all dead-end fracture segments and non-connected clusters. Once the hydraulic backbone is identified, flow and particle transport can then be computed for the network.
mathematically possible to use (2) for the stochastic generation of asymmetric deviates. The authors are currently developing a method to simulate asymmetric deviations from mean fracture set orientations.  Fig. 2. Stereonet plots of poles to fracture planes with contour plots of all poles (left) and identified fracture sets along with prior probability (right). From Reeves et al. (2010).

Spacing
Fracture spacing refers to the linear distance between fractures. This distance also provides a length scale for unfractured matrix blocks. Fracture spacings from a data set require a correction (Terzaghi, 1965): to convert the apparent spacing D ′ measured along a transect to true fracture spacing D.
Values of α denote the angle of the transect relative to the mean fracture orientation or a pre-determined reference direction ( Figure 3). Apparent spacing is equal to true spacing if the transect is perpendicular to the mean fracture orientation or reference direction. If α = 90 • the Terzaghi correction factor f = sin(α) reduces to 1.  as the total number of fractures along the distance of a transect. Average fracture spacing [units of length, L] is simply the inverse of fracture frequency, and defines the size of the unfractured matrix block. This metric and its influence on radionuclide retention will be discussed in a later section.
Spacing in natural fracture networks is most commonly an exponentially distributed random variable. This can be tested by plotting the inverse empirical cumulative distribution function of fracture spacing, also known as a survival function (Ross, 1985). If spacing is exponential, the probability decay of the tail of the empirical fracture spacing distribution will exhibit a straight line on a semi-log plot. A Poisson point process defined by independent and identically distributed uniform deviates provides easy generation of exponential spacing (Ross, 1985).
Other possible distributions of fracture spacing include uniform (Rives et al., 1992) and fractal clustering (Barton, 1995;Darcel et al., 2003), both of which are considered extreme end members. Uniform spacing may occur in thin geologic layers which restrict fracture growth in the vertical direction and promote long horizontal fracture growth with nearly constant spacing (Rives et al., 1992). Exact causes of fractal clustering are less known and may be related to the role of mechanical fracture interaction during propagation that likely controls fracture length and spacing (Ackermann & Schlische, 1997;Darcel et al., 2003;Olson, 1993;Segall & Pollard, 1983). Networks with fractal clustering can be generated via a multiplicative cascade process (Mandelbrot, 1974;Schertzer & Lovejoy, 1987).

Length
Fracture length denotes the horizontal trace length of a fracture. There is a consensus in recent literature that fracture lengths above a lower length cutoff, l min , are power-law: with a power law exponent, a, that ranges between 1 and 3 in natural fracture networks (Bonnet et al., 2001;Bour & Davy, 1997;1999;Renshaw, 1999). C is a constant based on l min and a. Though lognormal distributions of fracture length have been reported in the literature, they are a result of improper sampling of the largest fractures within a sampling window. Lognormal distributions easily arise in data sets with power-law tails if the largest values are censored.
Determination of the distribution of fracture length is similar to that of fracture spacing and involves the analysis of an inverse empirical cumulative distribution function. Fracture lengths that are power-law will exhibit linear trends on a log-log plot for the tail of the distribution. In this example, the tail of the distribution refers to the greatest 5-10% of length values. The slope of the power-law trend of the data is equal to −a.
Truncations can frequently occur in fracture length data due to constraints imposed by the finite scale of the sampling window. For example in Reeves et al. (2010), the longest fracture measured in a tunnel drift was parallel to the drift and was approximately two-thirds of the total drift length. Instead of choosing between a traditional power-law or lognormal distribution, an upper truncated Pareto (power law) model (Aban et al., 2006): was used to compute the power-law trend in the data, where L (1) , L (2) ,...,L (n) are fracture lengths in descending order and L (1) and L (n) represent the largest and smallest fracture lengths, γ and ν are lower and upper fracture length cutoff values, and a describes the tail of the distribution. Truncated power law models like (5) can also be useful for imposing an upper length scale to the generation of stochastic networks at the regional scale. Lacking evidence of domain-spanning faults (with the exception of bounding faults of the stock itself) f o ra5k mwide granitic stock, Reeves et al. (2010) assigned an upper limit of 1 km in the stochastic generation of fault networks.

Hydraulic conductivity
Boreholes are commonly used to characterize fractured media. Borehole geophysics can provide useful information about fractures within the rock including fracture frequency, orientation, aperture and mineral infilling. Fracture aperture, defined as the width of the void space normal to fracture walls, can be used to infer hydraulic properties of fractures. Fracture apertures at land surface have low confining stresses that are not representative of subsurface confining stresses within a rock mass. We therefore recommend that fracture aperture values used to compute flow are measured in boreholes where in-situ stress is preserved.
The cubic law, a solution to the Navier-Stokes equation for laminar, incompressible flow between two parallel plates, describes a general relationship between fluid flow and fracture aperture (Snow, 1965): where fluid discharge per unit width, Q [L 2 /t], is proportional to the cube of the hydraulic aperture, b. Similar to Darcy's Law, the cubic law (6) assigns discharge through a fracture as a linear function of the hydraulic gradient, ∇h. The relationship between hydraulic aperture and transmissivity (T), hydraulic conductivity (K) and permeability (k) is described by: 12µ b 2 and k = b 2 12 , respectively. Fluid-specific properties density, ρ, and viscosity, µ, allow for conversions between permeability and hydraulic conductivity/transmissivity. As a note of caution, the relationship between mechanical aperture, the physical distance between fracture walls, and hydraulic aperture, the equivalent aperture for a given flow rate, is unclear. As a general rule, hydraulic aperture is typically smaller than mechanical aperture (Chen et al., 2000;Cook et al., 1990;Renshaw, 1995;Zimmermann & Bodvarsson, 1996). Discrepancies between mechanical aperture and hydraulic aperture are attributed to surface roughness, flow path tortuosity, and stress normal to the fracture. Though empirical correction factors have been used to correlate mechanical and fracture apertures (Bandis et al., 1985;Cook et al., 1990;Renshaw, 1995), no method is reliable for a wide range of aperture values.
Hydraulic testing of boreholes yields reliable estimates of fracture T and K. While there are many different hydraulic testing techniques, the isolation of specific intervals during testing with the use of dual-packer systems provides the best data to characterize the distribution of transmissivity/hydraulic conductivity. These tests yield flow rate information for applied fluid pressures, which also allows for the inverse computation of hydraulic aperture using (6). These aperture values, in addition to T and K estimates, are useful for parameterizing flow and transport models.
Studies in highly characterized rock masses have shown that fracture K is extremely heterogeneous and may encompass 5 to 8 orders of magnitude (Andersson et al., 2002a;b; Guimerá & Carrera, 2000;Paillet, 1998). Often the distribution of K (and T) is thought to be lognormal: where x is the mean and σ is the standard deviation. Values of log(σ K ) are typically around 1 for fractured media (Andersson et al., 2002a;b;Stigsson et al., 2001). However, other studies suggest power law distributions (Gustafson & Fransson, 2005;Kozubowski et al., 2008) and that these lognormal distributions, similar to length, could be caused by censoring flow data possibly due to instrument limitations. Additionally, flow through rough-walled fractures can be non-Darcian (Cardenas et al., 2007;Qian et al., 2011;Quinn et al., 2011). This may further complicate the estimation of hydraulic conductivity in field hydraulic tests as flow is no longer linearly proportional to a pressure gradient as described by (6).

Density
Fracture networks consist of two-dimensional planes embedded within a rock matrix ( Figure  1). The lack of access to the total rock volume makes it impossible to directly measure the three-dimensional fracture density of a rock mass. Instead, three-dimensional density for discrete fracture networks is estimated from density measurements of lower dimensions, i.e., one-dimensional fracture frequency from boreholes and/or tunnel drifts or two-dimensional fracture density from outcrops and fracture trace maps. Definitions of fracture density according to dimension are: one-dimensional density (also known as fracture frequency), , is expressed as the ratio of total number of fractures, f i , to transect length, L: , is expressed as the ratio of the sum of fracture lengths, l, to area, A: ρ 2D = A −1 ∑ n i=1 l i ; and three-dimensional fracture density, ρ 3D [L 2 /L 3 ], is expressed as the ratio of the sum of fracture plane area, A i , to rock volume, V: ρ 3D = V −1 ∑ n i=1 A i . Numerical techniques can be used to upscale one-dimensional fracture frequency [L −1 ] estimates to a three-dimensional spatial density [L 2 /L 3 ] (Holmén & Outters, 2002;Munier, 2004). For example, one-dimensional transects can be used to upscale two-dimensional networks by adding fractures until the one-dimensional transect density is satisfied along several transects placed along the two-dimensional network. Three-dimensional networks can be generated in a similar fashion by either generating fractures until the frequency along one-dimensional boreholes is satisfied or by projecting fractures onto sampling planes (e.g., Figure 1) until a two-dimensional density criterion is satisfied.
Fracture density is highly dependent on the distribution of fracture lengths in a model domain where the density at the percolation threshold increases with increasing values of a (Darcel et al., 2003;Reeves et al., 2008b;Renshaw, 1999). This will become apparent in the fracture network examples in the next section.

Flow and advective transport properties of fracture networks
A central theme of this Chapter is that flow and radionuclide transport characteristics of fractured media can be inferred from fracture statistics. The previous section discussed in detail how fractured rock masses can be analyzed according to statistics of fracture orientation, spacing, length, hydraulic conductivity and values of fracture density. This section contains detailed explanations on how these attributes can provide a priori insight into flow and transport properties of a fractured rock mass. We begin by defining network structure, and then discuss how the structure relates to transport characteristics such as shape of breakthrough curves, potential for early arrivals, and variability of individual breakthroughs about the ensemble. Simulations of flow and transport in two-dimensional discrete fracture networks (DFN) with physically realistic parameters are used to illustrate specific concepts.

Network structure and transport
The structure of natural fracture networks is the end result of the complex interplay between stress fields and their anisotropy, mechanical properties of the rock, mechanical fracture interaction and distribution of initial flows in a rock mass. The reliance on probabilistic descriptions of fracture attributes reflects our lack of ability to accurately construct fracture networks based on mapping studies alone. The limited accessability to the network leaves an incomplete understanding of the patterns of fracturing within a rock mass that can often be improved through visual inspection of representative networks generated according to site-specific statistics.
A total of three different network types are generated from two fracture sets with power law distribution of lengths with exponent values in the range 1.0 ≤ a ≤ 3.0, fracture density, 1.0 ≤ ρ 2D ≤ 2.0 m/m 2 , orientations of ±45 • with variability described by κ=20 and a lognormal hydraulic conductivity distribution with log(σ K )=1 (Figures 4-5). Once a network is generated, the hydraulic backbone is identified by eliminating dead-end segments and isolated clusters. This is accomplished in our model through both geometric and flow techniques. The hydraulic backbone represents the interconnected subset of a fracture network that is responsible for conducting all flow and transport across a domain. Hence, analysis of backbone characteristics can provide insight into these processes. The generated networks in this study do not explore the full parameter space for fractured media. However, the wide range of fracture length exponents provides sufficient variability and produces three distinct types of hydraulic backbones. Networks generated with a = 1.0

Hydrogeologic Characterization of Fractured Rock Masses Intended for Disposal of Radioactive Waste
www.intechopen.com produce backbone structures dominated by long fractures (Figure 4), and networks with a = 3.0 produce backbone structures dominated by short fractures (Figure 5). Backbones with a mixture of short and long fractures are produced for networks generated with a = 2.0 ( Figure  6). Another feature of these networks is that density of the network increases from ρ 2D =1.0 m/m 2 to ρ 2D =2.0 m/m 2 as the value of a increases from 1.0 to 3.0. This increase in density is necessary to maintain a percolating backbone. For example, the density values assigned to a = 1.0 and a = 2.0 (ρ 2D =1.0 and 1.5 m/m 2 , respectively) result in a non-percolating networks if used with a = 3.0. Conversely, networks generated with a = 1.0 and ρ 2D = 2.0 m/m 2 (assigned to a = 3.0) produces unrealistically dense networks.  Computation of flow in two-dimensional discrete fracture networks involves solving for hydraulic head at all internal nodes (intersection point of two or more fractures) inside the domain according to Darcy's law (de Dreuzy & Ehrel, 2003;Klimczak et al., 2010;Priest, 1993).
In the simulations presented, boundary conditions are applied to all nodes on the domain to induce flow from top to bottom according to a linear hydraulic gradient of 0.01. Fluid flow through the backbone is then solved iteratively at each node via a biconjugate gradient method (Figure 7). A total of 500,000 conservative (non-sorbing) radionuclide particles are then introduced at the top of the domain, traced through the network, and recorded as they exit the model domain. The velocity along each fracture segment is constant, and radionuclide particles are moved by advection from node to node. Radionuclide particles change direction at nodes proportional to the flux of the down gradient fracture segments. This process is repeated to produce 50 statistically equivalent networks and resultant breakthrough curves for each of the three network types. The radionuclide breakthroughs reflect the characteristics of the fracture networks ( Figure  8). First, fracture length defines the correlation structure of the network. Lower values of a place more probability weight on extreme values of fracture length, and consequently these networks consist of very long fractures that often span the model domain. Radionuclides traveling through these networks only change direction at fracture intersections, and the combination of long fractures and sparse density promotes fast transport. This can be observed in Figure (8) where mean ensemble breakthrough occurs at 1.5 years for a transport distance of 100 meters. Networks with a = 2.0 assign less probability weight on extreme values and produce backbones with a mixture of short and moderately long fracture lengths. This backbone structure requires radionuclides to move through a mixture of short and long fractures and sample a wider variety of velocities in the process of migrating across the model domain. Consequently, mean ensemble breakthrough for networks with a = 2.0 occurs at 11 years, which is approximately 7 times slower than for a = 1.0 (Figure 8) transport occurs for networks with a = 3.0 which consist of backbones dominated by short fracture segments. The short segments and higher density causes radionuclide particles to sample a very large number of fracture segments and velocities in the process of migrating across the domain. The mean breakthrough for a = 3.0 networks occurs at 26 years, which is approximately 17 and 2.4 times slower than the a = 1.0 and a = 2.0 networks (Figure 8).
Second, the breakthrough curves for individual realizations all show a high degree of variability ( Figure 8). In fact, this variability is the reason why breakthroughs are presented in a cumulative form, as the shape of the breakthroughs can range from sharp peaks (more prevalent for networks with a = 1.0) to more broad "Gaussian-like" plumes (more prevalent for networks with a = 3.0). Networks with a = 2.0 represent an intermediate case where breakthroughs for individual realizations can encompass both scenarios. All of the networks exhibit heavy late-time tails which are briefly discussed in the following section. The variability in fracture networks is quantified by comparing deviations of individual realizations from the ensemble. Specifically, the standard deviation of breakthroughs is computed for cumulate probabilities ranging from 0.05 to 0.40 in intervals of 0.05 (Figure 9). This metric is then normalized by the ensemble mean for comparison to the other network types, and termed realization variability. Trends of realization variability are highest for a = 1.0 and lowest for a = 3.0. The largest realization variability values occur for early breakthroughs (0.05 to 0.10 of total mass) for a = 1.0 and a = 2.0. This can be attributed to the fact that these networks, especially a = 1, produce early arrivals that vary greatly in time. Reeves et al. (2008c) studied deviations in spatial distribution of individual solute plumes from the ensemble, and their results support our findings that networks with a = 1.0 are inherently unpredictable, and predictability in transport increases as mean fracture length decreases.

Radionuclide retention mechanisms
The previous section illustrates how statistics of fractured massses can be used to infer advective characteristics of radionuclide transport. In this section, we address the retentive properties of fractured rock masses which include adsorption of radionuclides onto fracture walls and within matrix blocks, and diffusional mass exchange of radionuclides between fractures and low-velocity zones. Adsorption is a chemical process that is related to the mineralogy of the host rock, and cannot be inferred from fracture network statistics. One must also consider that several radionuclides, such as 3 H, 14 C, 36 Cl, 39 Ar and 85 Kr, are conservative and not subject to sorption processes. The dominant retention mechanism for conservative radionuclides is diffusional mass exchange. The retention characteristics of diffusional mass exchange of radionuclides between fractures and matrix blocks can be inferred from fracture network properties and will receive the majority of attention in this section. A short discussion on sorption is included for completeness in the context of the screening framework. A network-scale retention mechanism arising from velocity contrasts among conductive fracture segments and is also introduced for additional consideration.

Diffusional mass exchange
Molecular diffusion promotes the exchange of radionuclides between active fracture flow paths and low-velocity zones such as stagnant zones within fractures, fault gauge and unfractured matrix blocks. Among these low-velocity regions, the influence of diffusional mass exchange of radionuclides in unfractured matrix blocks is much greater than diffusional  Fig. 9. Plot of realization variability in particle breakthroughs about ensemble breakthough for networks generated according to power law length distribution exponents of a = 1.0, a = 2.0, and a = 3.0. Note that realization variability decreases as the value of the power law length exponent decreases.
exchange within fracture stagnant zones and fault gauge (Andersson et al., 2002a;b;Grisak & Pickens, 1980;Neretnieks, 1980;Sudicky & Frind, 1982;Winberg et al., 2000). The reason for this is simple: the time a radionuclide molecule spends trapped in an immobile zone is proportional to the size of that zone, and the scale of matrix blocks is typically on the order of centimeters to meters while within fracture stagnant zones and fault gauge are on the order of micrometers to millimeters. Diffusion into matrix blocks was found to be a major retention mechanism at the time scales of the Tracer Retention Understanding Experiments conducted at the Äspö underground laboratory at both single fracture (Winberg et al., 2000) and block scales (Andersson et al., 2002a;b). For these in-situ experiments, the prevalence of molecular diffusion as a retention mechanism is exhibited by late-time breakthroughs with slopes equal to -3/2 for both conservative and non-conservative radionuclides.
Fracture spacing is a network parameter that denotes unfractured matrix block size and provides a length scale for molecular diffusion. The influence of fracture spacing on contaminant breakthroughs is illustrated in Figure 10 for 3 H. The figure is generated from an exact analytical solution to conservative transport through a series of 40 m parallel-plate fractures with 110 µm apertures subject to dual-domain mass transfer through fractures and matrix blocks (Sudicky & Frind, 1982). The mean age of 3 H breakthroughs versus fracture spacing (2B) is shown in Figure 11. Figures 10 and 11 show that increases in mean breakthrough time for conservative radionuclides are highly correlated to fracture spacing with a power-law (nearly linear) trend. These results indicate that the retention potential for a candidate repository rock mass can be inferred from fracture spacing, and that larger average spacing leads to greater retention due to diffusional mass transfer between fractures and matrix.  (Sudicky & Frind, 1982). Note the dramatic affect of fracture spacing on radionuclide breakthrough due to diffusional mass exchange between fractures and matrix blocks.  Figure 10 with best-fit linear (red) and power-law (green) trendlines. The slightly non-linear trend in mean age versus spacing is more apparent in log-log axes (b). The best-fit power-law trendline is: A = 971(2B) 0.9934 .
Although Figures 10 and 11 are helpful in illustrating and partially quantifying the influence of fracture spacing on radionuclide retention by diffusional mass transfer, it is important to acknowledge that fractured media is inherently much more complex and many simplifying assumptions were made in this analysis. First, cystalline rocks typically have very low porosities and this void space may become disconnected further into the matrix block, giving rise to the concept of maximum penetration depth which describes a cutoff to the distance radionuclides can diffuse into the matrix block (e.g., Neretnieks, 1980). Thus, for crystalline rock masses, average fracture spacing would have to reflect penetration depth which is usually unknown. Second, fracture networks are discontinuous in nature and exhibit power-law fracture lengths which promote a broad distribution of matrix block sizes. For these networks, a "characteristic" or average fracture spacing does not exist and the retention of radionuclides in matrix blocks of different sizes can result in a continuum of exchange rates that are both faster and slower than the what would be predicted using the arithmetically 365 Hydrogeologic Characterization of Fractured Rock Masses Intended for Disposal of Radioactive Waste www.intechopen.com computed average block size. Third, weathering processes at the interface of a fracture and the host rock can create higher porosity rim zones that enhance diffusion (Cvetkovic, 2010). Despite these complicating factors, it can be generally stated that rock masses with greater average fracture spacing have greater potential for radionuclide retention.

Radionuclide adsorption
There may be cases where potential repositories have similar fracture characteristics, the most important being the distribution of fracture length, network density and average fracture spacing, and the sorption capacity of the host rock may serve as a secondary screening criterion. Adsorption of radionuclides is a reversible chemical process that can include both cation exchange and surface complexation. It can often be difficult to distinguish in-situ between retention processes as molecular diffusion plays a role in introducing sorbing radionuclides to adsorptive sites within fractured rock environments. The use of multiple tracers, including conservative and non-conservative with varying sorbing rates, can be used to ascertain the relative contributions of retention from diffusional mass transfer and adsorption (e.g., Andersson et al., 2002b).
Adsorption in fractured rock can occur along fracture walls and within fault gauge and matrix blocks. Sorption along fracture walls is a surface process rather than a volumetric process, like sorption within matrix blocks. Similar to the discussion on molecular diffusion, the TRUE-1 and TRUE Block Scale tests conducted at Äspö Hard Rock Laboratory indicate that sorption in the unfractured rock matrix is a dominant process over adsorption onto fracture walls and fault gauge (Andersson et al., 2002a;b). This is because matrix blocks have a greater surface area and length scale, and hence, a greater number of available sorption sites than fracture walls and fault gauge.
Sorption is typically characterized using equilibrium isotherm models parameterized via laboratory batch experiments containing soil, sediment or crushed rock from a site of interest (e.g., Langmuir, 1997;Zheng & Bennet, 2002). This approach to modeling sorption assumes that rates of sorption/desorption between radionuclides and the sorptive medium are orders of magnitude more rapid than the advective groundwater velocity. However, recent studies have shown that desorption rates of U and 237 Np may yield a 3-4 order of magnitude range of desorption rate constants despite uniform flow fields approximating equilibirum conditions (Arnold et al., 2011;Dean, 2010). We recommend from a cost savings perspective that candidate rock masses be screened on the basis of sorption coefficients determined from laboratory batch experiments. However, the findings of large desorption rate constants questions the usefulness of characterization techniques relying on equilibrium assumptions for use in radionuclide transport models.

Network-scale retention
Numerical investigations have shown that velocity contrasts within fracture segments on the hydraulic backbone can lead to solute retention within the backbone itself (Berkowitz and Scher, 1997;Painter & Cvetkovic, 2002;Reeves et al., 2008c). Painter & Cvetkovic (2002) and Reeves et al. (2008c) found that distributions of inverse velocity, 1/v, a surrogate for retention, exhibited power law decay trends at late times. Values of γ measured from the tail of the distribution of 1/v versus time were in the range 1.1 ≤ γ ≤ 1.8 for Painter & Cvetkovic (2002) and approximately γ=1 for Reeves et al. (2008c). Analysis of the survival function of 366 Radioactive Waste www.intechopen.com ensemble breakthrough versus time for the a = 1.0 networks in this study yielded a value of γ = 1.5 (Figure 12). Values of γ for a = 2.0 and a = 3.0 networks were not computed due to incomplete breakthroughs that leave the tail of the distribution undefined. None of these studies found a correlation between distributions of fracture length and γ; however, all of these studies used transmissivity distributions encompassing several orders of magnitude. It is entirely possible that this type of retention results from the interaction of power law distributions of trace length and distributions of transmissivity that vary by several orders of magnitude, such as lognormal or power law. More work in this area is needed to provide further insight into this retention process in the context of fracture statistics. γ=1.5 Fig. 12. Survival function (inverse empirical cumulative distribution function) of particle breakthrough versus time for the a = 1.0 networks generated in this study (blue) along with best-fit tail trendline (black). The trendline scales according to −1 − γ where a slope of -2.5 corresponds to γ=1.5.

Geologic repository screening framework
The geologic barrier relies on the intrinsic ability of the host rock to limit the release of radionuclides to the biosphere. The discussion up to this point has included both advective and retentive characteristics of fractured media. Desirable qualities for a potential repository rock mass from a radionuclide transport perspective include: low potential for early arrivals, longer bulk breakthroughs, and low variability in predicted results. Simulated particle breakthroughs show a clear trend that rock masses with networks dominated by short fractures exhibit advective transport characteristics that are most desirable for geologic repositories. These networks have the lowest potential for early arrivals, breakthroughs occur at much longer timescales than networks consisting of longer fractures, and the variability of individual realizations about the ensemble are lower than for networks generated using a = 2.0 and a = 3.0. These conclusions are supported by the study of Reeves et al. (2008c).
Retention of radionuclides via molecular diffusion favors networks with large fracture spacings, as retention time for conservative radionuclides in a low-velocity matrix block is proportional to block size. This suggests that the sparsest fracture networks are the most desirable. Given the link between fracture density and the distribution of trace length, networks dominated by long fractures (e.g., a < 2.0) would be the optimum choice.
However, these networks produce breakthroughs that exhibit the earliest arrivals and bulk transport times, and have the greatest variability. It is our opinion that advective transport characteristics should be weighted higher than retention since slower advective transport rates through a network decreases the propensity for fast radionuclide transport and early arrivals. With these network types, retention processes such as molecular diffusion and adsorption will further retard the relatively slow migration rates of radionuclides. Note that molecular diffusion has the ability to significantly retard contaminant migration even for relatively small matrix block sizes (Figures 10 and 11). Networks with identical length distributions can be distinguished from another on the basis of average fracture spacing, a secondary criterion. The rock mass with larger values of average spacing (assuming similar host rock porosity) will have a greater retention potential via molecular diffusion. Sorption characteristics of a rock mass may serve as a tertiary criterion, and this consideration must take into the account the waste itself and the prevalence of non-conservative radionuclides.
This screening framework is intended to be cost effective as flow and transport characteristics of fractured rock masses are inferred solely on the basis of fracture statistics. The selection of an ideal geologic repository for high-level radioactive waste disposal is a complex endeavor that must include considerations of geologic stability (e.g., lack of recent seismic events, volcanic activity, high geothermal gradients), hydrology, rock mechanical strength, geochemistry and proximity to large population centers and natural resources (e.g., National Research Council, 1978). It is assumed that these factors will be considered at some point in the selection process.