Comparing the prospectivity of hydrogeological settings for deep radioactive waste disposal

Nuclear power has the potential to provide significant amounts of reliable electricity generation without carbon dioxide emissions. Disposing of radioactive waste is, however, an ongoing challenge, and if it is to be buried, the characterisation of the regional groundwater system is vital to protect the anthroposphere. This aspect is understudied in comparison to the engineered facility; yet, selecting a suitable groundwater setting can ensure radionuclide isolation hundreds of thousands of years beyond that provided by the engineered structure. This paper presents a multi-faceted scoping tool to quantitatively assess, and directly compare, the regional hydrogeological prospectivity of different groundwater settings for disposal at an early stage of the site selection process. The scoping tool is demonstrated using geological data from three distinct UK groundwater settings as a case study. Results indicate a significant difference in the performance potential of different regional groundwater settings to ensure long-term waste containment.


Introduction
Nuclear power provides a reliable and low carbon source of electrical energy, currently accounting for 21% of the total UK electricity production (BEIS 2017). An unfortunate byproduct of this process is the generation of highly radioactive and chemical-toxic waste (NDA 2017), which must be managed safely. By 2125 the UK is predicted to have 4.77 million m 3 of packaged radioactive waste requiring management (NDA 2017). The chosen international approach is to dispose of the most radioactive of this waste (approx. 650,000 m 3 in the UK (NDA 2014)) within a deep geological disposal facility (GDF) situated 200-1,000 m below the ground surface (Apted and Ahn 2017;NDA 2013NDA , 2014Streffer et al. 2011). Sweden and Finland are the only countries to have officially selected disposal locations for their higher activity civil waste legacy. The Waste Isolation Pilot Plant in the USA is accepting military-derived radioactive waste. All other waste producing nations (including the UK) remain in the search phase.
Geological disposal facilities must be designed as a series of complimentary yet independent physical and chemical barriers (IAEA 2009(IAEA , 2011a. These barriers are intended to contain and isolate the waste within the subsurface for many hundreds of thousands to millions of years (Environment Agency 2009;IAEA 2011a;RWM 2016a). Typical repository safety assessments are undertaken over a 1 million-year timeframe (Metcalfe et al. 2008), with the risk of death regulated at less than one-in-a-million per annum (Environment Agency 2009).
The engineered barriers will lose their integrity over timescales of tens to hundreds of thousands of years (RWM 2016a) due to chemical and mechanical degradation (Corkhill et al. 2013;Gin 2014;King et al. 2016;Metcalfe et al. 2008;RWM 2016a;Sharland et al. 2008), after which the natural barrier will become the main operational barrier (RWM 2016a). Although independent natural barrier assessment is vital to ensure overall repository performance (RWM 2016b), very little research has followed this approach. This is likely due to the geospatial and temporal complexity of subsurface characterisation.
The natural barrier comprises the rock mass and groundwater surrounding the engineered facility. Groundwater movement is an important component of a disposal facility as it is a primary transport mechanism through which released waste will be returned to the surface (RWM 2016b). It also controls the degradation rate of the engineered system (Corkhill et al. 2013;Gin 2014;NDA 2010a;RWM 2016a).
The protection of UK groundwater resources is driven by the European Union's Groundwater Directive (2006/118/ EEC) which states that the release of hazardous substances to groundwater should be prevented, and nonhazardous substances minimised. The approach to groundwater protection from developments and activities, e.g. infrastructure, landfills and liquid effluent discharges, is set out in domestic legislation, e.g. Environment Agency (2018).
The risk to groundwater is assessed through structured risk assessments. Within the UK, developments or activities deemed to pose a risk to groundwater are generally disallowed within 50 days (source protection zone 1; SPZ1) or 400 days (SPZ2) travel time to a potentially significant receptor, e.g. a potable well. Furthermore, developments or activities are tightly regulated within the same groundwater catchment (SPZ3), including when separated by a protective geological cover (Environment Agency 2018). Long groundwater pathway lengths (from source to receptor), slow groundwater velocities, and low-permeability geological barrier are all therefore important characteristics for groundwater protection.
Where disposal or storage of potential pollutants beneath the water table is concerned, e.g. some landfills, both engineered barriers and a low permeability geological barrier are required. Common engineered barriers for landfills include low permeability clay or geosynthetic liners (<1E-09 m/s, Landfill Directive, 1999/31/EC), and leachate extraction to ensure hydraulic containment.
Argillaceous (including clay), evaporitic and crystalline rocks have all been identified as potentially suitable lowpermeability GDF host rocks (NDA 2014). Research has found that within these types of formations, the groundwater flow rate can be low enough that diffusion dominates contaminant transport mechanisms (Colley and Thompson 1991;Falck and Hooker 1990;Miller et al. 1994;Wollenberg and Flexser 1984). Under these conditions, and in combination with strong sorption properties-argillaceous rocks and to a lesser extent crystalline rocks (Berry et al. 1999; NDA 2010b)-migrating radionuclides can move extraordinarily slowly. For example, a naturally occurring uranium ore surrounded by a clay rich halo (Cigar Lake, Canada) indicates containment of uranium elements in the near field for circa. 1.3 billion years (NDA 2010c). Because of the slow flow rates and strong sorption potential, clays have been the focus of host rock investigations in both France (ANDRA 2005) and Switzerland (Nagra 2002).
Although hydraulic containment can be ensured over the short-term using groundwater extraction techniques (like for near-surface landfills), long-term containment relies more heavily on passive measures (Streffer et al. 2011). A setting must therefore be selected in which pollutants are transported naturally away from resources, not towards them. This requires consideration of the regional groundwater system.
For clarification, this research considers groundwater as a pathway, and not a receptor. To ensure protection of potentially significant resources, e.g. aquifers or mineral seams (RWM 2016b), the groundwater pathway should exhibit certain 'protective' physical and chemical characteristics. Hydrogeological characteristics (HCs) identified by Chapman et al. 1986 as being of benefit for radioactive waste containment include: HC.1, slow regional groundwater movement; HC.2, long groundwater pathways prior to surface discharge; HC.3, groundwater progressively mixing with older deeper waters; HC.4, slow local groundwater movement; and HC.5, predictable groundwater pathways despite parameter and process uncertainties (Chapman et al. 1986). Separation of near-surface from deeper groundwater systems (HC.6) has also been listed as a beneficial hydrogeological characteristic (RWM 2016b).
These characteristics are exemplified through a series of idealised regional hydrogeological regimes (HRs) which include: inland basin, modified basin limb, seaward-dipping and offshore sediments, basement rocks under sedimentary cover, hard rocks in low-relief coastal environments, and small islands (Chapman et al. 1986; Fig. 1).
The HCs are controlled by regional geographical and geological features such as mountain chains, river basins and geological formations (Tóth 1963). The natural barrier should thus be assessed with respect to the regional setting, covering areas of tens of kilometres.
Consideration of the quality of the natural barrier and regional setting has been of variable importance within international site selection processes-for example, Forsmark (Sweden) was selected over Laxmar (both situated in low relief tectonic lenses of the Ferroscandianivian shield) due to natural barrier features which promote 'better prospects for achieving long term safety', including: 'fewer water-conducing fractures' and 'limited groundwater flow through repository' (SKB 2009). In contrast, the shortlisting and selection of Sellafield (UK) for an underground rock characterisation facility (RCF) (commonly perceived as a precursor to a GDF) (Letter from the Director of Infrastructure & Planning for the North West, dated 17 March 1997;UK Parliament 1997) was 'strongly influenced by non-technical [social, political and economic] factors' such as 'restricted to sites that were owned by central government, or by its nuclear industry shareholder' (Nirex 2005a) and because 'waste transport would be less of an issue' (Nirex 2002(Nirex , 2005a. The planning application for the RCF at Sellafield was refused in 1997 due, in part, to the site being 'geologically and hydrogeologically much less simple and more complex than would be expected of a choice based principally on scientific and technical grounds' (Letter from the director of Infrastructure & Planning for the North West, dated 17th March 1997;UK Parliament 1997). Lack of focus on the natural barrier characteristics was thus a contributing factor to the rejected application for the RCF.
Although demonstration of the natural barrier performance is a prerequisite for the safety assessment (IAEA 2011a), no obligatory method exists which makes comparison of the relative benefits of different regional settings difficult. This paper provides a solution to this challenge by proposing a method, based on the set of established and transferable groundwater characteristics defined by Chapman et al. 1986, to compare the prospective performance of different settings. Three coupled process models are developed based on UK data as a case study. Variable rock permeability scenarios are explored as key performance uncertainties (Domenico and Schwartz 1997). Settings are compared to an 'idealised' benchmark scenario. Results indicate that the natural barrier, and thus overall repository performance, varies up to a factor of 4 between assessed settings.

Methods
The method is undertaken in three stages: firstly, a benchmark scenario is defined against which the HCs of the chosen settings could be directly compared; secondly, three settings are selected for assessment; and thirdly, numerical regional groundwater flow models to represent each of the selected settings are developed.

The benchmark scenario
The benchmark scenario is defined based on the previously listed beneficial HCs from Chapman et al. 1986, which is aligned with a common resource protection strategy as discussed previously. The HCs are converted into quantifiable parameters (Ps), such that P.1 quantitatively represents HC.1, etc. The Ps have been developed so as to be applicable to any modelled setting The conditions chosen for P.1-P.4 (summarised in the following) are detailed further in Table S1 of the electronic supplementary material (ESM), illustrated in Fig. S4 Fig. 1 Conceptual illustration of idealised regional hydrogeological regimes which exhibit wide ranging groundwater characteristics considered to be of benefit for long-term radioactive waste containment. Image digitised and adapted from Chapman et al. (1986) with categorization as: a inland basin, a′ modified basin limb, b seawarddipping and offshore sediments, c basement rocks under sedimentary cover, d hard rocks in low-relief coastal environments, and e small islands P.4. (representing HC.4) is the distance released radionuclides travel over the first 10,000 years.
As shown in the preceding, HC.1-HC.4 have been integrated directly into the assessment criteria. HC.5 is accounted for by running two versions of each modelled setting -a most-likely permeability and a high permeability model-in order to account for key epistemic uncertainties. The difference in overall scores, therefore, reflects the predictability of each setting as a result of permeability uncertainty. Finally, HC.6 is already accounted for through HC.1-HC.4 as a coupled near surface to deep groundwater system would permit direct return of radionuclides to the surface (short and ascending groundwater pathways) and would be subject to faster groundwater velocities through connection to topographic drive. The findings of HC.6 are therefore discussed qualitatively for each modelled site. The benchmark scenario scores 0, with points accrued (up to 20) for negative characteristics (Fig. 2).

Site selection for assessment
To illustrate this method, three exemplar settings are selected (see Fig. 3). Site selection is based on perceived geological (and therefore hydrogeological) diversity, and historic literature indications of potential hydrogeological suitability as discussed in the following. The UK is used as a case study due to the availability of geological and hydrogeological data. This method could however be applied anywhere in the world as the assessment method is based on universal hydrogeological parameters. At the time of writing, no sites have been officially selected for GDF investigation within the UK.
Setting 1 will be based on rock geometry and properties from Sellafield (West Cumbria), where £400 million of drilling investigations were undertaken up to 1997 for the purpose of the RCF (Nirex 1997a, b, c, d). Comparison of the HCs of this location to other UK settings is thus anticipated to be of interest to the wider scientific community. The setting also provides an opportunity for model verification from the datasets obtained during the RCF investigations (Nirex 1997a, b, c, d), and from previous site-specific modelling (Fraser-Harris et al. 2015).
Setting 2 is located offshore and is based on rock geometries and properties from the Tynwald Basin (East Irish Sea Basin). Setting 2 was chosen due to laterally extensive low permeability sediments (Jackson et al. 1995(Jackson et al. , 1997, and a possible offshore dense brine formation (Barnes et al. 2005), indicative of a hydrogeologically stable environment over inter-glacial timescales (Park et al. 2009). This type of setting could be considered potentially analogous to the seaward dipping and offshore sediments type HR (HRs listed in section 'Introduction'). At the time of writing, no  Fig. 2 a Idealised benchmark outflow scenario, including hydrogeological parameters 1-4, against which all model scenarios are assessed. b Proposed bar-chart assessment criteria to be applied to all model scenarios in which a score of 0 indicates beneficial hydrogeological characteristics for long-term waste containment, and a score of 20 indicates poor groundwater characteristics country has disposed of its higher activity waste legacy offshore. Setting 2 is therefore also anticipated to be of wider scientific interest.
Setting 3 beneath Thetford (East Anglia), was chosen due to a laterally extensive sedimentary sequence overlying crystalline basement (Lee et al. 2015), which could be considered analogous to the basement rock beneath sedimentary sequence type HR. This setting therefore has the potential to exhibit advantageous HCs for long-term radioactive waste containment as defined by Chapman et al. 1986.

Model development
This section details the construction of numerical models to represent each of the three chosen settings.

Geometry
Model geometries (approx. 30 km length by 2-4 km depth) were obtained from publicly available 2D British Geological Survey geological cross-sections, and converted into 2D conceptual hydrogeological cross-sections by combining geological units considered to behave in a similar hydrogeological manner. Setting 1 was constructed using a combination of geological cross-sections from Gosforth (BGS 1999a), the lake district fells (Michie 1996) and the previously investigated RCF-Nirex 1997a, b; Fig. S1 of the electronic supplementary material (ESM). Setting 2 was constructed using a geological cross-section from Bootle (BGS 1997), supplemented by offshore borehole and seismic-derived stratigraphic data (Akhurst et al. 1997;Barnes et al. 2005;Jackson et al. 1987;Schlumberger 2016) (Fig. S2 of the ESM). Setting 3 was constructed using a geological cross-section from Swaffham (BGS 1999b), supplemented with additional stratigraphic information for East Anglia (Lee et al. 2015;Fig. S3 of the ESM). Fault zones were represented as 50-m-wide continuum porous m edia in keeping with previ ous hydrogeological modelling assessments (Fraser-Harris et al. 2015;McKeown et al. 1999;Nirex 1997c).

Mesh
Two-dimensional (2D) conceptual hydrogeological cross-section geometries were converted into a mesh of 2D triangular elements using the GMSH meshing software (Geuzaine and Remacle 2009). Element sizes ranged from 5 to 100 m, with higher mesh densities applied to thinner hydrogeological units. Setting 1 comprised 173,527 elements, setting 2 comprised 258,920 elements, and setting 3 comprised 117,461 elements.

Material properties
Models were populated with the material properties of porosity, permeability (both 'most-likely' and 'high' values as a key performance uncertainty), mass dispersion, heat dispersion, storativity, bulk density, thermal capacity, thermal conductivity, specific heat capacity, and specific heat conductivity (Tables S2-S4 of the ESM).
Material properties were applied to each individual hydrogeological unit using site-specific information where available, else generic data ranges were applied. When a lithologically layered hydrogeological unit was present, e.g. for setting 3, the effective horizontal (K x ) and vertical (K z ) permeabilities (m s −1 ) were calculated using Eqs. (1) and (2) respectively (Lee and Fetter 1994) where k is the hydraulic conductivity (m s −1 ), n is the number of layers within the aquifer, i is the layer, and b is the layer thickness (m).
Hydraulic conductivities (m s −1 ) were converted into intrinsic permeability (m 2 ) using Eq. (3) (Domenico and Schwartz 1997). K itself comprises the fluid properties of density ρ (kg m −3 ), acceleration due to gravity g (m s −2 ), and the dynamic viscosity μ (Pa s), and the geometric component of intrinsic permeability k (m 2 ) of the host rock formation (McDermott et al. 2006).
Mass and heat dispersion, which represent characteristic properties of the geological medium (Domenico and Schwartz 1997), but can also be used as a tool for numerical stability, were set as either 10 or 50 m, depending on which most closely represented half the average element length for that unit. This ensured numerical (Peclet) stability (Anderson and Woessner 1992) (Eq. 4). Recommendations are for dispersion to be ascertained at a field scale.
where Pe is the Peclet number (−), v a is the advective velocity of groundwater (m s -1 ), x is the grid size (m), D x is hydrodynamic dispersion in the x-direction (m 2 s −1 ), and α is dispersivity (m).
Specific heat capacity and conductivity were applied as constants across the entire modelled domain, since dynamic viscosity is a more influential fluid parameter when simulating coupled process fluid flow (Watanabe et al. 2010). Fluid density and viscosity were therefore represented as concentration and temperature dependent functions across the model domain.

Boundary conditions
Boundary conditions were applied for liquid flow, heat transport and mass transport (see section 'Processes') as constant pressure, temperature and mass concentration boundaries respectively (Figs. S1-S3 of the ESM). Onshore surface pressure boundaries were calculated using Eq. (5).
where P a is atmospheric pressure at sea level (101,325 Pascals, Lide 2004), and h is the elevation above sea level (m). Offshore seabed pressure boundaries were calculated using Eq. (6).
where P a and h are as already described for Eq. (5), ρ seawater is the density of seawater (1,025 kg m −3 ; Grasshoff et al. 1999), and g is acceleration due to gravity (9.81 m s −2 ). Subsurface pressure boundaries were calculated using Eq. (7).
where P onshore surface/offshore seabed (Pascal) is either Eqs. (5) or (6) depending on the location of the sub-surface pressure boundary, ρ water is the density of the water, and g and h are as described for Eq. (6). Subsurface temperature boundaries were calculated using Eq. (8).
where T onshore surface/ offshore seabed is the temperature (°C) at either the onshore surface (8.5°C for settings 1 and 2, and 9.78°C for setting 3) or offshore seabed (average 8.5°C, Joint Nature Conservation Committee 2003), ΔT is the geothermal gradient for a typical continental setting (0.025°C (Downing and Gray 1986)), and h is the distance below ground level (m).
Mass concentration (C m ) boundaries (kg m −3 ) were calculated using Eq. (9), where Cl is chloride concentration in mg L −1 (making up the majority of dissolved ions in groundwater by weight, SOEST 2015), and ρ water is the density of freshwater (1,000 kg m −3 ).

Initial conditions
Initial conditions for water pressure (liquid flow) and temperature (heat transport) were calculated for each modelled ð3Þ setting using Eqs. (7) and (8) respectively, and applied linearly to (1) the onshore sub-surface, and (2) the offshore subsurface (Tables S5-S7 of the ESM).
Initial conditions for mass concentration (mass transport) were applied to setting 1 using a mass concentration distribution across the model domain obtained from inverse distance weighted borehole data. Mass concentration was based on chloride values obtained from boreholes 2, 3, 4, 9A/B and 10 at Sellafield, West Cumbria (Bath et al. 2006). Initial mass conditions were applied to setting 2 as a linear gradient for the offshore quaternary, grading from seawater at the surface to the fully salt saturated underlying Mercia Mudstone Group (MMG). Constant mass densities were applied to the MMG assuming equilibrium, and to the underlying sedimentary units based on average offshore well data (Barnes et al. 2005;Bastin et al. 2003;Cowan and Boycott-Brown 2003;Yaliz and Chapman 2003;Yaliz and McKim 2003;Yaliz and Taylor 2003). A coastal freshwater interface was subsequently superimposed onto the northeast corner. Initial mass concentration conditions were applied to setting 3 as a linear gradient to (1) the sedimentary cover, and (2) the Silurian basement, based on mapped chloride concentrations for the areas (BGS 1976).

Processes
The models were run using the open-source finite element code 'OpenGeoSys', which has been extensively benchmarked and tested OpenGeoSys 2017). OpenGeoSys couples the processes of liquid flow, heat transport, and conservative mass transport ). The three-dimensional liquid flow (saturated), heat transport and mass transport equations are presented in Eqs. (10-12) respectively (Niemi et al. 2017).
where S s is the storage coefficient (Pa −1 ), P is the fluid pressure (Pascal), t is time (s), k is the intrinsic permeability (m 2 ), μ is the dynamic viscosity (Pa s), ρ is the fluid density (kg/m 3 ), g is the acceleration due to gravity (m/s 2 ), z is the elevation head (m), and Q is the source/sink term (m 3 /s).
where D T is the heat-diffusion-dispersion tensor for porous medium (W/mK), T is the temperature (K), c w is the specific heat capacity of fluid (Jkg K), ρ w is the fluid density (kg/m 3 ), v a is the advective velocity (m/s), ρ is the density of the saturated porous rock (kg/m 3 ), Q T is the heat source or sink (J/kg K), c m is the specific heat capacity of the saturated porous rock (kg/m 3 ), ρ m is the density of the saturated porous rock (kg/m 3 ), and t is time (s). D T comprises an effective heat diffusion coefficient (W/mK) and a heat dispersion coefficient (J/K m 2 ) due to advective velocity (m/s).
where C is the solute concentration (kg/m 3 ), t is time (s), D is the hydrodynamic dispersion tensor (m 2 /s), v a is the advective velocity of the fluid (m/s), R f is the retardation factor, C s is a concentration source term (kg/m 3 s), such as input from a chemical spill, C r is a concentration source term due to chemical reactions (kg/m 3 s), and C λ is a concentrations source term due to radioactive decay (kg/m 3 s). D is itself a function of the dispersivity (m), advective velocity (v a ), and the effective molecular diffusion coefficient (m 2 /s). The process equations (Eqs. 10-12) are coupled via the fluid material properties of density ρ (kg m −3 ) (Eq. 13) and dynamic viscosity μ (Pa s) (Eq. 14).
where v a is the advective velocity (m s −1 ), n e is the effective porosity (−), K is the hydraulic conductivity (m s −1 ) (Eq. 3), and h is the hydraulic head (m).

Model run
Liquid flow (pressure) was first solved for. Once achieved, the pressure distribution was used as an input for the coupled liquid flow, heat transport and nonreactive mass transport models. Models were run to steady-state conditions using a 100-year timestep. Simulation took place on four cores of the University of Edinburgh's Linux high performance compute cluster 'Eddie3' (University of Edinburgh 2018).

Calibration
Setting 1 was calibrated to measured field data for (1) freshwater head, (2) salinity and (3) temperature. Calibration of freshwater head, salinity and temperature was closest to measured field values for the most-likely permeability model scenario, in particular within the vicinity of the proposed rock characterisation facility where the greatest confidence in model results is required. Images of the calibration for setting 1 is provided within the (Figs. S6-S8 of the ESM). No calibration was possible for settings 2 and 3 due to an absence of sitespecific pressure, salinity and temperature data. After calibration the quasi-steady-state-coupled-process models were made transient through the addition of a storage term, obtained from publicly available information, and run for a further 10,000 years.

Implementation of the proposed method
To determine the percentage of the far-field domain with 'slow groundwater movement', a representative groundwater velocity required defining. Although this research aims to consider hydrogeological processes, rather than explicitly solute transportation processes, solute transport via diffusion is, over longer timescales, a relatively slow method of contaminant transportation. This is especially true compared to the transportation speed which can be achieved by advection (Domenico and Schwartz 1997). Therefore, if the advective groundwater velocity is less than that of the rate of solute transportation via diffusion, it is reasonable to also call that advective velocity 'slow'. Using Eq. (16) (Atkins 2001) where x is the distance travelled (m), D e is the effective diffusion coefficient (m 2 s −1 ) and t is time (s), 50% of the solutes can be expected to travel 77 m over 10,000 years (the chosen assessment period), based on a conservative effective diffusion coefficient of 9.31E-09 m 2 s −1 (Domenico and Schwartz 1997). Solutes traveling via advecting groundwater would be required to travel at a velocity of approx. 2E-10 m s −1 to achieve the same travel distance over the same timeframe. Therefore, for the purpose of this research, an advective groundwater velocity of ≤2E-10 m s −1 will be considered 'slow' although it is recommended that further research is undertaken to refine this value.
The percentage of the regional area with 'slow groundwater movement' i.e. an advective velocity ≤ 2E-10 m s −1 was determined by laying a 20 × 20-m grid over the model domain, and extracting the advective groundwater velocity at the grid nodes over a 2 × 20 km area around the repository. The purpose of using a grid, rather than extracting groundwater velocities from the mesh (see section 'Mesh') was to avoid bias of hydrogeological features with higher mesh densities such as faults. This data manipulation and extraction was achieved using the post-processing visualisation software 'Tecplot' (Tecplot 2018).
The total length of the quasi-steady-state groundwater pathway from repository top to discharge/exit point (P.2), and the depth of discharge relative to repository top (P.3) was determined by measurement of the groundwater pathway length within Tecplot. Parameters P.2 and P.3 are required to ensure maximum number of advantageous groundwater characteristics in the face of subsurface geometric and lithological uncertainties, and to ensure over reliance is not placed on a single safety barrier function.
The procedure to determine the distance that released particles travelled over 10,000 years (P.4) was undertaken by generating 10 evenly spaced streak-lines along the top of the repository in Tecplot, then tracking which travelled furthest within 10,000 years based on the advective groundwater velocity. See Fig. S5 of the ESM for illustration. When particles travel less than 77 m, i.e. solutes could theoretically be transported faster via diffusion than via advecting groundwater, the travel distance was corrected to 77 m for the purpose of scoring.

Method uncertainties
Two levels of uncertainty are inherent in this research: real world aleatoric and epistemic uncertainties in natural barrier performance; and uncertainties that arises from the modelling method. Uncertainty in natural barrier performance arises from aleatoric uncertainty (uncertainty about the occurrence of future events such as earthquakes, glacial events and human intrusion) and epistemic uncertainty (incomplete knowledge about the physical properties of a system such as the location or occurrence of a fault, i.e. geometric uncertainty, or permeability, i.e. lithological uncertainty; Apted and Ahn 2017). Although epistemic uncertainty can be reduced (e.g. through site investigations), aleatoric uncertainty is considered irreducible (Apted and Ahn 2017). Detailed understanding of the epistemic uncertainty of a setting is however beyond the scope of this research, which is instead to develop an assessment method that enables identification of potentially viable settings. Where epistemic uncertainty could dramatically affect the outcome of the proposed assessment method, namely permeability uncertainty, which spans several orders of magnitude (Domenico and Schwartz 1997), this has been accounted for by running two versions of each modelled setting: one populated with most-likely permeability values; and the other with high permeability values. If the modelled groundwater characteristics are not significantly changed by permeability variation, the setting is likely to perform in a more predictable manner.
Uncertainties also arise from the modelling method including: incorrect conceptualisation (i.e. key processes) of the problem; mathematical misrepresentation; inappropriate selection of a numerical method to solve the mathematical equations; and uncertainty in data used to populate the model, e.g. permeability uncertainty and unknown boundary conditions (Konikow and Bredehoeft 1992). Again, the latter is beyond the scope of this research, which is instead designed as a scoping tool (using pre-existing information) to identify settings that are worth further investigation to reduce epistemic uncertainties. Previous research identified thermal, hydraulic and nonreactive chemical processes as key processes in farfield natural barrier research (Fraser-Harris et al. 2015). The influence of mechanical processes in a steady-state system are considered less that the uncertainty generated from parameter upscaling Hudson et al. 2005), i.e. from a single well over hundreds of meters, and as such have been excluded from consideration. Finally, uncertainties that arise from the chosen numerical method are also considered small compared to epistemic geological uncertainties (Istok 1989). For settings where enough site-specific data currently exists-i.e. setting 1 (Nirex 1997a, b, c, d)-models have been calibrated to improve model confidence (Figs. S6-S8 of the ESM).

Results
Setting 1: key hydrogeological parameters score P.1. 57% of the 20 × 2 km far-field area (see section 'Methods') exhibits very low rates of advective (pressure driven, Domenico and Schwartz 1997) groundwater movement, with 'very low' being defined as <2E-10 m s −1 (see section 'Methods'). This reduces to 20% when simulated with high permeability values P.2. Both high and most-likely permeability model simulations indicate groundwater pathways <2,500 m in length. P.3. In both high and most-likely permeability models, pathways ascend (600 m) directly from the repository to the surface (Fig. 4a,b). P.4. Released radionuclides travel small distances over 10,000 years within the most-likely permeability scenario, i.e. slow enough to be via diffusion (defaulted to <77 m, see section 'Methods'). Radionuclides travel 602 m in the high permeability scenario, exiting the host rock within 10,000 years. Overall score. Setting 1 accrues 12 and 16 out of 20 negative points for the most-likely and high permeability modelled scenarios respectively (Fig. 4c).
Modelling of setting 1 shows that the pattern of regional groundwater behaviour matches previous model simulations of the area (Fraser-Harris et al. 2015;Mckeown et al. 1999;Nirex 1997c), providing confidence in the method of regional flow assessment for settings 2 and 3.
Setting 2: key hydrogeological parameters score P.1. 66% of the 20 × 2 km-far-field area exhibits very low rates of advective groundwater movement (see section 'Methods'). This reduces to 39% when simulated with high permeability values. P.2. Both high and most-likely permeability model simulations indicate long groundwater pathways of >13,700 m from repository to model boundary. P.3. In both high and most-likely permeability models, pathways descend >1.6 km from the repository, and discharge to the onshore (northeast) boundary. P.4. Both high and most-likely permeability model simulations indicate released radionuclides travel small distances over 10,000 years i.e. slow enough to be via diffusion (defaulted to <77 m, see section 'Methods'). Released radionuclides would take at least 35,000 years to leave the host rock formation. Overall score. Setting 2 accrues 3 and 4 out of 20 negative points for the most-likely and high permeability modelled scenarios respectively (Fig. 5c).
Setting 3: key hydrogeological parameters score P.1. 56% of the 20 × 2 km far-field area exhibits very low rates of advective groundwater movement (see section 'Methods'). This reduces to 47% when simulated with high permeability values. P.2. When modelled with most-likely permeabilities the setting shows long groundwater pathways (12,900 m). When modelled with high permeabilities, two pathways form. The first is similar to the most likely permeability model (13,000 m). The second however is shorter at 3100 m. P.3. When modelled with most-likely permeabilities the groundwater pathway descends (150 m) from the repository and discharge to the west-northwest boundary (Fig. 6a). When modelled with high permeabilities, the first pathway is similar to the most likely permeability model descending (183 m) from the repository and discharge to the west-northwest boundary (Fig. 6b), however, the second ascends 630 m directly to the surface. P.4. When modelled with most-likely permeabilities released radionuclides travel 375 m in 10,000 years i.e. remaining within the host rock (Fig. 6a). When modelled with high permeabilities, released radionuclides travel c Bar chart presenting beneficial groundwater characteristics scores for the most likely permeability (blue) and high permeability (yellow) modelled scenarios 5574 m (Fig. 6b) along the first pathway, remaining within the host rock formation, but 3100 m along the second pathway, discharging to the surface and exiting the model. Overall score. Setting 3 accrues 7 and 16 out of 20 negative points for the most-likely and high permeability modelled scenarios respectively (Fig. 6c). The scoring for the high permeability scenario is based on the second ascending pathway as this is deemed of greatest risk from a safety perspective.

Discussion
Setting 1: hydrogeological characteristics and performance potential The regional groundwater flow pattern of setting 1 is controlled by enhanced topographic elevation (Lake District Fells) to the east, driving less dense rainfall derived groundwater westwards. The westward groundwater progression is blocked by the offshore dense 'Irish Sea Brine' regime, formed over millions of years from the dissolution of offshore salt rich layers (Bath et al. 1996;Black and Brightman 1996), which block and then forces groundwater up through the vicinity of the repository. This pattern of regional flow causes short and undesirable groundwater pathways (HC.2), which progress directly to the surface (HC.3). The direct coupling between near-surface and deep groundwater (HC.6) is undesirable from a containment perspective and means that setting 1 cannot be considered analogous to the basement rock beneath sedimentary cover type HR (see Fig. 1). In addition, the direct coupling is likely to increase deep groundwater vulnerability to glacial flushing. This is supported by research from West Cumbria which attributes the high hydraulic heads presently observed in the Borrowdale Volcanic Group (basement rock) to be from the late Devensian glacial retreat (Black and Barker 2015), suggesting hydraulic coupling.
The host rock formation (lithological formation in which the engineered facility is constructed) is only likely to have wide spread diffusion-dominated (concentration driven, Domenico and Schwartz 1997) solute transport, which is desirable for radionuclide containment (Apted and Ahn 2017;Domenico and Schwartz 1997;Streffer et al. 2011), if field permeabilities are closest to the most-likely permeability range (HC.1 and HC.4). This indicates a lack of predictability (HC.5) in the potential performance of the setting (20% variation in overall score results), despite the extensive ground investigations previously undertaken for the RCF. Furthermore, the regional flow rates are fastest (~1.00E-07 Most-likely permeability score = 7 Ideal hydrogeological scenario score = 0 High permeability score = 16 c. Setting 3 groundwater characteristics scores a Setting 3 most-likely permeability regional groundwater system b Setting 3 high permeability Fig. 6 Regional groundwater system of setting 3 a when modelled with most-likely permeability values, and b when modelled with high permeability values. c Bar chart presenting beneficial groundwater characteristics scores for the most likely permeability (blue) and high permeability (yellow) modelled scenarios to 1.00E-08 m s −1 ) within the overlying sedimentary sequence through which radionuclides will ultimately ascend, speeding up final discharge times (HC.1). Based on the overall score, setting 1 cannot be expected to exhibit broad ranging beneficial hydrogeological characteristics deemed advantageous for long term waste containment and isolation.
Setting 2: hydrogeological characteristics and performance potential The regional groundwater flow pattern of setting 2 is primarily controlled by the groundwater density difference between the host rock formation and underlying units, and the constant pressure head from the sea, facilitating flow north-eastwards. These environmental conditions create long groundwater pathways from repository to model boundary (HC.2), which progress deeper into the earth rather than shallower (HC.3). These characteristics are advantageous from a performance perspective as they allow more time for radionuclide decay, and would also provide additional safety functions in the case of unknown geometric and lithological uncertainties, such as faults.
The high offshore chloride content (mass density) in setting 2 ( Barnes et al. 2005;Bastin et al. 2003;Cowan and Boycott-Brown 2003;Yaliz and Chapman 2003;Yaliz and McKim 2003;Yaliz and Taylor 2003) is purportedly from dissolution of halite rich Mercia Mudstone Group layers (Barnes et al. 2005;Bath et al. 1996). Groundwater residence times have been reported as 2 million years for this Irish Sea Brine formation (Bath et al. 2006), which exceeds performance assessment timescales (1 million years; NDA 2014). Although dense brine formations can occur onshore (Bein and Arad 1992;Fritz and Frape 1982), the proximity from the coastline supports formation. This is because with distance from the coast, pressure (topographic) driven flow reduces, and very long duration geological processes, such as compaction, diagenesis, cementation, and mechanical stress, begin to dominate fluid flow (Bjørlykke 1993(Bjørlykke , 1994Bjørlykke and Høeg 1997;Gluyas and Swarbrick 2003). Under these conditions, groundwater moves at velocity orders of magnitude slower than onshore (Bjørlykke 1993;Bjørlykke and Høeg 1997;Ge et al. 2003), with gases and waters diffusing extraordinarily slowly (Lu et al. 2009(Lu et al. , 2011. Dense brines have also been found to reduce the upwards vertical velocity of groundwater under flushing conditions (Johns and Resele 1997;Park et al. 2009), providing hydrogeologically stable environments over inter-glacial timeframes (Park et al. 2009).
In support of this, modelling of setting 2 suggests that regardless of permeability uncertainty, the host rock formation is likely to transport solutes via diffusion (HC.4), whilst advection could only dominate solute transport in the underlying sedimentary sequence (Fig. 5a,b; HC.1). The overall scores (only 5% variation in results) shows a degree of predictability in the potential performance (HC.5), despite the setting having never been drilled to depth and thus shows promise.
It would however be beneficial for any future research to focus on the uncertain impact of glaciation on this groundwater system, especially as some degree of near-surface to deep groundwater coupling occurs in the model simulations (HC.6). Finally, simulation of setting 2 shows flow lines to descend through the low permeability sedimentary layers, and not along them. Setting 2 cannot therefore be considered directly analogous to the previously hypothesised 'seaward dipping and offshore sediments' HR ( Fig. 1). However, based on the overall score, this location indicates broad ranging beneficial HCs for long-term waste containment and isolation, and thus warrants further investigation. Discussion of the wider implication of offshore deep geological disposal facility development is presented in section 'Comparison of performance potential and wider implications and considerations'.
Setting 3: hydrogeological characteristics and performance potential The regional groundwater flow pattern of setting 3 is sensitive to permeability uncertainty (45% variation in overall score results) (HC.5). In the most-likely permeability model, the overlying sedimentary layers behave as a low permeability seal preventing near-surface to deep groundwater coupling (HC.6). This creates a long horizontal groundwater pathway through the host rock (HC.2 and HC.3) and is thus desirable from a hydrogeological performance perspective. In the high permeability model however, the seal is ineffective and a direct coupling between near-surface and deeper groundwater occurs (HC.6), permitting radionuclide transport along a short ascending pathway to the surface (HC.2 and HC.3). This coupling is undesirable and again could create vulnerability to glacial flushing.
Regardless of permeability uncertainty, the models indicate released radionuclides to be transported via advection within the host rock (HC.4), and with diffusion-dominated transport only likely to become effective deeper within the host rock formation, or within the overlying sedimentary sequence (HC.1).
Because of the coupling, setting 3 can only be considered analogous to the basement rock beneath sedimentary cover HR (Fig. 1) if regional lithological permeabilities are closest to modelled most-likely values, but in comparison to setting 1, setting 3 does show a regional-scale resemblance to a basement rock beneath sedimentary cover regime.
The overall score reflects the uncertainty in site performance with setting 3 only exhibiting wide-ranging positive characteristics for radionuclide containment if the sedimentary sequence permeabilities are found to be closest to most-likely values. The permeability of the overlying sedimentary sequence should therefore be a key focus for any future investigation in this area.

Comparison of performance potential and wider implications and considerations
Modelling and assessment of the three exemplar groundwater settings shows the diverse quality of hydrogeological characteristics (i.e. groundwater speed and direction) available to be part of a comprehensive multi-barrier containment facility. For example, the overall scores indicate that setting 2 could exhibit hydrogeological characteristics (listed in section 'The benchmark scenario') that are 4 times more advantageous for long-term waste containment and isolation than setting 1 (see Fig. 7), i.e. when populated with most-likely permeabilities, setting 2 received an overall score of 3, which is four times more advantageous than that of setting 2 which received a score of 12. Similarly, when populated with high permeability values, setting 2 received a score of 4, which is four times better than the 16 received by setting 1.
The low score of setting 2 (close to the idealised scenario score of 0) is due to the settings' long groundwater pathways, and slow local flow rates. These hydrogeological characteristics are due, in part, to its offshore location, the characteristics of which are discussed in section 'Setting 2: hydrogeological characteristics and performance potential'. The apparent hydrogeological advantage of offshore settings for enhanced repository performance should be of interest to nation states with abundant offshore territory. Indeed, in 2016, RWM (the UK implementor) extended the search area for a GDF up to 20 km offshore the UK (RWM 2016b), which makes this research directly applicable.
Although the UK is the first country to publicly consider disposing of its high-level waste legacy offshore (both Sweden and Finland are constructing onshore GDFs (Posiva Oy 2012Oy , 2017SKB 2009)), the Swedish Final Repository (SFR) for short-lived intermediate and low-level waste was constructed under the Baltic Sea floor (SKB 2017(SKB , 2018). An extension of the SFR from 60 to 120 m below sea level is also planned (SKB 2018). The construction of the SFR, in addition to major international subseabed tunnels and mines sites including Boulby Potash Mine and the Channel Tunnel, indicate engineering feasibility. Offshore disposal could also speed up final disposal. This is because, in the UK, territorial waters are under the jurisdiction of the Crown Estate (Crown Estate 2017) who could become the sole party with whom approval for site investigations would be required. This could bypass the need for a volunteering host community which is stipulated by the NDA (2014). Technical challenges could however arise in designing appropriate ventilation systems (Parsons Brinckerhoff 2010), ensuring the retrievability of the waste if required, and public acceptability (Nirex 2005b).
Any site(s) identified for final disposal would be required to go through a rigorous safety assessment (IAEA 2009). This proposed assessment method is not intended as a safety assessment but is instead designed as a scoping tool for use at an early stage of site selection process to identify sites, with wide-ranging advantageous groundwater characteristics, determined by the proposed assessment method. By combining the scoping tool with pre-existing geological and hydrogeological data, the cost of site investigations can be minimised as settings with greatest performance potential (e.g. setting 2 which scored 3 and 4 respectively), and features within those settings most pertinent to safety (e.g. setting 3 cover rock permeability), can be focused upon.
It is understood that both technical and nontechnical factors will be involved in the site selection process; however, identification of advantageous natural settings should be at the forefront of the process because although the engineered barrier can be adapted for performance, the natural barrier cannot. Ultimately the natural barrier will be the main operational barrier preventing radionuclide return to the surface (RWM 2016a).
Timescales of natural barrier control have the potential to transcend interglacial timeframes (IAEA 2011b; greater than tens of thousands of years, Clark et al. 2012), which can exert significant forces on groundwater systems (Degnan et al. 2005;McEvoy et al. 2016;Tóth 1963;Tsang and Niemi 2013). Although this assessment method does not directly assess hydrogeological performance through glacial events,  Fig. 7 Comparison of the farfield hydrogeological characteristics of the three selected sites using the newly proposed method of assessment. Setting 2 exhibits significantly more advantageous hydrogeological characteristics for the long-term containment and isolation of radionuclides than settings 1 or 3, despite uncertainties in permeability it does identify sites that at present show promise, and thus warrant further investigation for future features, events and process (McEvoy et al. 2016) resilience. However, even with advanced site investigations, significant aleatoric and epistemic uncertainty (see section 'Method uncertainties') remains. The magnitude of uncertainty is such that a 'good enough' site would be difficult to define. The best opportunity for performance success, despite these uncertainties, lies with selection of a high-quality natural barrier setting (comprising multiple advantageous characteristics) which operate in conjunction with a complimentary engineered barrier system. It is this concept that underpins the multi-barrier safety philosophy (IAEA 2009(IAEA , 2011a, prevents overreliance on a single operational barrier function (IAEA 2011a), and ultimately drives this research.
This assessment method provides a high-level measure of the margin of hydrogeological advantage/disadvantage between settings. It is the belief of the authors that this type of quantitative measure would aid decision making by stakeholders and the public-with whom a final test of support is made (NDA 2014)-and reduce the risk of a 'poor quality', or 'marginal', site being selected.
Finally, the range of overall scores of the three exemplar settings are testament to the diversity of hydrogeological regimes present within the UK (Fig.7). This contrasts with Sweden and Finland which are both dominated by 'low-lying fracture crystalline basement' hydrogeological regimes (Posiva Oy 2012SKB 2009) of the Baltic Shield. The diversity of UK natural barrier systems should be considered a resource and an opportunity for improved long-term repository performance, and should be explored as such.

Conclusion
A method is presented to assess and score the likely performance of any regional groundwater setting, which is required as part of a comprehensive multi-barrier deep geological disposal facility. This paper demonstrates how the assessment method, which is based on individual groundwater characteristics (such as speed and direction), can be used in conjunction with publicly available geological and hydrogeological data to score settings for prospectivity at an early stage of the site selection process. The method also enables identification of hydrogeological features, within assessed settings, that are most pertinent to performance. Using three UK settings as a case study, the approach indicates a significant difference (quantified as a fourfold variation) in the performance potential of difference regional groundwater settings to ensure longterm waste containment and isolation. Further research should focus on the settings identified here as showing the greatest prospective performance potential. Highlighted is the broad-ranging advantageous hydrogeological characteristics exhibited by the exemplar offshore setting.
Funding information This research was funded through the Natural Environment Research Council's (NERC) E3 Doctoral Training Programme.

Compliance with ethical standards
Disclaimer The settings assessed within this research were selected as a theoretical exercise due to their perceived hydrogeological differences and were constructed using publicly available geological data. The settings selected here are not, to the knowledge of the authors, under official investigation by the UK implementer.
The Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/.