Groundwater flow into underground openings in fractured crystalline rocks: an interpretation based on long channels

Hydrogeology Journal, Mar 2017

Rethinking an old tracer experiment in fractured crystalline rock suggests a concept of groundwater flow in sparse networks of long channels that is supported by results from an innovative lattice network model. The model, HyperConv, can vary the mean length of ‘strings’ of connected bonds, and the gaps between them, using two independent probability functions. It is found that networks of long channels are able to percolate at lower values of (bond) density than networks of short channels. A general relationship between mean channel length, mean gap length and probability of percolation has been developed which incorporates the well-established result for ‘classical’ lattice network models as a special case. Using parameters appropriate to a 4-m diameter drift located 360 m below surface at Stripa Mine Underground Research Laboratory in Sweden, HyperConv is able to reproduce values of apparent positive skin, as observed in the so-called Macropermeability Experiment, but only when mean channel length exceeds 10 m. This implies that such channel systems must cross many fracture intersections without bifurcating. A general relationship in terms of flow dimension is suggested. Some initial investigations using HyperConv show that the commonly observed feature, ‘compartmentalization’, only occurs when channel density is just above the percolation threshold. Such compartments have been observed at Kamaishi Experimental Mine (Japan) implying a sparse flow network. It is suggested that compartments and skin are observable in the field, indicate sparse channel systems, and could form part of site characterization for deep nuclear waste repositories.

A PDF file should load here. If you do not see its contents the file may be temporarily unavailable at the journal website or you do not have a PDF plug-in installed and enabled in your browser.

Alternatively, you can download the file locally and open with any standalone PDF reader:

Groundwater flow into underground openings in fractured crystalline rocks: an interpretation based on long channels

Groundwater flow into underground openings in fractured crystalline rocks: an interpretation based on long channels John H. Black 0 1 Nicholas D. Woodman 0 John A. Barker 0 0 Department of Engineering and Environment, University of Southampton , University Road, Southampton SO17 1BJ , UK 1 In Site Hydro Ltd , East Bridgford, Nottingham NG13 8NW, UK 2 Rethinking an old tracer experiment in fractured crystalline rock suggests a concept of groundwater flow in sparse networks of long channels that is supported by results from an innovative lattice network model. The model, HyperConv, can vary the mean length of 'strings' of connected bonds, and the gaps between them, using two independent probability functions. It is found that networks of long channels are able to percolate at lower values of (bond) density than networks of short channels. A general relationship between mean channel length, mean gap length and probability of percolation has been developed which incorporates the well-established result for 'classical' lattice network models as a special case. Using parameters appropriate to a 4-m diameter drift located 360 m below surface at Stripa Mine Underground Research Laboratory in Sweden, HyperConv is able to reproduce values of apparent positive skin, as observed in the so-called Macropermeability Experiment, but only when mean channel length exceeds 10 m. This implies that such channel systems must cross many fracture intersections without bifurcating. A general relationship in terms of flow dimension is suggested. Some initial investigations using HyperConv show that the commonly observed feature, 'compartmentalization', only occurs when channel density is just above the percolation threshold. Such compartments have been observed at Kamaishi Experimental Mine (Japan) implying a sparse flow network. It is suggested that compartments and skin are observable in the field, indicate sparse channel systems, and could form part of site characterization for deep nuclear waste repositories. Conceptual models; Crystalline rocks; Hydraulic properties; Numerical modeling; Radioactive waste repositories - Many geoscientists are interested in how groundwater flows through sparsely fractured crystalline rock, none more so than those considering constructing nuclear waste repositories in such formations. During the 1970s and 1980s, several countries initiated programmes aimed at disposing of radioactive waste at depth, some of which included the construction of underground research laboratories (URLs; Witherspoon and Bodvarsson 2006). One of the major intended activities in these URLs was to study groundwater flow and radionuclide transport within poorly transmissive fractures at roughly repository depth (say ∼500 m below surface). Early tracer experiments at Stripa Mine, Sweden (Abelin et al. 1985, 1991a), attempted to inject tracers into the high head gradients induced by the presence of the mine and then recover them in the underground opening. This tended to yield low rates of tracer recovery so subsequent tracer experiments in URLs were usually located some distance laterally into the rock in artificially created ‘point sink’ configurations (Sawada et al. 2000; Elert and Svensson 2001; Andersson et al. 2002; Einsiedl and Maloszewski 2005). In addition, these borehole to borehole tests were always sited in major fractures and fracture zones with concomitantly higher values of transmissivity compared to surrounding rock: say ∼10−7 relative to ∼10−10 m2/s. Obviously, these experimental choices were made in order to ensure positive results since migration rates within low transmissivity fractures for even conservative tracers over small distances require many years. However, one could reasonably question how much results from fracture and fault zones reveal about flow and transport in what is often referred to as ‘background rock’ (i.e. the average rock between fault and fracture zones and referred to as ‘fracture domains’ by Follin and Hartley 2014). The origin of such zones ought to produce a morphology and mineral content that differs from ordinary fractures in background rock. Indeed, Andersson et al. (2007) concluded after the extensive series of TRUE tracer tests at spö Hard Rock Laboratory, Sweden, that uncertainty concerning the tracer flowpath was greater than uncertainty about the water −rock interaction. Also that background fractures were about ten times more retentive than the main fractures involved in the tracer experiments. Notwithstanding that the objective of the tracer tests has been to determine radionuclide transport properties, it is the premise of this paper that an early, possibly considered inconclusive, tracer test yielded a valuable insight into how groundwater flows through background rock. Based on observations from the aforementioned tracer test, a particular concept of long channels and sparse intersections is proposed. The concept is encapsulated within an innovative lattice network code that is shown to reproduce some of the features noted or measured in other URL experiments. The 3D (Migration) Experiment The ‘Three-Dimensional (Migration) Experiment’ was undertaken at the Stripa Mine between the end of 1984 and early 1987. The experiment was based on the assumption that a drift at atmospheric pressure, hundreds of metres below the water table, would provide a suitable sink for surrounding groundwater. Then, if different tracers were injected into the local groundwater system at various locations, they could all be recovered eventually as inflow to the drift (Abelin et al. 1991a, b). An experimental drift (roughly 75 m long and 4 m diameter) was located at least 100 m from existing mine openings and connected to the old mine galleries by a 140-m-long access drift. Two ‘side drifts’ about 10 m long, were also excavated about 10 m before the end of the drift, thus forming a cruciform-shaped experimental gallery (Fig. 1). A major innovation of the experiment was the development of a polythene-sheet catchwater system so that inflow could be identified as emanating from individual rectangles marked on the upper walls and ceiling of the experimental drift. Three hundred and seventy-five polythene sheets, each a rectangle of 2 m2, were attached to the upper rock surfaces using resin glue which enabled each area of the grid to yield a measurement of inflow (Fig. 1a). Inflows varied greatly from sheet to sheet; less than half of them gave any measurable inflow and 11 sheets accounted for 50% of inflow. About half of the grid-based inflow emanated from one of the side drifts. The results clearly showed that natural inflow is limited to a few areas with large dry areas between. For the next stage of the experiment, three 70-m-long tracer injection boreholes were drilled vertically upwards into the ceiling of the main drift (Fig. 1b). Nine different tracers were injected in short (2.5 m), discrete, packer-isolated zones at distances between 10 and 56 m directly above the crown of the drift. They were injected for periods between 60 and 95 weeks followed by 25 weeks of continued monitoring. Only seven tracers were recovered in the experimental drift, of which one, bromide, incurred analysis difficulties that negated attempts to calculate cumulative recovery. The distribution of tracer inflow was as noteworthy as the distribution of water inflow. In contrast to the water inflow that was highest into the NE arm of the crucifix-shaped excavation, all tracers that were recovered at all arrived in the middle of the experimental drift (Fig. 1b). They ignored the wettest part of the drift located in one of the ‘arms’ (Fig. 1a). Also surprising, was the recovery of all four tracers from borehole C given its proximity to the high-inflow ‘arm’. Perhaps equally surprising was the absence of recovery of the tracers from the two more distant injection points in borehole B bearing in mind that they both appeared to show a pronounced head response to the experimental drift. A final conundrum was posed by the subsequent discovery, in significant quantities, of the tracer Eosin Y from borehole A when excavating a connecting drift about 150 m away from the injection point. In summary, the 3D Migration Experiment revealed the extremely uneven distribution of inflow to underground excavations together with flow paths that ignore the concepts underlying groundwater flow in a porous medium. Abelin et al. (1991b) concluded that the results indicated that much flow takes place in channels that interact infrequently with neighbouring channels. Even though Abelin et al. (1991a, b) conclude that flow occurs via channels within fractures, they do not attempt to gain further insight by examining the pattern of mixing exhibited by their experiment. The pattern is enlightening because most catchwater sheets that yielded tracer, yielded mixtures from several injection points, and not always from the nearest injection points. Also, not all possible mixtures were represented. It seems reasonable to assume that mixing must occur at certain points where actively flowing channels meet others. Based on the recovery of the six tracers, it is possible to reproduce the observed mixtures based on only eleven ‘mixing nodes’ (Fig. 2). The incidence of mixtures containing tracers from injection points in boreholes A and C, many metres apart, suggests that the mixing nodes are relatively close to the catchwater sheets, possibly within an excavation damaged zone (EDZ). This in turn suggests that the channelized flowpaths from the injection points to the Fig. 1 The 3D (Migration) Experiment: a natural inflows of groundwater; b the measured distribution of tracer inflows from nine injection points in the three boreholes A, B and C (after Figs. 5 and 13 in Abelin et al. 1991a) The sparse channel network concept The concept of channel-based flow is not new. However, since Rouleau and Gale (1985) report that fractures at Stripa Mine are seldom more than 0.5 m apart, the concept of long channels proposed here means that channels in background rock cross many fracture intersections without branching. The concept of sparse long channels with few intersections agrees with the frequent observation in tunnels and drifts that inflows are sporadic and wet patches are separated by large regions of dry tunnel/drift surface (Tsang and Neretnieks 1998; Tsang et al. 2015). An implication of isolated points of drift inflow is that groundwater must converge to those points and that the convergence must cause head loss proportional to the ‘dimension’ of the convergence (Barker 1988). Figure 3a shows theoretical head distributions towards a drift based on a 4 m diameter drift in a groundwater environment with a 200 m head Fig. 2 Pseudo three-dimensional representation of the movement of tracers in the 3D Migration Experiment invoking imaginary ‘mixing’ and ‘dispersal’ nodes (constructed from information contained in Fig. 13 of Abelin et al. 1991a) Fig. 3 Head versus distance from the centre of a drift. a Calculated values and log linear regression lines for ideal systems with zero skin. b Results from packer isolated sections in the radial boreholes of the Macropermeability Experiment at Stripa (based on Wilson et al. 1981) boundary located 50 m from the drift centreline. Figure 3a uses a logarithmic x-axis so that, based on porous medium theory, cylindrically convergent flow should plot as a straight line. If the convergent flow is of a higher dimension than cylindrical (i.e. >2D), the line of log-linear regression projected to the drift wall will intersect it at positive values of head. Most hydrogeologists would interpret this configuration as a positive skin, a region of locally reduced hydraulic conductivity relative to ‘average’ rock further from the drift (Fig. 3a). The general concept and mathematics of skin and its application to a drift deep below the water table is summarized in the Appendix. Recognizing and estimating skin in fractured crystalline rocks is difficult because of the inherent variability of the transmissivity of adjacent fractures. Nevertheless, if enough points around a drift (actually short straddle packer intervals in nearby boreholes) are measured, then a distance–drawdown relationship is usually apparent. Figure 3b shows such a set of measurements from the Macropermeability Experiment at Stripa Mine (Wilson et al. 1981). Similar results have been reported for other experiments at Stripa, such as the Simulated Drift Experiment (Black et al. 2006), and for other URLs such as the KD90 drift at Kamaishi in Japan (Chijimatsu et al. 1996). The idea proposed here is that observations of isolated patches of inflow on drift surfaces, positive skin effect and channelling effects in tracer experiments in ‘background rock’ are all linked to flow in long channels with infrequent junctions. It was decided to validate the concept by reproducing the ‘apparent’ skin effect in an appropriate model. The intention was to construct a numerical model of flow to limited inflow points to a drift that would, when using parameter values constrained by the Stripa experiments, yield similar values of apparent skin effect to those that were observed. A novel lattice network model Choice of model type Studies of flow in fractured crystalline rocks have employed a wide variety of different conceptual and numerical models (Selroos et al. 2002). Apart from porous medium and discrete fracture network models, lattice network models have sometimes been used to understand flow and solute transport behaviour in three dimensions. Moreno and Neretnieks (1993) and Margolin et al. (1998) both used a cubic lattice model characterized by the simple probability of any single bond being present. This approach, here termed the ‘classical lattice model’, has the characteristic that long channels (continuous lines of bonds) are associated with short gaps and vice versa. Since the experiments at Stripa described earlier suggest long channels separated by large gaps, the classical lattice model is considered inappropriate. A novel lattice network concept has therefore been developed that allows the mean lengths of the channels and gaps to be independent of each other. This ‘two-parameter bond lattice model’ has been implemented primarily as a 3D steady-state flow code named HyperConv. This name derives from the term ‘hyper-convergent flow’ introduced by Black et al. 2006 to represent flow converging to a point sink within cylindrical boundary conditions. It was felt that important characteristics of the code should include the ability to calculate flow and head quickly in order to run multiple realizations. HyperConv does not attempt to reflect rock or fracture structure except in as much as it is three-dimensional and channels are organized on an orthogonal grid. Hundreds of realizations are generated for two fixed parameters (e.g. mean channel and gap lengths). For all realizations, flow through each bond is determined for one of two types of boundary conditions and the key measures of bulk hydraulic conductivity and skin effect are evaluated. The measures are then averaged over all percolating realizations of a particular ensemble (i.e. of 100 realizations). The network is formed by creating bonds between the nodes of a 60 × 60 × 60 orthogonal lattice. A channel is conceived as a group of bonds in a continuous row. Two probabilities are used to generate channels along lines of nodes in each coordinate direction: 1. PN (abbreviation of PNew) is the probability of starting a new channel if the previous bond (within a particular line of nodes) is absent. 2. PA (abbreviation of PAgain) is the probability of continuing a channel to a subsequent node. This method of network generation gives rise to channel lengths and gap lengths (in the x, y and z-directions) that both have a geometric distribution (Table 1). The mean and variance are given in Table 2. The lengths are all scaled to the length of a bond (single sub-channel) which is given unit length. It should be noted that any two of the parameter set {L, G, PON, PN, PA} are sufficient to characterize the new lattice model and that the ‘classical lattice model’ is a special case where PON = PN = PA or, equivalently, 1/L + 1/G = 1. Applying boundary conditions and assigning parameters The basic channel network is produced in the form of a cube. Once this has been generated, different boundary conditions can be applied, namely face-to-face, and cylindrical (Fig. 4). Face-to-face boundary conditions are used to evaluate percolation thresholds and ‘compartmentalization’, whilst cylindrical boundaries are used to evaluate skin (in the context of a tunnel or drift). For face-to-face flow, four of the six faces are set as no-flow boundaries confining flow between the remaining two opposing faces at fixed head. Cylindrical boundary conditions are achieved by removing surplus nodes in the region of four parallel edges to form a cylinder and setting an inner cylindrical boundary at (zero) fixed head (Fig. 4). The modelling was intended to simulate observations at Stripa URL so values were assigned to boundary conditions with the Stripa measurements in mind. The main choice concerns the ‘equivalent’ distance between adjacent nodes within the lattice and hence the length of a ‘sub-channel’. A value of 1.5 m was chosen so that each 2 m × 1 m catchwater sheet of the 3D Migration Experiment (Fig. 1) could receive inflow at maximum model density (PON = 1). This in turn meant that the external boundary of the cylindrical version of the network was 44.25 m from the centre of the network. A value of 200 m head was chosen for the external boundary as this appeared to be the undisturbed background value near to the 3D Migration Experiment (Abelin et al. 1991a). The inner boundary was a 4-m diameter central drift connected to the atmosphere at 0 m head. Inflow to the Macropermeability Experiment, in conjunction with the near-drift hydraulic gradient, yielded an estimated value of bulk hydraulic conductivity of about 2 × 10−10 m/s. Within the HyperConv code, channel conductances (Cs) are selected randomly from a lognormal distribution. All bonds (subchannels) comprising a channel have the same conductance (Cs). The bond conductance values are based on a constant base (geometric mean) value of conductance, Csg. So, where u is a randomly generated unit normal variable, and σ is the standard deviation of log10 (Cs). A value of Csg = 4.5 × 10−10 m3/s was used to yield a target hydraulic conductivity of 2 × 10−10 m/s for a fully occupied grid. By setting the standard deviation of the log normal L = 1/(1 − PA) G = 1/PN Distribution of lengths of: Prfl ¼ ng ¼ PnA−1ð1−PAÞ ¼ ðL−1Þn−1=Ln Prfg ¼ ng ¼ PNð1−PNÞn−1 ¼ ðG−1Þn−1=Gn Relationships between probability generators, mean lengths and the bond probability PON L ¼ 1−1PA PA ¼ 1− L1 distribution (σ) of log10(Cs) to zero, a constant conductance network can be produced. In this way it is possible to compare and therefore distinguish the significance of network structure on flow from the effect of varying channel conductance. Solving the flow equation The last step is to solve the flow equation for the thinned channel network. The network has mass balance applied at each node with a linear, Darcian, relationship between flow and head difference between each bond pair. The resulting linear equations are solved for the nodal heads using different routines depending on bond density. Where the network is sparse and symmetric, a pre-conditioned conjugate gradient approach is favoured. For poorly connected systems (less than 30% of the nodes connected) a general purpose sparse matrix solver with a partial Cholesky pre-conditioner is used on a system containing only the connected nodes. For fuller systems, the full lattice is solved using a specialized block-structured pre-conditioner (see Black et al. 2006 for more detail). Finally, the channel network is assessed to identify which nodes have significant flow in them and which contain negligible flow. All the nodes containing significant flow (i.e. values greater than 10−16 m3/s) are termed ‘flowing’ nodes. The rest are termed ‘connected’ nodes. The evaluation of skin is based on finding the best-fit loglinear regression line through the heads calculated at the nodes Fig. 4 Schematic of the sequence of operations involving the thinning of the HyperConv network during the application of different boundary conditions and the solution of the flow equation of the network in the manner identified in Fig. 3b. The calculation is of ‘skin effect’ (hskin), and has units of metres of water head. It comes in two versions, hsckin and hfskin that refer to the calculation based either on all connected nodes (superscript c) or only on the flowing nodes (superscript f). The values for all connected nodes, hsckin, is used in the results below because that is the more likely type of value measured in the field, albeit at less frequency at greater radial distance. HyperConv was designed to construct a relatively complex network, solve the flow equation and derive a value for skin effect quickly enough that it could be run many times. The outcomes of many realizations can then be averaged to yield a single value appropriate to the single set of parameters used to generate the ensemble (of realizations). A key question concerned how big the ensemble needed to be. Based on three runs (σ = 0, 1 and 2), each of about 4,000 realizations, and using the fact that the variance of a mean is inversely proportional to the sample size, it is estimated that for 100 realizations the relative error is typically between 2 and 5%. Therefore an ensemble size of 100 was chosen as a reasonable compromise between accuracy and computer run time (Black et al. 2007). A fundamental property of any network is whether or not it percolates (i.e. conductive elements connect across the network). Non-percolating realizations are generally of little interest for flow or transport. This study used a series of ensembles of 100 realizations for five values of mean channel length, but with different sizes of intervening gaps. Starting with large gaps (i.e. low density), only a very few realizations within an ensemble would percolate. As gap size was decreased increasing numbers of realizations would percolate until all percolated. The results show that the average number of percolating realizations is dependent on mean channel length as well as network density (Fig. 5). More importantly, systems of long channels continue to percolate at lower values of network density than short channels. Directly linked to this finding is that, of the systems examined using HyperConv, classical lattice networks require the highest values of density to percolate. Based on this outcome, the percolation behaviour of a much larger million-node cubic lattice designed to evaluate only connectivity (i.e. without the ability to solve for flow) was evaluated. The range of behavior of a ‘two-parameter bond lattice’ type model is shown in Fig. 6 and spans all Fig. 5 The dependence of percolation probability on channel length and network density. L = average number of bonds (sub-channels) per channel. The value of probability of percolation of 50% is commonly referred to as the average percolation threshold or critical density (Harter 2005) combinations of average channel and gap length (i.e. all PA and PN values). The boundary between percolating and nonpercolating combinations of channel and gap length can be approximated by the equation Fig. 6 The combinations of channel and gap lengths that percolate and do not percolate and their relationship to bond density (PON) Figure 6 shows the two regions where percolation does and does not occur, though the boundary is only sharp for infinite networks. It also indicates that any model with a bond density, PON, greater than about 0.3 should percolate and that, as PON falls below this value, percolation only occurs for increasingly long channels and gaps. Most of the important results from the HyperConv simulations occupy this low density region. It can be seen that the possible combinations of average channel and gap lengths for a classical lattice network are limited to the single line running diagonally across Fig. 6, corresponding to 1/L + 1/G = 1. The intention of the modelling was to investigate the effect of channel length and sparseness on the derivation of values of skin effect. Margolin et al. (1998) have previously shown that sparse flow systems can be produced by classical lattice networks based on variable bond conductance. In order to avoid confusion with effects caused by variable bond conductance, all the results reported here (except for some ensembles investigating variable conductance) concern systems comprising sub-channels that have the same value of conductance. The variation with density in Fig. 7 is surprising. A general trend from extreme values at very low bond density to zero skin when the lattice is full (i.e. high bond density) might have been expected. Instead, while the ensemble results of some shorter channels exhibit a reasonably gradual evolution from large negative values of skin at low bond density towards zero at high density, some longer channels exhibit a peak of positive skin at intermediate values of density. Fig. 7 The dependence of skin effect on channel length and network density. Each point represents the average value of an ensemble of 100 realizations or the number of realizations that percolated: reported as the number beside data points. Where numbers are absent, all 100 realizations of a given ensemble percolate Considering the results in more detail, it can be seen that all networks yielded negative skin effect over at least some range of density values; that is they indicated a region of apparently increased permeability in the near-tunnel region. However, only some more dense ensembles of longer channels yielded a positive skin effect (i.e. apparent reduced permeability near the tunnel). The division between networks with a region of positive skin effect and those without, appears to occur around networks of 6.7 bonds average channel length (10 m). It should probably be borne in mind that the modelled system includes a 4-m diameter drift as its inner boundary and that this is of similar size to the 10-m scale dividing the two types of behaviour. Another point of note from the results is that classical lattice networks show no tendency to yield positive skin effects in the conditions applied here. In line with the percolation results shown in Fig. 5, the expected trend of percolation percentage decreasing with bond density is maintained with a marked decline for PON less than 0.04. With decreasing numbers of realizations percolating, the variability of their average increases as evidenced by the overlapping data curves at large negative values of skin effect (Fig. 7). The relatively sudden transition from negative to positive skin effect exhibited by networks of longer channels is a region of particular interest. It was decided to investigate this transition by examining in more detail the distance/drawdown behaviour of three ensembles of 10bonds average length per channel (PON ≈ 0.02, 0.04 and 0.07 as identified in Fig. 7). Head was calculated as a function of distance from the inner boundary for all connected nodes within every percolating realization of each of the three ensembles. Values were placed in annular bins of equal thickness and an average Fig. 8 a Head versus distance from the centre of a drift of 4-m diameter for percolating ensembles based on 15-m average channel length and three values of density. b Section across the axis of the inner boundary of the HyperConv model showing the configuration of potential inflow points. Note that the lattice extends out in all directions beyond the inner region depicted derived (Fig. 8a). In all cases, the value of conductance of the bonds is kept constant (σ = 0) so that all outcomes are a function of geometry. It should be pointed out that these values of density represent average ‘gap’ distances between channels ranging between 125 and 500 bonds (i.e. beyond the boundaries of the model). The variation of head with the logarithm of distance from a drift centreline is plotted in Fig. 8a against a background of 1D, 2D and 3D flow in a porous medium with zero skin. The straight-line sections extending from the outer boundary for all three cases indicates reasonable adherence to cylindrically convergent flow essentially imposed by the outer boundary. However, as the three types of flow system approach the drift wall, the least dense system diverges from cylindrical behavior at greatest distance and yields an excellent fit to the model of Barker (1988) with a flow dimension of 1.6. The most dense system of the three (i.e. PON ≈ 0.07) adheres to cylindrically convergent flow until very close to the drift where almost half of the overall head loss occurs within 2 m of the drift wall. If the procedure summarized in the Appendix is applied, the data yields an apparent positive skin effect only for the densest of the three realizations. It is useful to consider some practical dimensions of the lattice network at the point of inflow. Eight bonds (potential inflow points) emanate from each lattice plane perpendicular to a 4-m diameter drift (Fig 8b), so the probability of one or more of these bonds being ‘on’ is (1 – PON)8 and from this it follows that the average distance between perpendicular planes with inflow points to the tunnel is For the results in Fig. 8a) this gives distances of: 2.1, 3.7 and 6.8 lattice intervals (3.3, 5.6 and 10.3 m), respectively, increasing with decreasing density. It is notable that for the two higher density results in Fig. 8a, these intervals are reasonably consistent with the radial distance at which flow ceases to be cylindrical, which is analogous to the flow pattern for a partially penetrating well. Although the objective of constructing a sparse long channel model was to simulate skin effect, it is valuable to examine the relationship between average channel length, channel density and hydraulic conductivity. The HyperConv model is capable of calculating overall network conductivity using different boundary conditions: face-to-face and radial flow to a line sink (Fig. 4). The results shown in Fig. 9 are much as expected: as bond density decreases, hydraulic conductivity decreases. In line with skin effect, the impact of density is most pronounced where average channel lengths are shortest. It is noteworthy that the percentage of percolating realizations drops below 100% at roughly the same value of hydraulic conductivity (about 2.5 orders of magnitude less than the full lattice value) for all lengths of channel. Effect of spatially variable channel conductance on skin effect and hydraulic conductivity In the investigations reported in the preceding sections, the conductance of the bonds of the networks in all the realizations has the same value by setting the standard deviation of the log normal distribution of log10 (conductance) to zero (σ = 0). This means that all the effects observed in the model results thus far are the product of network geometry alone. Obviously, natural systems cannot be expected to be composed solely of constant transmissivity channels, so variance in conductance was introduced by setting σ greater than zero on the set of ensembles based around an average channel length of 10 bonds (i.e. average length = 15 m and the same as the investigation of transitional behaviour above) (Fig. 10). Fig. 9 The dependence of hydraulic conductivity on channel length and network density. Each point represents the average conductivity of the percolating realizations from an ensemble of 100 realizations The skin effect values (Fig. 10a) are displaced upwards but by a relatively small amount given that the range of conductances for one standard deviation differ by a factor of 30. This suggests that the results for equal conductance channels (σ = 0) are of reasonably general applicability. The outcome for hydraulic conductivity can be understood by reference to the analytical result obtained by Rubin (2003) for an uncorrelated, spatially stochastic, lognormal permeability field: where Kg is the geometric mean of K, m is the geometrical (network) dimension and σy2 is the variance of ln(K). At high density, the network dimension, m, can be expected to be at a maximum of 3 which yields a value for Keff of roughly 4 × 10−10 m/s via Eq (4) given Kg = 2 × 10−10 m/s. In other words, flow occurs preferentially via the more conductive channels. Thus, it must be assumed that, for any given density, the number of channels containing the bulk of flow is fewer when conductance varies than when it is held constant. This may also be the cause of the increase in skin effect at maximum density because flow continues to converge to the more conductive channels intersecting the inner boundary. Based on Eq. (4), when the network dimension (m) = 2, Keff = Kg which agrees well with Fig. 10b) for PON ≈ 1. The transition from positive to negative values of skin effect (at PON ≈ 0.04) also marks the boundary between all Fig. 10 Results for networks averaging 10 bonds per channel but with variable conductance: a skin effect versus density, and b hydraulic conductivity. Numbers beside data points indicate percolating percentage and apply to both figures (a). A ‘spanning’ channel is one that extends the entire distance between the inner and outer boundaries of the model realizations percolating and decreasing percentages percolating (Fig. 10). At network densities less than PON ≈ 0.04, results are best understood in the context of a single individual channel extending the entire distance between inner and outer boundaries. The parameters appropriate for a single line of ‘flowing’ nodes (i.e. 29 bonds) are: −62.27 m (log-linear regression to 29 equally spaced head steps of 7.1 m) It can be seen in Fig. 10b that once node density decreases to 0.01 and percolation frequency is at 21%, the value of conductivity of the equal conductance networks decreases to the value for a single (‘spanning’) channel intercepting the inner boundary. For the spatially variable networks where σ = 0.75, each single spanning channel can vary within the limits of the log normal distribution and an effective hydraulic conductivity given by the average conductance can be expected. For a lognormal distribution the average is Kg exp(σy2/2) giving Kmin = 2.5 × 10−13 m/s for σ = 0.75 which agrees reasonably well with Fig. 10b). Bounding the model behaviour A ‘semi-generic’ approach It is possible to derive semi-generic sets of results by linking the hand-calculated minimum results (for a single spanning line of bonds) to those calculated by the lattice network model using estimated results (Fig. 11). Some of the lattice model results in the region where only small percentages of ensembles percolate (and therefore results are variable) have been manually smoothed (compare Figs. 7 and 9 with Fig. 11a, b). The resulting diagrams, Fig. 11a, b are not generic since the axes are not dimensionless. They are termed ‘semi-generic’ because it is considered that realizations based on different parameters from the set calculated will exhibit similar relative behaviour: in particular, the peak values of positive skin just above the percolation threshold for longer channels. Fig. 11 The relationship between a skin effect and b hydraulic conductivity, and bond (network) density shown in semi-generic form to indicate generalized behaviour. Key applies to both parts of the figure The semi-generic form of the skin effect results (Fig. 11a) indicates a division into two groups, those short channels that always yield a negative skin effect and the longer ones that exhibit an interval of positive skin effect. It is surmised that the value of this division is linked to the particular size of the inner boundary (4 m diameter) relative to the average length of the channels within each set of results (10 m in this case). The envelope of possible values of hydraulic conductivity (Fig. 11b) is bounded by the density of a single line of bonds between the inner and outer boundaries (i.e. 1.88 × 10−4) and the density of a full lattice (i.e. 1). The upper edge of the envelope is formed by the density–conductivity relationship of multiples of single lines of connected channels. One of the curiosities of the results set is that a single line of channels does not yield the minimum value of hydraulic conductivity. Instead, the minimum values of hydraulic conductivity are produced by a small number of very tortuous flow paths. It is not clear whether short channel networks (i.e. less than three bonds average length and including classical lattice networks) exhibit such behavior because the probability of percolation declines so rapidly with decreasing bond density that numerical results were unobtainable in this study. The skin effect results for head versus distance (Fig. 8) and the explanation of the results for spatially variable hydraulic conductivity Fig. 10b) indicate the possible usefulness of a dimension-based interpretation. A conceptual model involving three different dimensions is suggested. The channels of the percolating lattice form a network which has a network dimension, Dnet (‘m’ in Eq. 4) which decreases from 3 at the highest densities (PON approaching unity) to 1 when the density is below about 1%. As density decreases, the number of inflow points to the tunnel decreases, which increases the flow convergence close to the tunnel wall. This convergence yields positive values of skin and is here termed the ‘boundary-condition dimension’, Dbc, close to the tunnel. The skin factor closely relates to the flow dimension, Dflow (Fig. 10a). Relating these three dimensions to two parameters of the lattice model, PON and L, requires some speculative assumptions. If it assumed that the network dimension is linearly related to the number of intersections per channel, B, which, for the two-parameter lattice network model is: Then, combining Eq. (5) with relations in Table 2, and bearing in mind that a 1D system will have zero intersections, Dnet can be defined in terms of PON and L: where Fnet is a constant. Similarly, it is proposed that the boundary condition dimension is linearly related to the average interval between inflow planes along the tunnel (Δinflow as defined by Eq. 3). Since, if each plane of nodes perpendicular to the axis of the inner boundary has an inflow point to the tunnel and results in 2D flow, then a value of Δinflow greater than 1 should yield a flow dimension greater than 2; hence, Eq. (7) imposes cylindrical flow when the inflow interval is equal to the lattice grid spacing. 1−ð1−PONÞ where Fbc is a constant. A reasonable supposition is that the flow dimension is limited by both the network dimension and the boundarycondition dimension, so: Since the intention is to relate apparent skin effect to flow dimension, Dflow, the results from HyperConv of skin effect for varying density, PON, and average channel length, L, are used for appropriate scaling. It seems reasonable that, based on Fig. 7, the maximum value of skin effect for L = 20 at PON = 0.03 should equate to a flow dimension of 3. This yields values for Fnet of 1.5 and for Fbc of 0.3 and gives the variation of flow dimension with density and channel length shown in Fig. 12. A better correspondence could not be expected given the simplicity of the model, the use of averaged parameters and functions, and the fact that the results in Fig. 7 are averages over percolating realizations alone. Although Eqs. (6)–(8) do not include the size of the inner boundary directly, the outcome involves Eq. (7) which is based on eight possible inflow points connecting the lattice to the inner boundary. This means that the channel lengths are a function of tunnel diameter and the interpretation may have general application. The implications of a long-channel concept The idea of long channels between bifurcations does not, at first, appear to be a particularly radical concept; however, it does have many implications that are not immediately obvious. Firstly, the HyperConv modelling clearly shows that long-channel systems are able to percolate at lower values of network density than networks consisting of shorter channels. Assuming that once water finds the shortest pathway between regions of higher and lower pressure, it will not seek extra pathways unless there is a change in boundary conditions or the original pathways alter due to external processes (e.g. seismicity, mineral precipitation, solution, etc.). This should leave most groundwater in background crystalline rock flowing within a channel network at a density just above the percolation threshold for an average channel length available within Fig. 12 Variation of flow dimension with density, PON, and channel length based on Eqs. (6)–(8) that particular rock. This should mean that the network containing flowing groundwater is sparse. Channel sparseness has its own implications. HyperConv simulations show how, close to the percolation threshold, flow within a sparse long channel network converges to pass through a few key conductors that can be regarded as ‘chokes’. The same choke effect is reported in previous variable conductance lattice network studies (Margolin et al. 1998) and in a fracture network study by Robinson (1984). Within such a flow system containing chokes, the bulk of head loss within the system will occur within the chokes, leaving the network between the chokes with little head loss as a consequence. These patches of similar head adjacent to patches of quite different head are termed ‘compartments’. These ideas of sparseness, chokes and compartments are discussed in the following. The main reason that channel networks such as those modelled here using the HyperConv lattice network model, emulate the skin effect behaviour observed at Stripa is their sparseness. It is difficult to appreciate the impact of sparseness without a visual cue. Figure 13 shows a single representative 2D slice through two HyperConv-like 3D lattice networks. Figure 13a is of a classical network at critical density (i.e. 50% of realizations percolate), whilst Fig. 13b is of a long channel network, also at critical density. The channels in Fig. 13b are 5× longer than those in Fig. 13a and the bond density is also 5× less. Neither slice has a direct connection across the network from face-to-face so relies for percolation on connections perpendicular to the 2D slice. There is ample field evidence of sparseness in the form of isolated patches of inflow to drifts at Stripa URL—e.g. the Buffer Mass Test area of Pusch et al. (1985) and the site characterization and validation drift of Harding and Black (1992)—into the Kymmen (Moreno and Neretnieks 1989) and Bolmen tunnels in Sweden (Neretnieks 2006), and into tunnels in the Aar and Gotthard massifs of Switzerland (Masset and Loew 2010). At Pocos de Caldas uranium mine in Brazil, patches of oxidized pyrite associated with groundwater flow are observed together with channels of 100 m or more in length (Chapman et al. 1992). Chokes and compartments The model of Margolin et al. (1998) was quite similar to HyperConv, though being of the classical lattice network type (see Fig. 6) it lacked the independent spatial correlation of bonds and gaps that HyperConv is able to simulate. Even so, Margolin et al. (1998) found that by either reducing the density of channels or increasing the variation of channel conductance they could reduce the number of percolating flow paths in a region to a small number. They also reported that chokes have a greater effect on overall network conductivity as network density decreases towards the percolation threshold (see Fig. 8 of Margolin et al. 1998). Hunt (2005), in his review of percolation theory applied to hydrogeology, stated that, Bcritical resistances on the optimal path control the entire field of potential drops.^ In other words, a sparse channel network should be associated with regions (i.e. clusters of channels) of roughly equal head (termed ‘compartments’) separated by head drops that reflect the head losses through the dominant chokes in the flowing system. Fig. 13 2D slices from two HyperConv-type lattice networks at critical bond density: a a classical network (average channel length (L) = 1.331 bonds, average gap length (G) = 4.02 bonds), in which lines A and B identify two of the many continuous lines of gaps across the entire network; b a long channel network (L = 10, G = 200 bonds) Almost all of the results reported here are based on constant conductance (σ = 0) networks, where any choke-like behaviour must result from the convergence of several channels towards a single channel. Whilst that behaviour appears to be exhibited by flow in long channels converging to a central opening (Fig. 8a), it does not demonstrate chokes and compartments within a more general flow regime. It was decided therefore to search for chokes using the faceto-face boundary conditions shown in Fig. 4 and networks of long channels close to critical density. The approach required the development of an add-on code (POOLS) that searches for clusters of nodes, between which no adjacent pairs have a head drop greater than a pre-determined value (Δh). This process is repeated a number of times using increasing values of Δh until there is a small irreducible number of clusters, sometimes a single cluster. An initial investigation of this type was performed on three realizations based on the network generation parameters used to investigate skin around tunnels. That is, a mean channel length of 10 (sub-channels/bonds) and average sub-channel densities (PON) of 0.02, 0.04 and 0.08 where Pc is just less than 0.02. The densest network did not yield any specific chokes and the number of clusters decreased gradually with increasing Δh. The two networks closer to the percolation threshold displayed behaviour associated with chokes, and heads and flows within the denser one (i.e. PON = 0.04) are shown in Fig. 14. The system contains five clusters of wideranging sizes, namely 5,301, 418, 289, 71 and 21 bonds. One of the clusters (289) is entirely separate from the rest of the network and contains flowpath 1 linking top and base directly. The largest cluster virtually coincides with the high head channels (coloured red) located in the top half of Fig. 14a and connects via several chokes to the remaining three clusters. The contrast in density between Fig. 14a, b shows how few channels provide the dominant flow pathways. Given that Fig. 14 is of a cubic system of 90 m edge, it shows how unlikely it is for boreholes in real-life tracer tests to intercept a flowing groundwater channel network. In summary, this limited investigation suggests that only networks close to the percolation threshold are likely to Fig. 14 Heads and flows in a single realization of a HyperConv network where L = 10, G = 240 bonds and PON = 0.04: a head variation [note that the chokes and flowpath 1 have been visually enhanced for diagram clarity]; b Location of flowing channels (thickness of line is proportional to the logarithm of flow rate) exhibit choke-and-compartment behavior. Hydrogeological compartments (of head) have only recently begun to be identified in fractured crystalline rocks since it requires timeconsuming groundwater-based investigations (Sawada et al. 2000). Notwithstanding this difficulty, two projects around the world seeking sites for deep repositories for nuclear waste have begun to describe fractured crystalline rocks in terms of compartments (e.g. Kamaishi research mine in Japan (Sawada et al. 2000) and Bátaapáti possible repository site in Hungary (Benedek and Dankó 2009). Other repository investigations have recorded abrupt head changes in deep borehole profiles but possibly only the UK Sellafield investigation had sufficient density to warrant a compartment-based interpretation (Black and Barker 2016). The investigation at Kamaishi research mine in granitic rock in Japan was probably the best-reported example of head compartments (Sawada et al. 2000). Seven sub-horizontal boreholes (totalling 560 m) about 350 m below surface intersected more than 3,000 fractures of which about 600 were considered open by borehole TV logging. Some of the fractures were extrapolated between boreholes based on their location, similarity of orientation and other properties such as alteration and mineralization (Fig. 15a). Hydrogeological testing for a tracer test revealed a complex arrangement of interborehole connections that did not seem to reflect the preceding structural analysis; however, cross-hole testing did agree with pre-existing knowledge about the disposition of head and a system of compartments was hypothesized (Fig. 15b). Sawada et al. (2000) conclude, Bthe cause of the compartmentalization of flow at this site is not well understood.^ They also concluded that, Bconductive fractures are heterogeneous internally, and a conductive fracture may not produce flow at all points of borehole penetration.^ The key aspect of the compartment definition at Kamaishi is the form of the boundary of the high head compartment with its complicated intrusion into its neighbours. It is obviously impossible to employ essentially planar features to define compartment boundaries of such complexity; however, it is exactly how one would expect to observe the complex interactions between the clusters of channels of Fig. 14a). Fig. 15 Plan views of the hydrogeological and tracer testing area at Kamaishi Mine: a boreholes, open fractures and extrapolated fractures; b head compartments as derived from cross-hole testing and environmental head measurements. The compartments are coloured Physical appearance of long-channel networks There were three phases of tracer testing performed at Stripa Mine: small-scale double fracture testing (Abelin et al. 1985), small-scale single fracture testing (Abelin et al. 1994) and the medium-scale 3D experiment described in the preceding (Abelin et al. 1991a, b). Channelling was seen to underly the Fig. 16 Schematic of how the sparse channel network appears at Stripa: a artist’s impression by the author, b a ‘perspectivized’ version of Figure 19 of Abelin et al. (1994) with dimensions added; c a section of borehole wall, a ‘perspectivized’ version of Fig. 9 of Abelin et al. (1994) from dark equals high head to lightest equals low head. Numbers in compartments are values of head in metres above the access drift floor (figure adapted from Figs. 3, 7 and 8 of Sawada et al. 2000) results from the outset and a detailed picture emerged from the hole-to-hole measurements in the single fracture experiment (reproduced here as Fig. 16b, c). The intertwined braided channels in the artist’s impression in Fig. 16b are important to the overall long channel concept because they help to explain certain phenomena and are similar to the sample-based simulations of Hirano et al. (2010). The arrangement satisfies the conclusions about dead-end channels containing water that is readily accessible to tracer in actively flowing channelized groundwater (Neretnieks 2006). It also accords with the observation of ‘channel swapping’ in the Site Characterization and Validation Experiment at Stripa (Harding and Black 1992), since it places ‘semi-separate’ channels of likely similar flow resistance in close proximity. This means that if a channel becomes clogged with either precipitate (Milodowski et al. 1998) or exsolved gas (Long et al. 1992), then there is a nearby alternative flowpath readily available. Over long periods of geological time, precipitates of very different geological age could end up in close proximity within the same fracture such as was found in the Borrowdale Volcanic Group rocks beneath West Cumbria, UK (Milodowski et al. 1998). It should be noted that channelized flow has been proposed many times before, though mainly at the scale within a single fracture (Neretnieks et al. 1982; Hakami and Larsson 1996; Brown et al. 1998 and Vilarrasa et al. 2011). Some authors have proposed multi-fracture channel models, though most, like Cacas et al. (1990), envisage channel bifurcation at every fracture intersection. Only the model of Tsang et al. (1991) explicitly conceptualizes channels crossing fracture intersections without bifurcating. The concept visualized in Fig. 16 only differs from that of Tsang et al. (1991) in the extensiveness of the non-bifurcating channels such that the distance between major channel junctions—i.e. ignoring channel choices within clusters (Fig. 16b)—is shown to be 10–20 m. A review of a decades-old tracer experiment in Stripa URL coupled with the widely noted patchiness of inflows to tunnels prompted the concept of sparse networks of long channels— that is, networks of channels that do not bifurcate at every fracture-to-fracture intersection but form long sinuous channels with only occasional junctions. The concept has been investigated using an innovative 3D two-parameter bond lattice network model named HyperConv where the mean length of ‘strings’ of connected bonds (channels) and the gaps between them can be controlled by separate probability functions. It is found that networks of long channels are able to percolate at lower values of (bond) density than networks of short channels. A general relationship between mean channel length, mean gap length and probability of percolation has been developed which incorporates the well-established result for classical lattice network models as a special case. Positive skin effects were common in the underground experiments at Stripa and were interpreted as regions of reduced hydraulic conductivity. HyperConv simulations were performed in this study using cylindrically convergent boundary conditions and yielded apparent positive skin under certain circumstances. More specifically, only realizations of channels longer than the diameter of the inner boundary (e.g. the tunnel) exhibited positive skin and reached maximum values just above critical density. Based on these numerical model results and calculations of a single channel configuration, a semi-generic pattern of behaviour is postulated that indicates that networks of short channels (relative to the size of the underground opening) are unlikely to yield positive values of apparent skin. A few ensembles based on variable conductance bonds (all other realizations used constant conductance) had only a minor effect on this unexpected behaviour. A relatively simple analytical relationship in terms of overall channel density and the frequency of channels intersecting the inner boundary is developed that largely replicates the behaviour seen in the numerical experiments. Several researchers have suggested in the past that flow in sparse networks will be dominated by a few key connecting channels (or fractures) termed ‘chokes’. A few numerical experiments were undertaken where it was found that chokes only become significant close to the percolation threshold. Values of channel density greater than about 4× critical density seem to rule out the creation of recognisable chokes. The occurrence of chokes is linked to the creation of head ‘compartments’ which, although the HyperConv simulations represent conditions at the Stripa URL, appear very similar to actual head compartments reported at Kamaishi experimental mine in Japan. It is concluded that, based on this relatively limited numerical study, the long-channel concept is capable of reproducing many groundwater-related effects observed in underground research laboratories in fractured crystalline rocks. Networks of long channels are intimately linked to sparseness which is itself linked to ‘apparent skin’, chokes and head compartments. Where these effects are observed would seem to indicate that the groundwater is flowing in a long-channel network close to critical density. These links between sparse longchannel networks and observable features in the field present the possibility of recognising and characterising such networks in future site characterization projects for nuclear repositories in fractured crystalline rock. Acknowledgements The model ‘HyperConv’ was developed by Peter Robinson of Quintessa Ltd. This study was supported in part by the Swedish Nuclear Fuel and Waste Management Company, SKB. Appendix: the interpretation of (hydraulic) skin around tunnels based on head measurements ‘Skin’ is a concept used to describe poor performance of groundwater and oil production boreholes relative to what would be expected under ‘ideal’ conditions. It is a catchall term for effects due to permeability variations close to the borehole wall. Within nuclear repository projects, the concept has been adapted to apply to tunnels and drifts (Rhén et al. 1997), even though tunnels are horizontal rather than vertical and timescales, length-scales and rocks are significantly different. In groundwater hydraulics, skin is usually quantified in terms of a dimensionless ‘skin factor’, sf, additional to any well function describing head loss involved in flow through an aquifer towards a sink. Normally, this factor would be added to the well function evaluated at the well, but to adopt the convention of the tunnel being at zero head, the skin factor is added to the head in the formation outside the well: where the mathematical terms are as follows: Using the well function representing cylindrically symmetrical flow in an infinite system, Eq. (9) can be re-arranged to give: where an effective radius is introduced given by: reff ¼ rte−s f Thus, a positive skin factor (equating to a skin of reduced permeability) corresponds to an effective radius of the tunnel which is less than the actual radius, and vice versa (see Fig. 17e, f). A somewhat more real-world solution than Eq. (9) is obtained by assuming a horizontal water table (at a constant head) combined with constant pressure conditions at the tunnel wall (like an air-filled tunnel) to yield a solution of more general application. Equations for the head distribution were derived by El Tani (2003). The full solution is too complex to reproduce here but a key result is the flux to the tunnel: ς ¼ rzt − t where zt = depth of the tunnel centre below the water table. When zt/rt is large, the El Tani (2003) solution differs from the Thiem equation (with zero head at radius 2rt) by less than 1% when zt/rt > 7.2. As can be seen in Fig. 17a, the system is symmetrical about the vertical plane intersecting the axis of the tunnel but not about the horizontal plane, since all the inflow is sourced from the water table above. As the ratio of depth below the water table relative to tunnel diameter increases, the head distribution becomes more radially symmetrical about the tunnel axis. Cylindrically symmetrical flow at great depth is the assumption underlying the well-known approximate solution of Goodman et al. (1965). The asymmetry of the solution is seen as divergent heads dependent on direction from the tunnel centreline (Fig. 17c). The maximum amount of divergence occurs between vertically upward and vertically downward directions and can be generalized by expressing depth and distance in terms of tunnel diameters (Fig. 17d). For example, a 4-m-diameter tunnel 400 m below the water table (i.e. 100 tunnel diameters) would have upward and downward measurements of head that differed by 4 m (i.e. 1 tunnel diameter) at 6.1 tunnel diameters distance from the drift centreline (i.e. 24 m; from Fig. 17d). It can be concluded that for either a water-filled tunnel or a tunnel at atmospheric pressure, groundwater flow can be approximated as cylindrically convergent flow to a line source when at distances of more than a few times the tunnel diameter. This can be represented by the Thiem equation with skin: r ln r þ s f t Hence, a plot of head versus the log of distance from the tunnel centreline should yield a straight line passing through zero at the tunnel wall (i.e. Fig. 17c). The deviation of real measurements from this ideal gives rise to interpretations of skin for tunnels. Most commonly, head variation around tunnels is measured in boreholes drilled ‘out’ from the tunnel itself so that measurements are readily expressed in terms of distance from the tunnel (Fig. 17b). In nuclear repository projects, boreholes are instrumented with packers to hydraulically isolate individual sections so that resultant head measurements apply to specific borehole intervals (e.g. as in Fig. 17b, e, f). All of the formulae given earlier are based on the assumption of flow in a homogeneous isotropic porous medium. In real fractured crystalline rock, measurements usually vary much more widely than that caused by the asymmetry of the theoretical solution. This variation is usually interpreted as resulting from the basic heterogeneity of fractured rock and the regression line is assumed to be indicative of the average behaviour of the groundwater system interacting with the underground-opening. In many cases, the regression line does not intercept the tunnel wall at zero head (Fig. 17e, f). Open Access This article is distributed under the terms of the Creative C o m m o n s A t t r i b u t i o n 4 . 0 I n t e r n a t i o n a l L i c e n s e ( h t t p : / /, which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. Abelin H , Neretnieks I , Tunbrant S , Moreno L ( 1985 ) Final report of the migration in a single fracture: experimental results and evaluation . OECD/NEA, Stripa Project TR 85-03, Swedish Nuclear Fuel and Waste Management Co ., Stockholm, 103 , pp Abelin H , Birgersson L , Gidlund J , Moreno L , Neretnieks I ( 1991a ) A large-scale flow and tracer experiment in granite: 1. experimental design and flow distribution . Water Resour Res 27 : 3107 - 3117 . doi:10.1029/91WR01405 Abelin H , Birgersson L , Moreno L , Widén H , Ågren T , Neretnieks I ( 1991b ) A large-scale flow and tracer experiment in granite: 2. results and interpretation . Water Resour Res 27 : 3119 - 3135 . doi:10.1029/91WR01404 Abelin H , Birgersson L , Widén H , Ågren T , Moreno L , Neretnieks I ( 1994 ) Channeling experiments in crystalline fractured rocks . J Contaminant Hydrol 15 : 129 - 158 . doi:10.1016/ 0169 -7722(94 )90022- 1 Andersson P , Byegård J , Winberg A ( 2002 ) TRUE Block Scale Project Final Report 2: tracer tests in the block scale . Technical report TR02-14, Swedish Nuclear Fuel and Waste Management Co ., Stockholm Andersson P , Byegård J , Billaux D , Cvetkovic V , Dershowitz W , Doe T , Hermanson J , Poteri A , Tullborg E-L , Winberg A ( 2007 ) TRUE Block Scale Continuation Project final report . TR-02-14, Technical report, Swedish Nuclear Fuel and Waste Management Co ., Stockholm Barker JA ( 1988 ) A generalised radial flow model for hydraulic tests in fractured rock . Water Resour Res 24 : 1796 - 1804 . doi:10.1029 /WR024i010p01796 Benedek K , Dankó G ( 2009 ) Stochastic hydrogeological modelling of fractured rocks: a generic case study in the Mórágy Granite Formation (South Hungary) . Geol Carpath 60 /4: 271 - 281 . doi:10.2478/v10096- 009 - 0019 -y Black JH , Barker JA ( 2016 ) The puzzle of high heads beneath the West Cumbrian coast , UK: a possible solution . Hydrogeol J . 24 : 439 - 457 . doi: 10.1007/s10040- 015 - 1340 -4 Black JH , Robinson PC , Barker JA ( 2006 ) A preliminary investigation of the concept of hyper-convergence . Research report R-06-30, Swedish Nuclear Fuel and Waste Management Co ., Stockholm, 66 pp Black JH , Barker JA , Woodman N ( 2007 ) An investigation of 'sparse channel networks': characteristic behaviours and their causes . Research report R-07-35, Swedish Nuclear Fuel and Waste Management Co ., Stockholm, 123 pp Brown S , Caprihan A , Hardy R ( 1998 ) Experimental observation of fluid flow channels in a single fracture . J Geophys Res 103 (B3): 5125 - 5132 . doi:10.1029/97JB0354 Cacas MC , Ledoux E , de Marsily G , Tillie B , Barbreau A , Durand E , Feuga B , Peaudecerf P ( 1990 ) Modeling fracture flow with a stochastic discrete fracture network: calibration and validation-the flow model . Water Resour Res 26 /3: 479 - 489 . doi:10.1029 /WR026i003p00479 Chapman NA , McKinley IG , Penna Franca E , Shea MJ , Smellie JAT ( 1992 ) The Pocos de Caldas project: an introduction and summary of its implications for radioactive waste disposal . J Geochem Explor 45 : 1 - 24 . doi:10.1016/ 0375 -6742(92) 90120 - W Chijimatsu M , Fujita T , Sugita Y , Ishikawa H ( 1996 ) Coupled thermohydro-mechanical experiment at Kamaishi Mine: initial data around T-H-M experiment area . PNC report 03-95-03 , JAEA, Ibaraki, Japan Einsiedl F , Maloszewski P ( 2005 ) Tracer tests in fractured rocks with a new fluorescent dye: pyrene- 1 , 3 , 6 , 8 - tetra sulphonic acid (PTS). Hydrol Sci J 50 ( 3 ): 543 - 554 . doi:10.1623/hysj.50.3.543.65026 El Tani M ( 2003 ) Circular tunnel in a semi-infinite aquifer . Tunn Undergr Space Technol 18 : 49 - 55 . doi:10.1016/S0886-7798(02)00102- 5 Elert M , Svensson H ( 2001 ) Evaluation of modelling of the TRUE-1 radially converging tests with conservative tracers . The spö Task Force on Modelling of Groundwater Flow and Transport of Solutes, Tasks 4E and 4F . Technical report TR-01-12, Swedish Nuclear Fuel and Waste Management Co ., Stockholm Follin S , Hartley L ( 2014 ) Approaches to confirmatory testing of a groundwater flow model for sparsely fractured crystalline rock, exemplified by data from the proposed high-level nuclear waste repository site at Forsmark , Sweden. Hydrogeol J 22 : 333 - 349 . doi:10.1007/s10040- 013 - 1079 -8 Goodman RE , Moye DG , van Schalkwyk A , Javandel I ( 1965 ) Ground water inflows during tunnel driving . Eng Geol 2 ( 1 ): 39 - 56 Hakami E , Larsson E ( 1996 ) Aperture measurements and flow experiments on a single natural fracture . Int J Rock Mech Min Sci Geomech Abstr 33 ( 4 ): 395 - 404 . doi:10.1016/ 0148 -9062(95 )00070- 4 Harding W , Black J ( 1992 ) Site characterization and validation: inflow to the validation drift . OECD/NEA technical report 92-18 of the Stripa Project, Swedish Nuclear Fuel and Waste Management Co ., Stockholm, 84 pp Harter T ( 2005 ) Finite-size scaling analysis of percolation in threedimensional correlated binary Markov chain random fields . Phys Rev E 72:026120 . doi:10.1103/PhysRevE.72.026120 Hirano N , Ishibashi T , Watanabe N , Okamoto A , Tsuchiya N ( 2010 ) New concept discrete fracture network model simulator, GeoFlow, and three dimensional channeling flow in fracture network . Proceedings World Geothermal Congress 2010 Bali, Indonesia, 25 - 29 April 2010 Hunt A ( 2005 ) Percolation theory and the future of hydrogeology . Hydrogeol J 13 : 202 - 205 . doi:10.1007/s10040- 004 - 0405 -6 Long JCS , Olsson O , Martel S , Black JH ( 1992 ) Effects of excavation on water inflow to a drift . In: Myer LR, Tsang C -F, Cook NGW , Goodman RE (eds) Fractured and jointed rock masses . Proc. of Conf. on Fractured and Jointed Rock Masses , Lake Tahoe , California. Balkema, Rotterdam, The Netherlands, pp 543 - 549 Margolin G , Berkowitz B , Scher H ( 1998 ) Structure, flow, and generalized conductivity scaling . Water Resour Res 34 : 2103 - 2121 . doi:10.1029/98WR01648 Masset O , Loew S ( 2010 ) Hydraulic conductivity distribution in crystalline rocks, derived from inflows to tunnels and galleries in the Central Alps , Switzerland. Hydrogeol J 18 ( 4 ): 863 - 892 . doi:10.1007/s10040- 009 - 10569 -1 Milodowski AE , Gillespie MR , Naden J , Fortey NJ , Shepherd TJ , Pearce JM , Metcalfe R ( 1998 ) The petrology and paragenesis of fracture mineralization in the Sellafield area, west Cumbria . Proc Yorkshire Geol Soc 52 ( 2 ): 215 - 242 Moreno L , Neretnieks I ( 1989 ) Channeling in fractured zones and its potential impact on the transport of radionuclides . Symposium on Scientific Basis for Nuclear Waste Management XII , Berlin, 10 - 13 October 1988, pp 779 - 786 Moreno L , Neretnieks I ( 1993 ) Fluid flow and solute transport in a network of channels . J Contaminant Hydrol 14 : 163 - 192 . doi:10.1016 / 0169 -7722(93) 90023 - L Neretnieks I ( 2006 ) Channeling with diffusion into stagnant water and into a matrix in series . Water Resour Res 42 , W11418, 15 pp. doi:10.1029/2005WR004448 Neretnieks I , Eriksen T , Tähtinen P ( 1982 ) Tracer movement in a single fissure in granitic rock: some experimental results and their interpretation . Water Resour Res 18 : 849 - 863 . doi:10.1029/WR018i004 p00849 Pusch R , Nilsson J , Ramqvist G ( 1985 ) Final report of the Buffer Mass Test, vol I: scope preparative fieldwork and test arrangement . Tech. Report no . 85- 11 , Report of the Stripa Project, Swedish Nuclear Fuel and Waste Management Co ., Stockholm Rhén I , Gustafson G , Wikberg P ( 1997 ) Äspö HRL : geoscientific evaluation 1997/4-results from pre-investigations and detailed site characterization-comparison of predictions and observations-hydrogeology, groundwater chemistry and transport of solutes . Technical report no . 97- 05 , Swedish Nuclear Fuel and Waste management Co ., Stockholm, 334 pp Robinson P ( 1984 ) Connectivity, flow and transport in network models of fractured media . PhD Thesis , Oxford Rouleau A , Gale JE ( 1985 ) Statistical characterization of the fracture system in the Stripa Granite , Sweden. Int J Rock Min Sci Geomech Abstr 22 /6: 353 - 367 . doi:10.1016/ 0148 -9062(85)90001- 4 Rubin Y ( 2003 ) Applied stochastic hydrogeology . Oxford University Press, Oxford, 390 pp Sawada A , Uchida M , Shimo M , Yamamoto H , Takahara H , Doe TW ( 2000 ) Non-sorbing tracer migration experiments in fractured rock at the Kamaishi Mine, northeast Japan . Eng Geol 56 ( 1-2 ): 75 - 96 . doi:10.1016/S0013-7952(99)00135- 0 Selroos J-O , Walker DD , Strom A , Gylling B , Follin S ( 2002 ) Comparison of alternative modelling approaches for groundwater flow in fractured rock . J Hydrol 257 : 174 - 188 . doi:10.1016/S0022- 1694(01)00551- 0 Tsang C-F , Neretnieks I ( 1998 ) Flow channeling in heterogeneous fractured rocks . Rev Geophys 36 ( 2 ): 275 - 298 . doi:10.1029/97RG03319 Tsang C-F , Tsang YW , Hale FV ( 1991 ) Tracer transport in fractures: analysis of field data based on a variable aperture channel model . Water Resour Res 27 ( 12 ): 3095 - 3106 . doi:10.1029/91WR02270 Tsang C-F , Neretnieks I , Tsang Y ( 2015 ) Hydrologic issues associated with nuclear waste repositories . Water Resour Res 51 : 6923 - 6972 . doi:10.1002/2015WR017641 Vilarrasa V , Koyama T , Neretnieks I , Jing L ( 2011 ) Shear induced flow channels in a single rock fracture and their effect on solute transport . Transp Porous Media 87 : 503 - 523 . doi:10.1007/s11242- 010 - 9698 -1 Wilson CR , Long JCS , Galbraith RM , Karasaki K , Endo HK , DuBois AO , McPherson MJ , Ramqvist G ( 1981 ) Geohydrological data from the Macropermeability Experiment at Stripa , Sweden. Report no . LBL-12520, Lawrence Berkeley Natl. Lab., Berkeley, CA, 194 pp Witherspoon PA , Bodvarsson GS ( 2006 ) Geological challenges in radioactive waste isolation, fourth worldwide review , Rep. LBNL-59808 , Lawrence Berkeley Natl. Lab., Berkeley, CA

This is a preview of a remote PDF:

Groundwater flow into underground openings in fractured crystalline rocks: an interpretation based on long channels, Hydrogeology Journal, 2017, DOI: 10.1007/s10040-016-1511-y