P1Freely available online through the publisher-supported open access option.

Abstract

The Nuclear Decommissioning Authority Radioactive Waste Management Directorate have been participating in the current DECOVALEX-2011 project (development of coupled models and their validation against experiments) one task of which has been examining the Mont Terri Ventilation Experiment (VE). This long-term (>9 years), field-scale experiment in the Opalinus Clay near the Swiss–French border, was designed to examine the coupled hydraulic–mechanical–chemical changes caused in the tunnel and in the surrounding geology, by the controlled ventilation of a 1.65 m diameter micro-tunnel.

In contrast to many conventional benchmarking and validation exercises, a key aspect of the VE as examined in DECOVALEX was that some data were held back and participants were required to make predictions of key metrics for the future evolution of the system. This paper presents an overview of the work conducted by the Quintessa and University of Edinburgh team including selected results. The coupled models developed include multiphase flow, elastic deformation and chemical processes in both detailed and upscaled geometries.

The models have been able to replicate the observed desaturation around the tunnel, tunnel deformation and localized failure, vapour migration in the tunnel, and the transition in redox conditions into the host rock.

Introduction

The performance and safety assessments of a disposal facility in a deep geological formation are supported to a large extent by the ability of the proponent to demonstrate a good fundamental understanding of the evolution of the system under construction, disposal and closure conditions. Demonstration and development of this understanding is heavily based on experimental investigations, which are conducted in the field, in dedicated underground research laboratories (URLs), and on physical models and computer codes that describe and predict the outcomes of such experiments. Among different possible host rocks, claystone is being investigated by several countries. In Switzerland, the Mont Terri URL was created in 1995 in an argillite layer (Opalinus Clay) to characterize and study the geological, hydrogeological, geochemical and geotechnical properties of this clay formation. Many instrumented experiments have been designed and installed, thus supplying a large amount of data which can be used for model validation. Among these experiments, the ventilation experiment (VE) has been carried out in a micro-tunnel to study the desaturation process which may happen because of the forced ventilation in galleries and drifts, during the construction and operation phases of a repository. In particular, this desaturation process can change the hydromechanical properties of the rock, thus influencing the performance and safety of the repository.

A specific task of the international DECOVALEX (development of codes and validation against experiments) project, in which the Nuclear Decommissioning Authority Radioactive Waste Management Directorate (NDA RWMD) is participating, relates to predicting phenomena observed during this experiment. A preliminary task considered modelling a laboratory drying test, in order to check the ability of the computer codes to reproduce the main phenomena involved in the ventilation experiment, and to identify a first set of hydraulic parameters relevant for the Opalinus Clay. Then in the following part of the task, the different phases of the ventilation experiment have been simulated: a first period (‘phase 0’) following the excavation of the micro-tunnel in 1999, up to the installation of doors that allow controlled ventilation of a portion of the micro-tunnel, and then a single resaturation–drying cycle (‘phase 1’) finishing in January 2004. This was then followed by further desaturation–resaturation cycles to the end of 2010, although data were only made available for this analysis (both descriptive and predictive) from May 2003 to April 2007. A summary of the analysis of the drying test and ventilation experiment using hydro-mechanical–chemical process sets is presented in the remainder of this paper, however available publication space does not permit an exhaustive discussion. Full details of the work undertaken under Task A including details of parameterization and formulation are discussed by Garitte et al. (2012), Bond et al. (2012a,b) and Millard et al. (2012).

Ventilation experiment and supporting experiments at Mont Terri

A short description of the ventilation experiment

The Mont Terri URL is located near a security gallery of a motorway tunnel in northern Switzerland (Bossart and Nussbaum, 2007). It is at a depth of about 400 m in Opalinus Clay, which is a stiff layered Mesozoic clay of marine origin. After the excavation of niches in 1996, a new gallery was excavated in 1998, followed by a micro-tunnel of 1.3 m in diameter in early 1999. The ventilation experiment took place in a 10 m long section of this micro-tunnel (Fig. 1).

After its excavation, the micro-tunnel was left without control of the ambient relative humidity for about 3.4 years. Then, doors were installed in order to allow control of a 10 m section, where the air inflow and the relative humidity could be imposed (Fig. 2). The micro-tunnel has then been subjected to two wetting–drying cycles. The first cycle, which is addressed in this paper, lasted from 8th July 2002 to 29th January 2004 (phase 1). After a prescribed period of 100% relative humidity (RH) of the inflowing air, a desaturation period started, corresponding to 2% relative humidity of the inflowing air. This first cycle was then followed by a second cycle, and a final resaturation which continued until 2010 (phase 2), although data are only available until April 2007. The corresponding total sequence of prescribed relative humidity is illustrated on Fig. 3 (curve RH-in, shown in red).

The micro-tunnel has been instrumented with relative-humidity sensors, pore-pressure sensors and displacement sensors along the test section. Moreover, two water pans have been installed to record the evolution of the mass loss due to the ventilation. Their location is indicated on Fig. 2. The variation of the relative humidity with time at different points along the micro-tunnel (indicated in Fig. 2), is shown in Fig. 3.

A short description of the drying test

This complementary experiment, which was used as a precursor modelling exercise to the main ventilation experiment to build confidence in codes, process models and hydraulic parameterization, has been discussed in some detail by Garitte et al. (2010) and Floria et al. (2002). The experiment also provided a useful role in understanding the possibilities for upscaling processes and parameters from the laboratory scale to the field scale. The experiment is illustrated schematically in Fig. 4. Three cylindrical samples of Opalinus Clay (101 mm radius, 254 mm in height) were placed in a controlled drying chamber along with an evaporation pan, with the axial direction oriented vertically. Chamber relative humidity and airflow was monitored continuously throughout the 142 day experiment. The samples were covered so that only the upper circular surface could lose water through evaporation. The samples and evaporation pan were also weighed continuously so that water loss could be monitored. Samples were removed and dissected at 21, 99 and 142 days such that the water content profile vertically from the evaporation surface could be monitored.

From this combination of data a continuous record of water loss for each sample as a function of the chamber conditions could be established, along with the final sample water content profile.

Numerical analysis of the drying experiment

The Quintessa multi-physics code QPAC (Quintessa, 2012) was used to simulate the drying test. The model consisted of a one-dimensional cylindrical grid discretized in the axial (vertical) direction using 12 finite volume compartments, a discretization which has been shown to be numerically sufficient for accurate representation of this case. All boundaries were set as no-flow, except for the top boundary where an evaporation boundary condition was applied. The model included full multi-phase flow of water, air and water vapour with equations for gas (g) and water (w) as follows:  
t(θρiSi)=(ρiui)+qiui=kiμi(ρi+ρigz)ki=kr,i(Si)kφ=pgpw
(1)

Here, the subscript i denotes the phase (g or w), ρ is the density (kg m−3), u is the volumetric flux (m s−1), q is an external mass source (kg s−1), S is the saturation (dimensionless), k is the intrinsic permeability (m2), kr is the relative permeability (dimensionless) and φ is the capillary pressure (Pa).

The vapour mass fluxes (kg s−1 m−2) for diffusion and advection in bulk gas are:  
uv,diff=Dv(Sg)ρvuv,adv=ugρv
(2)
where Dv is the effective diffusivity of water vapour (m2 s−1), which is assumed to be a function of bulk gas saturation. For this case the effective vapour diffusivity was assumed to be the intrinsic vapour diffusivity multiplied by the total gas saturation.

The reference model used an effective liquid pressure boundary condition (pl) based on Kelvin's Law in the chamber (equation 3) and produced results that were a good fit to the observed water contents, taking into account heterogeneity between the samples (Fig. 5).

 
pl=paRTMwρlln(RH)
(3)
where pa is the air pressure (Pa), R is the ideal gas constant (J mol−1 K−1), T is the temperature (K), Mw is the molar mass of water (kg mol−1), ρl is the density of water (kg m−3) and RH is the relative humidity in the air expressed as a dimensionless fraction.

Hydro-mechanical analysis of the ventilation experiment

Conceptual Model

The hydro-mechanical response of the system is described in detail by Garitte et al. (2010, 2012). However, in general terms the evolution of the system can be described as follows:

  1. A known flux of air with a defined relative humidity enters the sealed section of the tunnel.

  2. Interaction between the water vapour in the tunnel and the unlined tunnel host rock results in vapour exchange between the tunnel and the host rock. Evaporation of liquid water from the tunnel surface may also occur depending on local tunnel relative humidity.

  3. Water vapour leaves the tunnel via a measurement gauge which determines relative humidity and air flow rate. The difference between entries 1 and 3 above constitutes the tunnel water balance.

  4. Loss of water from the host rock to the tunnel as vapour causes a reduction in water pressure and saturation as air invades the formation from the host rock.

  5. The reduction in liquid pressure and relative humidity around the tunnel causes liquid water and water vapour (where present) to migrate towards the tunnel.

  6. Desaturation and reduction in fluid pressure causes reduction in pore volume and limited shrinkage of some of the rock skeleton, causing a local net drop in volume of the host rock.

  7. The volume change of the host rock causes localized stress changes and coupling with the hydraulic evolution through a reduction in porosity, which creates a coupling with fluid pressures and intrinsic permeability.

The dominant processes in the model are therefore: vapour diffusion in, and advection by, air in porous media and engineered volumes; viscous dominated multi-phase flow of air and water in porous media; and poro-mechanical deformation of the host rock. These are illustrated in Fig. 6.

Numerical model

Numerical modelling of the hydro-mechanical system was done in four steps. The first used a restricted dataset corresponding to phase 1 of the data to permit limited calibration and comparison against the drying test with the primary boundary condition being the observed relative humidity in the tunnel. The second step required a limited predictive hydraulic analysis whereby only the relative humidity and inflow rate of air into the tunnel was provided with no additional porous media or tunnel data and hence a prediction of behaviour focussing on tunnel water mass balance and tunnel relative humidity was required. Following a limited data release to review the mass-balance and tunnel relative humidity prediction of the second step, these hydraulic predictive models were then used for a third step, modelling reactive and non-reactive transport with comparison against geochemical measurements (section 5). For the final stage, fully blind predictions of water contents, water pressure and relative humidity were made for data collected up until 2010. Due to the complexity of this final step, the results cannot be presented here, but are discussed in Garitte et al. (2012).

The results of the first step showed that extending the model to include hydro-mechanics using a simple poro-elastic effective stress formulation (Bishop, 1959) and the hydraulic parameterization from the drying experiment matched the available measurements well. A 1D cylindrical discretization of 45 compartments using an approximately geometrical increase in compartment size from the tunnel wall was found to be sufficient given the available data. To permit solution of the second step, a simplified representation of the tunnel was created which was fully coupled with the porous medium model through the evaporation condition and diffusive vapour exchange. The tunnel model took a distribution of air velocities and assumptions regarding turbulent mixing as input to a simple free air water vapour model (equation 2). Given the low air flow velocities, and hence low Reynolds number of a maximum of 750 (Douglas et al., 2005), a parabolic air velocity profile and no enhanced radial turbulent mixing was adopted, consistent with the dominantly laminar flow characteristics for air. Predictions of tunnel relative humidities and water mass balances were very good over the whole model (Figs 7 and 8) as revealed by the data release between and the second and third steps.

The results suggest that the hydro-mechanical model represented the available data well, and moreover included a representation of the tunnel which allows the true boundaries of the experiment to be modelled directly.

Non-reactive and reactive transport analysis of the ventilation experiment

The details of the analysis are discussed in Bond et al. (2012b), and only a summary is provided in the following sections.

Data and conceptual model

Fernandez et al. (2007a,b) discuss the analysis and interpretation of transport data in full detail. However the following key points can be established:

  1. The Opalinus Clay water is marine in origin and has a background salt concentration of approximately 0.35 mol l−1.

  2. From supporting work cited in Fernandez et al. (2007b), chloride is expected to be present in between 0.62 and 0.55 of the free pore space. This is termed the ‘chloride porosity’ and was assumed to occur primarily from anion exclusion processes.

  3. There is a considerable increase in concentration, and in the absolute amount of chloride salts, close to the tunnel wall. Although the boreholes are constructed at different locations in the tunnel, at each time, the profile of chloride is quite consistent.

    Fig. 8.

    Comparison between the calculated and observed outflow water balance in the ventilation experiment tunnel.

    Fig. 8.

    Comparison between the calculated and observed outflow water balance in the ventilation experiment tunnel.

  4. There is a significant increase in sulfate close to the tunnel wall. However, because the ratio of sulfate to chloride is constant in water saturated samples at approximately 0.05 (the value for seawater), it is clear that the increase in sulfate is in excess of the relative increase seen in the chloride concentration.

It should also be noted that the geochemical analysis was revised significantly by Fernandez (2007a,b) from earlier analysis, adding to and significantly changing the experimental interpretation and stated measurements especially for the observed sulfate concentrations.

Non-reactive (chloride) analysis

The mathematical model comprises two processes; advection and diffusion. The equations governing the basic model are:  
gdiff,i=ADEff,iCigadv,i=qliquidwaterCiDEff,i=Dw,iSwθτ
(4)
where gdiff,i is the diffusive flux of species i (mol/s), gadv,i is the advective flux of species i (mol s−1), Ci is the concentration of species i in water (mol m−3), DEff,iis the effective diffusion coefficient (m2 s−1), A is the total area of porous medium transfer (m2), qliquidwater is the volumetric flux of liquid water (m3 s−1), θ is the available porosity (including any geochemical porosity effects), τ is the tortuosity (dimensionless), Sw is the water saturation (dimensionless) and Dw,i is the free water diffusivity of species i (m2 s−1).

The mathematical equations above were directly implemented in QPAC taking into account the definition of the chloride porosity (Fernandez et al., 2007a,b), and hence used the liquid water flow rate and water saturation from the hydro-mechanical analysis directly. Illustrative results are shown in Figs 9 and 10 and compared with experimental data at two different times. The consistency between the numerical results and the experimental data is clearly very good, and it is encouraging that features such as the localized peak seen on 9/10/2006 is replicated by the model. This result lends considerable support to the hydro-mechanical results, and provides a high degree of confidence that the process representations are well captured by the numerical tool.

Reactive (sulfate) analysis

A prototype model was developed using the QPAC (Quintessa, 2012) reactive transport module, designed to investigate the plausibility of a redox pyrite–gypsum reaction in describing the observed evolution of the sulfate. The reactive transport module has been designed to model interactions of groundwater and other subsurface fluids with rocks and man-made materials. In an open system (i.e. systems other than ‘closed-box’ batch-type systems) the interactions in the water–rock system are represented by non-linear reactive transport equations, which couple the fluid flow and advective-diffusive transport equations to equations representing the geochemical reactions between the pore water components and the solid materials. The model is ‘fully coupled’ (rather than being implemented as a two-step process as in some other modelling codes) and allows alteration processes in the rock and man-made materials to feed back into the fluid flow equations through variations in porosity and other material properties such as permeability and tortuosity (although this option is not used in the current calculations). Geochemical porosity (e.g. ‘chloride porosity’) was neglected in this case because data on the probable geochemical porosity for all the aqueous species were not available.

For these initial calculations it was assumed that the feedback of the chemical processes on the hydraulic calculations is weak and that the rate at which oxygen is dissolved into the water is not sufficient to significantly alter the O2 partial pressure in the pores due to a relatively rapid resupply of O2 through the connected gas (air) body compared to other processes in the system. Darcy pore water fluxes are obtained from the hydro-mechanical flow modelling, and use the same flows as for the chloride calculations.

The geochemical subsystem was deliberately kept simple for these prototype calculations. The key process of interest is the oxidation of pyrite (FeS2), as this is likely to exert the strongest control on sulfate concentrations in the pore water. The pyrite oxidation reaction can be represented as:  
FeS2(s)+3.75O2(aq)+0.5H2OH++2SO42+Fe3+
(5)

Overall, pyrite oxidation will release sulfate to the pore water and reduce pH. Sulfate concentrations in the pore water will vary according to the rate at which sulfate is supplied from oxidation of pyrite and the rate at which it is transported through the system by the bulk movement of water through the pore space and diffusion in the bulk water. Sulfate concentrations are expected to be highest near the tunnel wall where there is a strongly connected gas phase. Dissolved O2 concentrations will be transported similarly but are also subject to a source of oxygen from dissolution wherever free gas is present in the pore space.

In addition to pyrite, calcite (CaCO3(s)), quartz (SiO2(s)) and albite (NaAlSi3O8) were included in the model to represent the solids present in the host rock, with dissolution/precipitation modelled kinetically. Calcite tended to buffer pH in the model by dissolving. Quartz and albite were expected to be relatively inert.

Although the intention was to keep the geochemical system simple, in order to properly represent the redox system, around 50 aqueous geochemical species were included per compartment as listed below. A temperature of 25°C was assumed throughout to ensure availability of thermodynamic data, noting that this temperature is not very different from the 15°C observed at the tunnel wall, and in the context of all the other uncertainties regarding the geochemical system this is assumed to be a relatively minor issue. The aqueous species are included in the model are shown in Table 1. LogK data for all redox and complex species are taken from thermo.com.V8.R6.230 (Johnson et al., 2000) at 25°C.

The key output of interest is the sulfate content in the porewater, as this can be directly compared to measured data reported by Fernández and Melón (2007a,b). Measured and simulated sulfate concentrations (mg kg−1 rock) are shown in Figs 11, 12 and 13. Given the scatter in the measured data, the fit provided by the numerical model is reasonably good and is well within the range of the measured data noting that the geochemical modelling indicated a ‘background’ sulfate concentration of 2.58 mg kg−1 rock in water saturated clay. There are few data points at 5/7/2002 with the simulated sulfate content underestimating concentrations at the tunnel boundary but matching reasonably well in the 0.3–0.8 m range. At 26/1/2004 the model over-predicts concentrations with respect to the bulk of the data points in the 0.2–0.9 m range, but falls about midway between the measurements in the first 0–0.2 m. By 9/10/2006 the spread in the measured data is significant.

The reasonable fit despite ignoring the concept of a reduced geochemical porosity may lead to the conclusion that multi-component transport processes are occurring, with different aqueous species being advected at different rates. However, sulfate concentrations are likely to be affected mostly by mineral interactions (notably pyrite dissolution) with transport being a secondary effect. This can be seen in the concentration profiles in Figs 11, 12 and 13, where the ‘hump’ in the profile is mostly a consequence of localized pyrite dissolution with the localized peak near the tunnel boundary being primarily a consequence of advective transport. If the assumption were made that sulfate behaved like chloride (i.e. it is not present in all calculated pore space) it would be likely to most markedly affect results where advection processes dominate (i.e. near the tunnel boundary). It would have the effect of reducing the local peak, which might lead to a better fit to the non-peak sulfate concentrations that were measured, but would not be expected to have a large effect on the main sulfate peak close to the tunnel wall.

Conclusions

A hydro-mechanical representation of the Mont-Terri ventilation experiment and the supporting drying test has been developed which replicates all the major characteristics of the system, and has a parameterization which is fully consistent between the laboratory and field scales. Furthermore, in order to be able to replicate the true boundaries of the ventilation experiment, a simplified representation of the tunnel was implemented which was found to be remarkably effective and robust in replicating the interaction between the tunnel and the surrounding Opalinus Clay.

In addition, a simple representation of the key geochemical properties that are likely to have operated has been developed, namely, kinetic treatments of pyrite dissolution (releasing dissolved sulfate and reducing pH), calcite dissolution (pH buffering); and solute transport of reactive and non-reactive species. The model includes a simplified representation of mineralogy and porewater compositions and a relatively simple treatment of reactive mineral surface area. However, based on the comparison of measured data with simulated data, the models provide a good representation of key processes. The fit of simulated dissolved chloride concentrations data to measured values is especially encouraging, given the simplification of the system and this gives additional confidence in the basic hydraulic modelling approach.

Acknowledgements

The authors would like to acknowledge the support of Cherry Tweed and Simon Norris at NDA RWMD. In addition the authors wish to recognize the technical input of Benoit Garitte at UPC, Spain and Alain Millard at CEA, France in developing the understanding of many aspects of the Ventilation Experiment presented in this paper.