Search

We conduct coupled thermo-hydro-mechanical modeling of a KBS-3V repository design in crystalline rocks, using data and conditions from the Forsmark in Olkiluoto repository sites in Sweden and Finland. The study focuses on repository performance related to the impact of thermal and hydraulic evolution on the potential for thermal–mechanical damage to underground repository excavations. For the designs and conditions considered at the Forsmark and Olkiluoto repository sites, the simulations show a peak temperature well under the adopted performance target of a 100°C maximum temperature, whereas there is still a high potential for thermal–mechanical damage to the KBS-3V waste deposition holes. The thermal–mechanical damage is much more likely if rock permeability is so low that it delays saturation and swelling of bentonite-clay-based backfill beyond the time for the thermal–mechanical peak, which occurs 50 to 100 years after nuclear waste deposition. We also found that sidewalls of the KBS-3V emplacement tunnels are vulnerable to tensile fracturing due to the combined effect of thermal stressing and backfill swelling. The study highlights a strong interaction between bentonite-based backfill and host rock through capillary suction along with induced rock desaturation. A careful design and selection of the bentonite-clay-based backfill materials for KBS-3V tunnels and deposition holes can facilitate a timely saturation and backfill swelling that in turn can minimize thermal–mechanical damage.

Understanding the impact of repository gas, generated from degradation of waste and its interaction with the host rock, is essential when assessing the performance and safety function of long-term disposal systems for radioactive waste. Numerical models based on conventional multi-phase flow theory have historically been applied to predict the outcome and impact of gas flow on different repository components. However, they remain unable to describe the full complexity of the physical processes observed in water-saturated experiments (e.g., creation of dilatant pathways) and thus, the development of novel representations for their description is required when assessing fully saturated clay-based systems. This was the primary focus of Task A within the international cooperative project DECOVALEX-2019 (D-2019) and refinement of these approaches is the primary focus of this study (Task B in the current phase of DECOVALEX-2023). This paper summarises development of enhanced numerical representations of key processes and compares the performance of each model against high-quality laboratory test data. Experimental data reveals that gas percolation in water-saturated compacted bentonite is characterised by four key features: (i) a quiescence phase, followed by (ii) the gas breakthrough, which leads to a (iii) peak value, which is then followed by (iv) a negative decay. Three models based on the multiphase flow theory have been developed. These models can provide good initial values and reasonable responses for gas breakthrough (although some of them still predict a too-smooth response). Peak gas pressure values are in general reasonably well captured, although maximum radial stress differences are observed at 48 mm from the base of the sample. Here, numerical peak values of 12.8 MPa are predicted, whereas experimental values are about 11 MPa. These models are also capable of providing a reasonable representation of the negative pressure decay following peak pressure. However, other key specific features (such as the timing of gas breakthrough) still require a better representation. The model simulations and their comparison with experimental data show that these models need to be further improved with respect to model parameter calibration, the numerical representation of spatial heterogeneities in material properties and flow localisation, and the upscaling of the related physical processes and parameters. To further understand gas flow localisation, a new conceptual model has been developed, which shows that discrete channels can possibly be induced through the instability of gas-bentonite interface during gas injection, thus providing a new perspective for modeling gas percolation in low-permeability deformable media.

Controlling brine flow is important for the safe disposal of radioactive waste in salt rocks. Thermal and mechanical processes can have a significant impact on brine flow, although this has never been thoroughly investigated. In this study, we conduct fully coupled THM modeling for analyzing brine flow in rock salt, considering non-isothermal two-phase flow through deformable porous media. To rigorously represent rock salt behavior, we incorporate suitable phenomenological models for creep and shear and tensile-induced dilatancy/damage, and their effect on the flow properties. To validate such a complex model, we use it to analyze an experiment on a meter-scale salt block subjected to multistage heating and cooling under controlled laboratory conditions. The experimental data and our model predictions of temperature and brine inflow show good agreement. Our modeling shows that it is important to consider the coupling between the heating- and cooling-induced damage and the flow properties for rigorously estimating brine inflow. Based on our modeling, a cooling-induced brine inflow spike is caused by a permeability increase due to tensile dilatancy. Thus, our modeling will be useful for a proper design of a cooling phase in the salt repositories to avoid or minimize damage and the resulting brine inflow.

A key feature of unconventional shale reservoirs is low permeability. On average, only 10% of original oil in place and 25% of original gas in place is produced from these reservoirs. Understanding the mechanisms behind fracture conductivity past the first couple of years of shale production is still an ongoing research topic, as most such reservoirs report fast and massive production declines. Hydraulic fracturing is indispensable for stimulating low-permeability rocks for economical oil and gas production from unconventional reservoir rocks such as shales. In contrast, for subsurface sequestration and storage of fluids including supercritical CO2 (for carbon capture, utilization, and storage) and hydrogen gas, inadvertent creation of fractures can lead to breaching of shale caprock. The Caney Shale in southern Oklahoma is being evaluated as both a potential hydrocarbon-producing reservoir and a caprock for fluid sequestration. The Caney Shale is a Mississippian-age, organic-rich mudrock with intermittent calcareous laminae, that produces. But currently there is limited data reported regarding productivity of horizontally drilled wellbores. While many scholars have investigated fracture conductivity in shale reservoirs, the mechanisms of proppant embedment in relation to lithology are still under investigation. This study implemented a multi-scale approach towards investigating proppant embedment and fracture conductivity, from nano-scale instrumented micro/nano-indentation to millimeter (mm) scale, mono-layer propped fracture flow at reservoir temperature and pressure, and American Petroleum Institute (API)-RP19D conductivity tests using inch/cm-scale shale platens. At each scale, various material characterization tools were utilized, including Focused Ion Beam-Scanning Electron Microscopy (FIB-SEM), Energy Dispersive Spectroscopy (EDS), Raman spectroscopy, Laser profilometry, Computed X-ray Microscopy and X-ray Diffraction (XRD). The outcomes of the micro-indentation revealed variations in mechanical properties attributed to changing mineral composition and microstructures. Outcomes from API fracture conductivity testing and flow-through testing using a monolayer of proppant demonstrated: (1) a 50% impact on proppant embedment compared to ductile samples, (2) a significant decline in fracture conductivity with increasing stress and temperatures, and (3) conductivity that was also influenced by organic and inorganic content as well as internal sample architecture. Mineralogically, ductile samples contain about 25% more clay minerals compared to brittle regions. Heterogeneity in mineral composition causes the Caney Shale to show different responses at reservoir temperature and pressure, particularly in relation to creep and proppant embedment. Geomechanical and fluid-flow modeling have been conducted at several scales to help interpret laboratory results and apply findings to field-scale operation. The outcomes from this integration should aid geological and engineering predictions for the hydrocarbon production and caprock integrity of the Caney Shale.

The effectiveness of Carbon Capture and Storage (CCS) as an imperative decarbonization technology relies on the sealing capacity of a fine-grained caprock to permanently store CO2 deep underground. Uncertainties in assessing the caprock sealing capacity increase with the spatial and temporal scales and may delay CCS deployment at the gigatonne scale. We have developed a computationally efficient transport model to capture the dynamics of basin-wide upward CO2 migration in a multi-layered setting over geological time scales. We find that massive capillary breakthrough and viscous flow of CO2, even through pervasively fractured caprocks, are unlikely to occur and compromise the storage security. Potential leakage from the injection reservoir is hampered by repetitive layering of overlying caprocks. This finding agrees with geologic intuition and should be understandable by the public, contributing to the development of climate policies around this technology with increased confidence that CO2 will be indefinitely contained in the subsurface.

This study investigates the effect of shale creep on the stress evolution of a generic nuclear waste repository. Thermohydromechanically (THM) coupled simulations were carried out with the TOUGH-FLAC simulator to model temperature/pore pressure/stress changes due to the decay heating from nuclear waste. The Norton-Bailey creep model, which was calibrated against published creep experiments on shale, was used to simulate different creep properties of shale (i.e., no creep vs. creep). Results show that shale without creep generated and maintained high stress concentrations (~40 MPa) near the disposal tunnel for an extremely long period (~10, 000 years), whereas shale with creep dissipated the stress concentrations and developed uniform and isotropic stress field of roughly 12 MPa (in compression), which corresponded with the overburden stress level at the tunnel depth, by as early as 100 years since the emplacement of nuclear waste. Also, it was found that shale with creep increased compressive stresses in the bentonite backfill compared to shale without creep by as much as 130% at 10, 000 years due to creep-induced tunnel volume reduction.

Oil and gas shales are a system where understanding the mechanical properties at the microscale is of paramount importance, e.g. to better understand the behavior of proppant-shale contacts and their role in the evolution of propped fractures in unconventional reservoirs. This work shows for the first time an experiment coupling indentation testing with in situ X-ray imaging in a Green River shale sample. A full compliance curve has been measured with the sample in water, allowing to visualize the indentation of the sample in function of axial load, in a purpose-built system for combined in situ indentation and X-ray imaging. A series of 3D datasets were used for a digital volume correlation study to obtain local strain fields. This analysis has been complemented with the analysis of cracks. Finally, geomechanical modeling has been carried out to replicate and generalize the observed behavior in the shale. This study validated this experimental approach, providing a breakthrough in understanding micro-mechanics in shales, and demonstrates how this class of experiments can be important for studies involving the prediction of the evolution of propped fractures in shale reservoirs, with possible applications in a much larger number of application fields (geothermal, materials science, etc.)

The complex geometric features of subsurface fractures at different scales makes mesh generation challenging and/or expensive. In this paper, we make use of neural style transfer (NST), a machine learning technique, to generate mesh from rock fracture images. In this new approach, we use digital rock fractures at multiple scales that represent ‘content’ and define uniformly shaped and sized triangles to represent ‘style’. The 19-layer convolutional neural network (CNN) learns the content from the rock image, including lower-level features (such as edges and corners) and higher-level features (such as rock, fractures, or other mineral fillings), and learns the style from the triangular grids. By optimizing the cost function to achieve approximation to represent both the content and the style, numerical meshes can be generated and optimized. We utilize the NST to generate meshes for rough fractures with asperities formed in rock, a network of fractures embedded in rock, and a sand aggregate with multiple grains. Based on the examples, we show that this new NST technique can make mesh generation and optimization much more efficient by achieving a good balance between the density of the mesh and the presentation of the geometric features. Finally, we discuss future applications of this approach and perspectives of applying machine learning to bridge the gaps between numerical modeling and experiments.

The coupled Thermo-Hydro-Mechanical (THM) behavior of the Callovo-Oxfordian claystone (COx) is of great importance for the design and safety calculations of the high-level radioactive waste disposal project in this potential host rock in France. The heat emitted by the waste causes a pore pressure increase within the surrounding rock essentially due to the differential thermal expansion of the pore water and the solid skeleton. The low permeability of the COx and its relative rigidity inhibits the discharge of the induced pressure build-up. Moreover, thermal loading may provoke thermo-mechanical stresses within the formation due to mechanical confinement by the rigidity of the surrounding host rock. An important research program has been conducted by the French National Radioactive Waste Management Agency (Andra) since 2003 in order to investigate the THM response of the COx under thermal loading, through laboratory tests, in-situ experiments, model development and numerical modeling. Within Task E of the DECOVALEX-2019 project, five research teams investigated upscaling THM modeling from a small-scale in-situ experiment (TED) to a full-scale in-situ experiment (ALC). The upscaling modeling started with a verification test to validate the numerical codes. Then, an interpretative modeling of the TED experiment was performed to calibrate the THM parameters of the COx. Finally, the calibrated THM parameters were used for a blind prediction of the ALC experiment. The modeling teams each adopted a thermo-poro-elastic approach which yielded satisfactory results. The blind prediction of the temperature field showed an overestimation of less than 2 °C which was considered acceptable. On the other hand, pore pressure was well predicted only in the direction parallel to the bedding whereas the slow dissipation of the pore pressure in the direction perpendicular to the bedding was not captured by any of the modeling teams – which remains an open question of the present study.

This paper is devoted to the study of the Thermo-Hydro-Mechanical (THM) responses of a porous rock with low permeability under thermal loading in the context of deep geological disposal of radioactive waste. To this aim, numerical simulations of a benchmark exercise of a hypothetical high-level radioactive waste (HLW) repository were performed. This benchmark exercise considered as a host formation the Callovo-Oxfordian claystone (COx), which has been selected for a deep geological disposal in France. Within the framework of the DECOVALEX-2019 project, five modelling teams (Andra, LBNL, NWMO, Quintessa, UFZ/BGR) adopted a thermo-poro-elastic approach and proposed different 3D representations of the HLW repository. The differences between the teams consisted mostly in the simplification of the geometrical model and the interpretation of the boundary conditions. Numerical results for temperature, pore pressure, and effective stress evolution in the far field (i.e., at the mid-distance of two HLW cells) were compared between the teams, to quantify the impact of modelling simplifications/assumptions for the assessment of the HLW repository. The THM behaviour of the COx formation in the near field (i.e., excavation damaged zone around the HLW cells) is not the objective of this study. Moreover, plane strain conditions were considered and evaluated in comparison to 3D modelling. Key parameters influencing the THM responses of the HLW repository were assessed by both mono- and multi-parametric analyses. Spatial variability analyses of THM parameters were also carried out to study the influence of the spatial correlation length on the Terzaghi effective stress and to estimate its probability distribution. The conclusions of this study provide reliable numerical techniques for modelling large-scale deep geological disposals and deduce the main behavior of the HLW repository.

Experimental studies show that compacted bentonite used as a backfill material for nuclear waste repository experiences strong coupling between chemical and mechanical processes. In this paper, we use a dual-structure expansive soil model, referred to as the Barcelona Expansive Model (BExM), to predict the behavior of bentonite buffer in subsurface emplacement tunnels for high-level radioactive waste. After emplacement, the bentonite is subjected to complex and coupled Thermal-Hydraulic-Mechanical-Chemical (THMC) processes. The BExM constitutive model is implemented in a multi-phase reactive transport and geomechanics simulator, TOUGHREACT-FLAC3D, and the model is verified with one benchmark test on FEBEX bentonite. We utilize a one-way Chemo-Mechanical (C-M) coupling approach, in which chemical changes affect the mechanical behavior of bentonite through the BExM linked with the evolution of mass fraction of smectite, exchangeable cation concentration, and ionic strength via osmotic suction. The parameters of the new coupled model for FEBEX bentonite are calibrated against a series of laboratory experiments with various salinity solutions. Finally, coupled THMC modeling is conducted for a generic argillite repository with bentonite buffer under high temperature, focusing on the long-term chemical change and its effect on the mechanical process.

In a repository for radioactive waste hosted in a clay formation, hydrogen and other gases may be generated due to the corrosion of metallic materials under anoxic conditions, the radioactive decay of waste and the radiolysis of water. If the gas production rate exceeds the gas diffusion rate within the pores of the clay, a discrete gas phase will form and accumulate until its pressure becomes large enough to exceed the entry pressure of the surrounding material, at which point dilatant, advective flow of gas is expected to occur. The purpose of Task An under DECOVALEX-2019 is to better represent the processes governing the advective movement of gas in both low-permeability argillaceous repository host rocks and clay-based engineered barriers within numerical codes. In this paper special attention is given to the mechanisms controlling gas entry, flow and pathway sealing and their impact on the performance of the engineered clay barrier. Previous work suggests gas flow is accompanied by the creation of dilatant pathways whose properties change temporally and spatially within the medium. Thus, four new types of approaches have been developed: (i) standard two-phase flow models (continuous techniques) incorporating a range of different mechanical deformation behaviours, (ii) enhanced two-phase flow models in which fractures are embedded within a plastic material (continuous techniques) or incorporated into the model using a rigid-body-spring network (discrete approaches), (iii) a single-phase model incorporating a creep damage function in which only gas flow is considered, and (iv) a conceptual approach used to examine the chaotic nature of gas flow. The outputs from these different approaches are compared. This is an essential step as the choice of modelling approach strongly impacts the representation and prediction of gas flow in a future repository. In addition, experience gained through this task is of direct relevance to other clay-based engineering issues where immiscible gas flow is a consideration including hydrocarbon migration, carbon capture and storage, shale gas and landfill design. This paper summarises the outcomes of work in Task A conducted between May 2016 and May 2019 and provides a brief overview of the experimental data and a synthesis of the work of the participating modelling teams.

The greatest challenges of rigorously modeling coupled hydro-mechanical processes in fractured rocks at different scales are associated with computational geometry. In addition, selections of continuous or discontinuous models, physical laws, and coupling priorities at different scales based on different geometric features determine the applicability of a numerical model for a certain type of problem. In this study, we present our multi-scale modeling capabilities that have been developed based on the numerical manifold method for analyzing coupled hydro-mechanical processes in fractured rocks. Based on their geometric features, the fractures are modeled as continua—finite-thickness porous zones, and discontinua—discontinuous interfaces and microscale asperities and granular systems. Different governing equations, physical laws, coupling priorities, and approaches for addressing fracture intersections and shearing are then applied to describe these. We applied these models to simulate coupled processes in fractured rocks using realistic geometry obtained from rock images at different scales. We first calculated shearing of a single fracture with different models and demonstrated the impacts of asperities on shearing. We then applied the continuous and discontinuous models to simulate a network of rough fractures, demonstrating that contact dynamics contribute significantly to the geometric, multi-physical evolution of systems where rough fractures are not mineral filled. For a discrete fracture network, our coupled processes modeling demonstrates that shearing of the discrete fractures can have a major impact on stress and pore pressure distribution. Lastly, we applied the discontinuous granular model to simulate evolution of a complex granular system with a deformation band, demonstrating that the deformation band can dominate contact dynamics, the structural and the stress evolution of the granular system.

Geologic disposal is a promising solution for a safe permanent isolation of accumulated high-level nuclear waste from nuclear powerplants. The behavior of host rock is highly coupled thermally, hydromechanically, and chemically. Numerical simulations of such coupled phenomena for the extremely long term (>100,000 years) and large length scale (>1 km) of geologic disposal remain to be computationally challenging. In this study, a methodology has been developed to approximate stress and stress-induced permeability change in host rock using only thermo-hydrological (TH) variables. A coupled thermo-hydromechanical (THM) simulation is carried out using TOUGH-FLAC simulator to model THM behaviors of a generic nuclear waste repository, in order to evaluate the performance of the developed methodology, which is implemented in a coupled TH simulation. Results show that stress and permeability change estimated by the developed methodology in the TH-coupled simulation match those calculated in the THM-coupled simulation over the simulated timeframe of over 10,000 years. Details about errors in stress and permeability estimates accrued by the developed methodology are discussed in this paper. The developed methodology will help incorporate stress-induced permeability change into existing TH simulators for the long-term radionuclide transport in geologic disposal.

A coupled multiphase fluid flow and discrete fracturing model is applied to simulate bench-scale gas migration experiments on compacted bentonite. The numerical modeling is based on the linking of the multiphase fluid flow simulator TOUGH2 with a Rigid-Body-Spring Network model, which enables a discrete (lattice) representation of elasticity and individual fractures. The evolution of a complex network of dilatant flow paths is modeled through opening and breakage of lattice interface bonds between porous-elastic matrix elements. To achieve a good match with the experimental results, including an abrupt gas breakthrough along with pressure and stress responses, it was necessary to calibrate model parameters for (1) air-entry pressure, (2) shear and tensile failure of lattice interface bonds, (3) moisture swelling/shrinkage effects on stress, and (4) aperture-dependent permeability of dilatant flow paths. Our best-fit conceptual model considers a pervasive network of discrete flow paths propagating from the gas injection point, whereas some of the experimental data indicate the potential for heterogeneous and unstable flow paths.

This study presents comprehensive coupled hydro-mechanical numerical modeling of the first and second hydraulic stimulations at the Pohang enhanced geothermal system (EGS) site in order to improve the understanding on the key stimulation mechanisms in the fractured reservoir. Two models for PX-2 and PX-1 wells were developed including fracture zones of different orientations and permeability. Shear dilation and frictional plastic strain-softening were implemented in a non-linear stress dependent fracture model. A history matching of the wellhead pressure curves was carried out in the early days of each stimulation. The PX-2 model with a fracture jacking mechanism achieved the reproduction of highly reversible non-linear aperture changes observed in situ. The PX-1 model successfully simulated a wellhead pressure peak and drop associated with shear slip and dilation by applying frictional plastic strain softening. Simulated critical wellhead pressure for shear slip near the PX-1 well was greater than the prediction by a simple slip potential analysis due to local stress changes by the poroelastic effect. The numerical modeling confirmed that the combination of shear dilation and jacking was adequate to capture the pressure evolution during the increasing step-rate test at PX-1. Besides, spatio-temporal changes in pressure, total stress, and permeability evaluated in the fracture zones greatly enhanced the understanding on the coupled behavior caused by the hydraulic stimulations. The possibility of zones with lower permeability away from the PX-1 well was suggested by an alternative model. The current numerical study demonstrates that the key hydro-mechanical processes of shear slip and dilation and hydraulic jacking observed in the fractured reservoir can be successfully reproduced.

Geodetic observations, often in conjunction with other data, provide a cost-effective means for identifying and characterizing geothermal resources. A review of the various methods reveals how the technology for measuring deformation has advanced considerably in the past few decades. Currently, interferometric synthetic aperture radar is the method of choice for monitoring deformation at a geothermal field. A discussion of geodetic monitoring at The Geysers geothermal field, California, illustrates some of the progress made and the challenges that remain.

We present results from an international model comparison study involving a series of controlled fault activation experiments in Opalinus Clay at the Mont Terri Laboratory, Switzerland. The fault activation experiments were conducted in situ by water injection at variable pressure from boreholes targeting different parts of the Main Fault crossing several tunnels and galleries of the Mont Terri Laboratory. The model simulations focused on (1) an experiment activating a discontinuity (fracture or minor fault) within the damage zone of the main fault, and (2) an experiment activating discontinuities close to the core of the Main Fault. The experimental data consist of coupled hydraulic and mechanical responses monitored at an injection borehole and a monitoring borehole located a few meters away. After overcoming several modeling issues along with necessary model developments, a reasonably good agreement was achieved between the modeling results and the field observations. The fault activation experiments displayed an abrupt flow rate increase associated with a sudden fracture opening and rupture propagation after the injection pressure reached above the estimated normal stress on the fracture. This was followed by an abrupt flow rate decrease, indicating hydraulic closing, once the injection pressure decreases to an estimated 1–2 MPa below the stress normal to the opened fracture. The models were able to capture this abrupt hydromechanical behavior, including an observed dominant opening behavior along with the rupture propagation, while the spatial extent of the shear rupture and the quantity of peak injection flow were the most challenging to predict.

We present coupled Thermo–Hydro–Mechanical (THM) modeling of geologic nuclear waste disposal in argillaceous claystone, focusing on thermally-induced pressure changes and the potential for such pressure changes to induce hydro-fracturing. To investigate this possibility, we first conduct a three-dimensional repository scale model, with host rock properties, repository design, and nuclear waste decay heat functions derived from the French concept of geologic disposal in argillaceous claystone. The model simulations show that the highest potential for hydro-fracturing occurs between emplacement micro-tunnels (cells) at the center rather than at the edge of the repository. We further investigate the use of a two-dimensional single cell model as a simpler surrogate for a full three-dimensional model. Our results reveal that such geometric simplification is reasonably accurate for modeling coupled THM processes at the center of the repository, though it overestimates the likelihood for hydro-fracturing and substantially overpredicts ground surface uplift. A parameter study shows the importance of adjacent higher permeability geological layers that play a significant role in dissipating overpressure in the host rock layer. Finally, the study shows the importance of the spacing between emplacement tunnels, which if too short results in a higher fluid pressure and a strongly increased potential for hydro-fracturing. Overall, the study suggests, the limiting factor in the thermal management and design of a repository is not necessarily the maximum temperature in the engineered barrier system near the waste packages, but rather the more modest host rock temperature between emplacement tunnels, due to the potential for thermal damage.

We investigate the influence of foliation orientation and fine-scale folding on the static and dynamic elastic properties and unconfined strength of the Poorman schist. Measurements from triaxial and uniaxial laboratory experiments reveal a significant amount of variability in the static and dynamic Young’s modulus depending on the sample orientation relative to the foliation plane. Dynamic P-wave modulus and S-wave modulus are stiffer in the direction parallel to the foliation plane as expected for transversely isotropic mediums with average Thomsen parameters values 0.133 and 0.119 for epsilon and gamma, respectively. Static Young’s modulus varies significantly between 21 and 117 GPa, and a peculiar trend is observed where some foliated sample groups show an anomalous decrease in the static Young’s modulus when the symmetry axis (x3-axis) is oriented obliquely to the direction of loading. Utilizing stress and strain relationships for transversely isotropic medium, we derive the analytical expression for Young’s modulus as a function of the elastic moduli E1, E3, ν31, and G13 and sample orientation to fit the static Young’s modulus measurements. Regression of the equation to the Young’s modulus data reveals that the decrease in static Young’s modulus at oblique symmetry axis orientations is directly influenced by a low shear modulus, G13, which we attribute to shear sliding along foliation planes during static deformation that occurs as soon as the foliation is subject to shear stress. We argue that such difference between dynamic and static anisotropy is a characteristic of near-zero porosity anisotropic rocks. The uniaxial compressive strength also shows significant variability ranging from 21.9 to 194.6 MPa across the five sample locations and is the lowest when the symmetry axis is oriented 45° or 60° from the direction of loading, also a result of shear sliding along foliation planes during static deformation.

Micromechanical modeling of geomaterials is challenging because of the complex geometry of discontinuities and potentially large number of deformable material bodies that contact each other dynamically. In this study, we have developed a numerical approach for micromechanical analysis of deformable geomaterials with dynamic contacts. In our approach, we detect contacts among multiple blocks with arbitrary shapes, enforce different contact constraints for three different contact states of separated, bonded, and sliding, and iterate within each time step to ensure convergence of contact states. With these features, we are able to simulate the dynamic contact evolution at the microscale for realistic geomaterials having arbitrary shapes of grains and interfaces. We demonstrate the capability with several examples, including a rough fracture with different geometric surface asperity characteristics, settling of clay aggregates, compaction of a loosely packed sand, and failure of an intact marble sample. With our model, we are able to accurately analyze (1) large displacements and/or deformation, (2) the process of high stress accumulated at contact areas, (3) the failure of a mineral cemented rock samples under high stress, and (4) post-failure fragmentation. The analysis highlights the importance of accurately capturing (1) the sequential evolution of geomaterials responding to stress as motion, deformation, and high stress; (2) large geometric features outside the norms (such as large asperities and sharp corners) as such features can dominate the micromechanical behavior; and (3) different mechanical behavior between loosely packed and tightly packed granular systems.

We apply thermoporoelasticity and a sequentially coupling technique for modeling thermally-driven coupled Thermo-Hydro-Mechanical (THM) processes in tight claystone. A THM benchmark case with a corresponding analytic solution for thermoporoelasticity under a constant heat loading verifies the model. Thereafter, two in situ heating experiments are simulated for model validation: a smaller-scale heating experiment (TED experiment) and a larger-scale experiment (ALC experiment) in Callovo-Oxfordian (COx) claystone at the Meuse/Haute-Marne underground research laboratory in France. The model exhibits good performance to match the observed temperature and pore pressure evolution for the smaller-scale TED experiment. For the larger-scale ALC experiment, general trends of thermal-pressurization are captured in the modeling, but pressure is underestimated at some monitoring points during cool-down. This indicates that the THM response in the field may be affected by the variability of rock's properties or irreversible or time-dependent mechanical processes that are not included in the current thermoporoelastic model. The main contributions of this work are as follows: (1) we verify and validate the numerical simulator, TOUGH-FLAC, to be a valuable coupled THM modeling tool; (2) prove that the laboratory determined material parameters can be used as reference values for upscaling experiments. However, to better identify and quantify THM processes with modeling of in situ tests, more emphasize should be dedicated to obtaining high-quality mechanical deformation data.

The greatest challenges of rigorously modeling coupled hydro-mechanical (HM) processes in fractured geological media at different scales are associated with computational geometry. These challenges include dynamic shearing and opening of intersecting fractures at discrete fracture scales as a result of coupled processes, and contact alteration along rough fracture surfaces that triggers structural and physical changes of fractures at micro-asperity scale. In this paper, these challenges are tackled by developing a comprehensive modeling approach for coupled processes in fractured geological media based on numerical manifold method (NMM) at multiple scales. Based on their distinct geometric features, fractures are categorized into three different scales: dominant fracture, discrete fracture, and discontinuum asperity scales. Here the scale is relative, that of the fracture relative to that of the research interest or domain. Different geometric representations of fractures at different scales are used, and different governing equations and constitutive relationships are applied. For dominant fractures, a finite thickness zone model is developed to treat a fracture as a porous nonlinear domain. Nonlinear fracture mechanical behavior is accurately modeled with an implicit approach based on strain energy. For discrete fractures, a zero-dimensional model was developed for analyzing fluid flow and mechanics in fractures that are geometrically treated as boundaries of the rock matrix. With the zero-dimensional model, these fractures can be modeled with arbitrary orientations and intersections. They can be fluid conduits or seals, and can be open, bonded or sliding. For the discontinuum asperity scale, the geometry of rough fracture surfaces is explicitly represented and contacts involving dynamic alteration of contacts among asperities are rigorously calculated. Using this approach, fracture alteration caused by deformation, re-arrangement and sliding of rough surfaces can be captured. Our comprehensive model is able to handle the computational challenges with accurate representation of intersections and shearing of fractures at the discrete fracture scale and rigorously treats contacts along rough fracture surfaces at the discontinuum asperity scale. With future development of three-dimensional (3D) geometric representation of discrete fracture networks in porous rock and contacts among multi-body systems, this model is promising as a basis of 3D fully coupled analysis of fractures at multiple scales, for advancing understanding and optimizing energy recovery and storage in fractured geological media.

The present study is aimed at developing a numerical model to reproduce coupled hydro-mechanical processes associated with fault reactivation by fluid injection in low permeability rock, as part of the DECOVALEX-2019 project Task B. We proposed a modeling approach for simulating the processes using the TOUGH-FLAC simulator, and modeled a fault reactivation experiment conducted at Mont Terri Rock Laboratory in Switzerland. The first step of the study involved benchmark calculations considering a simplified fault plane and geometry. Fluid flow along a fault was modeled using elements of aperture-sized thickness on the basis of Darcy's law and the cubic law in TOUGH2, whereas the mechanical behavior of a single fault was represented by zero-thickness interface elements in FLAC3D upon which a slip and/or separation is allowed. A methodology to connect a TOUGH2 volume element to a FLAC3D interface element was developed for handling the hydro-mechanical interactions on the fault during fluid injection. Two different fault models for describing the evolutions of hydraulic aperture by elastic fracture opening and failure-induced aperture increase were considered in the benchmark calculations. In the coupling process, the changes in geometrical features and hydrological properties induced by mechanical deformation were continuously updated. The transient responses of the fault and host rock to stepwise pressurization were examined during the simulation. The hydro-mechanical behavior, including the injection flow rate, pressure distribution around the borehole, stress conditions, and displacements in normal and shear directions were monitored in the surrounding rock and along the fault. The results of benchmark calculations suggest that the developed model reasonably represents the hydro-mechanical behavior of a fault and the surrounding rock. This modeling approach was applied to the fault reactivation experiment of the Mont Terri Rock Laboratory. In this interpretive modeling, a parametric study was conducted to examine the effects of input parameters regarding in situ stress and fault properties on the hydro-mechanical responses of the fault to water injection. Then, an optimal parameter set to reproduce the field experiment results was chosen by trial-and-error. The injection flow rate and pressure response during fault reactivation closely matched those obtained at the site, which indicates the capability of the model to appropriately capture the progressive pathway evolution during fault reactivation tests at the site. The anchor displacements were overestimated by the model, but a fair agreement was obtained in terms of the order of magnitude and the variation tendency.

Coupled thermal-hydro-mechanical (THM) numerical modeling is conducted to study thermal management associated with geologic disposal of spent nuclear fuel (SNF) in large dual-purpose canisters (DPCs). DPCs, are containers designed for SNF storage and transportation and if determined to be feasible for permanent geological disposal could provide a cost effective disposal solution. However, one of the challenges to direct disposal of DPCs is thermal management to avoid overheating of the Engineered Barrier System (EBS), including bentonite-backfill used as a protective buffer. The model simulations show that the use of a backfill that is thermally engineered for high heat conduction can reduce the EBS temperature to acceptable levels for disposal of large waste canisters in backfilled tunnels. On the other hand, the use of high heat conduction backfill will not reduce the far field rock peak temperature that can occur several thousand years after closure of the repository. This longer term host rock peak temperature generates thermal-poro-elastic stress and geomechanical changes that must be considered in the thermal management and design of a repository.

In this study, the combined effects of geometrical distribution and geomechanical deformation of fracture networks on fluid flow through fractured geological media are investigated numerically. We consider a finite-sized model domain in which the geometry of fracture systems follows a power-law length scaling. The geomechanical response of the fractured rock is simulated using a hybrid finite-discrete element model, which can capture the deformation of intact rocks, the interaction of matrix blocks, the displacement of discrete fractures and the propagation of new cracks. Under far-field stress loading, the locally variable stress distribution in the fractured rock leads to a stress-dependent variable aperture field controlled by compression-induced closure and shear-induced dilatancy of rough fractures. The equivalent permeability of the deformed fractured rock is calculated by solving for the fracture-matrix flow considering the cubic relationship between fracture aperture and flow rate at each local fracture segment. We report that the geometrical connectivity of fracture networks plays a critical role in the hydromechanical processes in fractured rocks. A well-connected fracture system under a high stress ratio condition exhibits intense frictional sliding and large fracture dilation/opening, leading to greater rock mass permeability. However, a disconnected fracture network accommodates much less fracture shearing and opening, and has much lower bulk permeability. We further propose an analytical solution for the relationship between the equivalent permeability of fractured rocks and the connectivity metric (i.e. percolation parameter) of fracture networks, which yields an excellent match to the numerical results. We infer that fluid flow through a well-connected system is governed by traversing channels (forming an “in parallel” architecture) and thus equivalent permeability is sensitive to stress loading (due to stress-dependent fracture permeability), whilst fluid flow through a disconnected system is more ruled by matrix (linking isolated clusters “in series”) and has much less stress dependency.

This paper presents modeling of mechanical anisotropy in argillaceous rocks using an irregular lattice modeling approach, namely the rigid-body-spring network. To represent the mechanical anisotropy, new schemes are implemented in the modeling framework. The directionality of elastic deformation is resolved by modifying the element formulation with anisotropic elastic properties. The anisotropy of strength and failure characteristics is facilitated by adopting orientation-dependent failure criteria into the failure model. The verification of the improved modeling procedures is performed against theoretical model predictions for unconfined compression tests with various bedding orientations. Furthermore, excavation damage and fracturing processes in rock formations are simulated for different geomechanical configurations, such as rock anisotropy and tectonic heterogeneity. The simulated excavation damage characteristics are realistic and comparable with the actual field observation at a tunnel located in an argillaceous clay formation. The simulation results provide insights into the excavation damage zone phenomena with an explicit representation of fracturing processes.

The hydration of a bentonite barrier in the early stage of a geologic nuclear waste repository with a bentonite buffer is a critical issue for its long-term performance and safety because bentonite might be permanently altered and subsequently affect the function of bentonite barrier. Large scale in situ testing integrated with modeling analysis is an effective way to study the key processes affecting the hydration of a bentonite barrier. In this paper, through the comparison between coupled thermal, hydrological, mechanical, and chemical (THMC) models and data from a long term in situ test, we attempt to pinpoint the importance of non-Darcian flow, thermal osmosis, and hydro-mechanical coupling (porosity and permeability change due to swelling) to the hydration rate of the bentonite barrier under heating conditions. We found that a TH model equipped with non-Darcian flow severely underestimates the relative humidity and water content measured in the bentonite. Calibration of the parameters associated with relative permeability overshadows the contribution of non-Darcian flow, and non-Darcian flow under unsaturated conditions is not yet fully understood. An empirical relationship between saturated permeability and dry density was found to work better than a saturated permeability that is the function of effective stress in matching the relative humidity, water content data, and the chloride concentration in pore water. We also found that chemical data are actually helpful in calibrating the THM model. A question regarding the relevance of thermal osmosis to the hydration process, in terms of matching models and data, remains unanswered. Although a THMC model with thermal osmosis matches all THMC data nicely, similar goodness-of-fit can also be achieved by a THMC model without thermal osmosis but with lower permeability. We learned that the robustness of the model could be increased if the model is tested against long-term data and multiple types of data, and given that non-uniqueness is inevitable, more independent measurements of key parameters and multi-scale and multi-physics tests may help approximate the right model for evaluating the safety of the repository.

Establishing the geomechanical stability of marine sediments in the vicinity of a production well is one of the key design considerations in planning offshore gas production from marine hydrate reservoirs. This paper presents an assessment of the sediment stability at India's National Gas Hydrate Program, Expedition 2 (NGHP-02) Site 16 Area B offshore eastern India, for which gas production is to be carried out by depressurization. One important feature of the study is that extensive calibration of constitutive model parameters has been conducted based on laboratory test data from pressured core samples. From analysis perspective, the site is challenging because the hydrate reservoir consists of thin layers of hydrate-bearing sands interbedded with mud. Moreover, depressurization at the depth of a reservoir more than 2750 m below sea surface will lead to a pore pressure drop, and accordingly an effective confining stress increase as high as 25 MPa. In dealing with thin interbedded hydrate-bearing strata, meshing requirements for flow and geomechanical analysis are quite different from those for reservoirs with thicker massive layers, An axisymmetric model and one-way coupling simulations were thus adopted for this study, in which the geomechanical study utilizes pore pressure and hydrate saturation output from the flow study, but the flow study does not takes the porosity changes from the geomechanical analysis. Instead, the reduction of porosity due to sediment deformation in the flow study is based on a pressure-dependent pore compressibility relationship derived from geomechanical modeling. The rationality is validated through back computing the pore compressibility from the geomechanical deformation results. The study shows that large compression in the reservoir will result in movement of the sediments from above and below, as well as laterally in smaller magnitudes; and the sediment is deemed stable during the gas production period.

The increase of temperature of a low-permeability, fluid-saturated media may trigger significant thermal pressurization, where the expansion of pore fluid cannot be accommodated by the thermal expansion of the pore space. With the aid of a coupled thermohydromechanical numerical simulator, we investigate the possible impact of thermal pressurization during the life of a deep geological repository (DGR) for high-level radioactive waste, characterized by the emplacement of radioactive material-filled canisters in a series of parallel tunnels, excavated in a low-permeability clay formation. We represent the fault as a planar structure embedded in a thermoporoelastic material and shear activation evaluated by a strain-softening Mohr-Coulomb failure criterion. The results show that stress changes caused by temperature and thermal pressurization of a rock mass around the emplacement tunnels may trigger a slip event on a fault plane in proximity of the geological disposal site: Rupture nucleates at depth, hundreds of meters below the DGR. Stress transfer plays a key role, while a direct hydraulic connection between the repository and the nucleation zone is not necessary in order to trigger rupture. A low stress ratio may favor the occurrence of slip up to a distance of 600 m of the fault from the outermost tunnel. These results highlight the need of investigating hydromechanical properties and local stress conditions at depth to characterize the geomechanical response of weak planes located in the surroundings of a high-level radioactive waste repository and to provide sufficient knowledge for the safe development of the DGR site.

Geologic carbon storage, as well as other geo-energy applications, such as geothermal energy, seasonal natural gas storage and subsurface energy storage imply fluid injection and/or extraction that causes changes in rock stress field and may induce (micro)seismicity. If felt, seismicity has a negative effect on public perception and may jeopardize wellbore stability and damage infrastructure. Thus, induced earthquakes should be minimized to successfully deploy geo-energies. However, numerous processes may trigger induced seismicity, which contribute to making it complex and translates into a limited forecast ability of current predictive models. We review the triggering mechanisms of induced seismicity. Specifically, we analyze (1) the impact of pore pressure evolution and the effect that properties of the injected fluid have on fracture and/or fault stability; (2) non-isothermal effects caused by the fact that the injected fluid usually reaches the injection formation at a lower temperature than that of the rock, inducing rock contraction, thermal stress reduction and stress redistribution around the cooled region; (3) local stress changes induced when low-permeability faults cross the injection formation, which may reduce their stability and eventually cause fault reactivation; (4) stress transfer caused by seismic or aseismic slip; and (5) geochemical effects, which may be especially relevant in carbonate-containing formations. We also review characterization techniques developed by the authors to reduce the uncertainty in rock properties and subsurface heterogeneity both for the screening of injection sites and for the operation of projects. Based on the review, we propose a methodology based on proper site characterization, monitoring and pressure management to minimize induced seismicity.

Faults in the host rock that might exist in the vicinity of deep geological repositories for radioactive waste, constitute potential enhanced pathways for radionuclide migration. Several processes might trigger pore pressure increases in the faults leading to fault failure and induced seismicity, and increase the faults’ permeability. In this research, we developed a mathematical model to simulate fault activation during an experiment of controlled water injection in a fault at the Mont-Terri Underground Research Laboratory in Switzerland. The effects of in-situ stress, fault shear strength parameters and heterogeneity are assessed. It was shown that the above factors are critical and need to be adequately characterized in order to predict the faults’ hydro-mechanical behaviour.

A numerical simulator entitled TOUGH-UDEC is introduced for the analysis of coupled thermal-hydraulic-mechanical processes in fractured porous media. Two existing well-established codes, TOUGH2 and UDEC, are coupled to model multiphase fluid flows, heat transfers, and discontinuous deformations in fractured porous media by means of discrete fracture representation. TOUGH2 is widely used for the modeling of heat transfers and multiphase multicomponent fluid flows, and UDEC is a well-known distinct element code for rock mechanics. The two codes are solved sequentially, with coupling parameters passed to each equation at specific intervals. After solving thermal-hydraulic equations within the TOUGH2 code, pressure and temperature information is imported into the UDEC model. After solving the mechanical equation within the UDEC code the calculated fracture apertures are converted to the equivalent permeability and porosity values for a TOUGH2 flow analysis. The solution is calculated by iteratively following an explicit sequence for numerical efficiency. Verifications are presented to demonstrate the capabilities of the coupled TOUGH-UDEC simulator. Three application examples of (1) shear dilation due to increased pore pressure, (2) thermal stress and (3) CO 2 injection, show that the new simulator can be an effective tool for geoengineering applications involving shear activation of fractures and faults.

Free gas migration through the gas hydrate stability zone (GHSZ) and subsequent gas seepage at the seabed are characteristic features in marine gas hydrate provinces worldwide. The biogenic or thermogenic gas is typically transported along faults from deeper sediment strata to the GHSZ. Several mechanisms have been proposed to explain free gas transport through the GHSZ. While inhibition of hydrate formation by elevated salinities and temperatures have been addressed previously in studies simulating unfocused, area-wide upward advection of gas, which is not adequately supported by field observations, the role of focused gas flow through chimney-like structures has been underappreciated in this context. Our simulations suggest that gas migration through the GHSZ is, fundamentally, a result of methane gas supply in excess of its consumption by hydrate formation. The required high gas flux is driven by local overpressure, built up from gas accumulating below the base of the GHSZ that fractures the overburden when exceeding a critical pressure, thereby creating the chimney-like migration pathway. Initially rapid hydrate formation raises the temperature in the chimney structure, thereby facilitating further gas transport through the GHSZ. As a consequence, high hydrate saturations form preferentially close to the seafloor, where temperatures drop to bottom water values, producing a prominent subsurface salinity peak. Over time, hydrates form at a lower rate throughout the chimney structure, while initial temperature elevation and salinity peak dissipate. Thus, our simulations suggest that the near-surface salinity peak and elevated temperatures are a result of transient high-flux gas migration through the GHSZ.

In this paper, a new numerical simulation tool named TOUGH-FEMM is presented and applied to model three-dimensional (3D) hydraulic fracturing in porous rock. The fluid flow in both fractures and porous rock is modeled using TOUGH2, which is a well-established code for analysis of multiphase and multi-component fluid flow. Rock deformations associated with fracture propagation are modeled using finite element-meshfree method (FEMM). FEMM is an approach to simulate fracture propagation without remeshing, in which the fracture path does not need to be predetermined. Fracture mechanics with mixed-mode stress intensity factors are employed to detect fracture instability and determine the direction of fracture propagation. TOUGH-FEMM is verified for modeling fluid-driven fracture propagation in 3D through a number of simulation examples, including modeling of hydraulic fracturing laboratory experiments and by comparison to independent numerical simulation results for multiple interacting hydraulic fractures at ten to hundred meter scale.

The center of the San Joaquin Valley is one of the most productive agricultural regions in the US. Farmers rely heavily on surface-water diversions to meet irrigation water demand. However, the 2007–2010 and 2012–2017 droughts have caused strong increases in groundwater pumping causing land subsidence with strong variability in location, magnitude (total subsidence) and rate of subsidence. In this study, we try to understand what caused these variations. We focus our analyses on three areas: (i) the Westland water district where subsidence was very small during the two drought periods, (ii) ‘El Nido’ area where the greatest subsidence rate was monitored from 2008 to 2010, and (iii) the’ Kings-Tulare counties’ area where the subsidence was small during the 2008–2010 and the largest during 2015–2017 droughts. Our main finding is that land subsidence is located in areas where the water demand for agriculture and the density of groundwater wells is the highest, whereas the rate of subsidence is strongly affected the amount of local and imported surface water and by groundwater resources. Based on these simple observations, we propose using continuous satellite-based ground deformation monitoring and geomechanical modeling to (i) localize areas most prone for future subsidence and (ii) to estimate and manage groundwater resources.

Numerous analytical solutions have been developed for modeling thermal perturbations to underground formations caused by deep-well injection of fluids. Each solution has been derived for a specific boundary value problem and a simplified flow network with one set of parallel fractures. In this paper, new generalized solutions G*(x, s) are developed using (existing) global transfer functionsG*0 (x,s) and a new memory function g*(s), where x and s are the space and Laplace variable. The memory function represents the solutions of conductive heat exchange between fractures and matrix blocks and between fractured aquifers and unfractured aquitards. The memory function is developed to account for multirate exchange induced by different shapes, sizes, properties, and volumetric fractions of matrix blocks bounded by multiple sets of orthogonal fractures with different spacing. The global transfer functions represent the fundamental solutions to convective, convective-conductive, and convective-dispersive heat transport in fractures (or aquifers) without exchange and are available for various (1-D linear, 1-D radial, 2-D dipole, and single-well injection-withdrawal) flow fields. The new solutionswith exchange are developed usingG*(x,s) = B*(s)G*0{x, s[1 + Θg*(s)]}, thereby greatly simplifying solution development in a novel way, where Θ and B*(s) are a fracture-matrix scaling factor and the boundary condition function. The new solutions are applied to several example problems, showing that heat transport in fractured aquifers is significantly impacted by (1) thermal dispersion in fractures that is rarely considered, (2) multirate heat exchange with a wide range of size and anisotropy of rectangular matrix blocks, and (3) heat exchange between aquifers and aquitards.

In this study, we analyze a deep hydraulic stimulation of a fracture zone that was conducted as part of the classical Fenton Hill Hot Dry Rock Program in the 1970s. At the time, it was suggested that a pre-existing fracture or multiple fractures within the fracture zone were jacked open by injection-induced increase in pressure. In this study, we analyze the same stimulation experiment, but investigate the possibility of an alternative mechanism of shear reactivation of pre-existing fractures. We conduct modeling that accounts for both jacking (or elastic-fracture opening) and shear-slip dilation and demonstrate that injection-induced shear reactivation (or hydroshearing) could have occurred simultaneously with seismic events of magnitudes lower than what can be felt by humans. In fact, simulations considering shear reactivation seem to better match observed fluid recovery after multiple injection cycles. Shear reactivation and shear dilation results in locked-open fractures, especially near the injection well that provides permeability of higher flow recovery. We then investigate the sensitivity of the proposed model by varying some of the critical parameters such as maximum aperture, dilation angle, as well as fracture density. Interestingly, none of the simulated cases resulted in a large event that could have been felt by humans, but did result in a cumulative seismic magnitude of less than 1 for each given stimulation step. These results suggest that a permanent irreversible permeability increase of several orders of magnitude can be obtained by hydroshearing in a “seismically” safe manner.

This paper presents numerical modeling of coupled thermal, hydraulic and mechanical processes in rock salt and crushed salt considering halite solubility constraints. The TOUGH-FLAC simulator is used, with a recently enhanced Equation-Of-State module that includes the thermodynamic properties of aqueous fluids of variable salinity. Laboratory and field scale tests performed on rock salt and crushed salt under temperature gradients are modeled first to evaluate the capabilities of the simulator to reproduce important features, such as porosity changes induced by halite dissolution/precipitation, and brine and heat migration. Since the results are quite satisfactory, the simulator is used to predict the long-term response of a generic salt repository for heat-generating nuclear waste. To evaluate the impacts of halite solubility on the predictions, two simulations that respectively consider or neglect solubility constraints are performed. In the scenario studied, the results are not significantly affected by dissolution/precipitation, and only some differences are observed due to changes in porosity, but the dominating processes remain the same. With the new provisions, TOUGH-FLAC is more complete in terms of processes occurring around a heat-releasing nuclear waste package and can therefore provide more accurate predictions of the long-term performance of a nuclear waste repository in salt formations.

In the paper already published, there is an error in the acknowledgements. The grant number concerning EU project is 636811 instead of 640979. There was a misunderstanding with the EU project number and we only found that after those papers are published. We confirm that the correction is to be made online and we also want to publish a corrigendum. The corrected text should be: The authors regret . “The authors gratefully acknowledge the Swedish Geological Survey (SGU), grant number 1724, and the EU project Fracrisk grant number 636811, for providing financial support to research reported in this paper. Additional support was provided by the U.S. Department of Energy under contract No. DE-AC02-05CH11231.” The authors would like to apologise for any inconvenience caused.

The authors regret the error regarding the grant number in the acknowledgements in the paper already published. The grant number concerning EU project is 636811 instead of 640979. There was a misunderstanding with the EU project number and we only found that after those papers are published. The corrected text should be: “The authors gratefully acknowledge the Swedish Geological Survey (SGU), grant number 1724, and the EU project Fracrisk grant number 636811, for providing financial support to research reported in this paper. Additional support was provided by the U.S. Department of Energy under contract No. DE-AC02-05CH11231.” The authors would like to apologise for any inconvenience caused.

For the deep geological disposal of high-level radioactive waste in argillaceous rocks, the heat production of the waste is an important driver for thermal–hydraulic-mechanical (THM)-coupled processes. These THM processes influence the properties and conditions of the near field that in many repository designs contains bentonite as a clay buffer. One task in the DECOVALEX-2015 (DEvelopment of COupled models and their VALidation against Experiments) project was the modelling of a heated bentonite column (Villar et al. in Long-term THM tests reports: THM cells for the HE-E test: update of results until February 2014. Deliverable-no: D2.2-7.3. CIEMAT Technical Report IEMAT/DMA/2G210/03/2014, 2014) in preparation for the in situ heater experiment HE-E at the underground rock laboratory Mont Terri. DECOVALEX is an international cooperative project that focuses on the development and validation of mathematical models for simulating such coupled processes associated with disposal in deep geological repositories. Eight modelling teams developed their own THM-coupled models for the bentonite column experiment, using six different simulation codes. Each of the teams individually calibrated the THM parameters for the bentonite material. The eight resulting parameter sets agree well and allow a satisfactory reproduction of the TH measurements by all models. The modelling results for the evolution of temperature and relative humidity over time at three sensors in the bentonite column are in good agreement between the teams and with the measured data. Also, changes of the temperature due to modifications of the insulation and the adjustment of the heating power during the course of the experiment are well reproduced. The models were thus able to reproduce the main physical processes of the experiment, both for vapour-dominated diffusion during the heating phase and combined liquid and vapour transport during a subsequent heating and hydration phase. Based on the parameter sets, the teams predict a penetration of the water infiltration front in the 48-cm column filled with bentonite pellets to a depth between 25 and 35 cm over the 15,000 h (i.e. over 20 months) of the hydration phase of the experiment.

We studied the relation between rupture and changes in permeability within a fault zone intersecting the Opalinus Clay formation at 300 m depth in the Mont Terri Underground Research Laboratory (Switzerland). A series of water injection experiments were performed in a borehole straddle interval set within the damage zone of the main fault. A three-component displacement sensor allowed an estimation of the displacement of a minor fault plane reactivated during a succession of step rate pressure tests. The experiment reveals that the fault hydromechanical (HM) behavior is different from one test to the other with varying pressure levels needed to trigger rupture and different slip behavior under similar pressure conditions. Numerical simulations were performed to better understand the reason for such different behavior and to investigate the relation between rupture nucleation, permeability change, pressure diffusion, and rupture propagation. Our main findings are as follows: (i) a rate frictional law and a rate-and-state permeability law can reproduce the first test, but it appears that the rate constitutive parameters must be pressure dependent to reproduce the complex HM behavior observed during the successive injection tests; (ii) almost similar ruptures can create or destroy the fluid diffusion pathways; (iii) a too high or too low diffusivity created by the main rupture prevents secondary rupture events from occurring whereas “intermediate” diffusivity favors the nucleation of a secondary rupture associated with the fluid diffusion. However, because rupture may in certain cases destroy permeability, this succession of ruptures may not necessarily create a continuous hydraulic pathway.

Subsurface manipulations such as those expected from the disposal of heat-emanating radioactive waste in deep repositories can induce strongly coupled Thermal (T), hydrological (H), mechanical (M) and chemical (C) processes. Adequate coupled THMC models are highly desirable or even indispensable for performance assessment of such repositories, for examples for the analysis of bentonite-based engineered barrier system (EBS) surrounding the emplaced waste. In this study, we present coupled THMC model simulations of a generic nuclear waste repository in an argillite with a bentonite-based buffer. The objective is to evaluate the chemical changes in the EBS bentonite and their effect on mechanical behavior under high temperature, attempting to shed light on whether EBS bentonite can sustain temperatures higher than 100 °C without significant impact on barrier's performance. Two scenarios were simulated for comparison: a case in which the temperature near the waste canister peaks at 200 °C and a case in which the temperature at the same spot culminate with about 100 °C. Simulations for a generic case with Kunigel-VI bentonite as backfill and Opalinus Clay as host rock were conducted for 1000 years and reported in the previous study (Zheng et al., 2015). In this paper, simulations for 100,000 years have been done for two types of bentonite-based buffer materials: Kunigel-VI and FEBEX bentonite. This enables us to evaluate how different types of bentonite behave in terms of the illitization and its impact on swelling stress and whether we can generalize these results to support decision making. The simulations show the occurrence of illitization in the bentonite buffer and the enhancement of illitization under high temperature. However, FEBEX bentonite undergoes less illitization mainly due to the higher ion concentration in pore water and the lower content of K-feldspar in the bentonite mineral composition. Moreover, the reduction of swelling stress by chemical changes is more pronounced for Kunigel-VI bentonite than for FEBEX bentonite. Overall, the results of our model simulations suggest that an argillite repository with a bentonite-based EBS that is similar to FEBEX bentonite could sustain temperatures much higher than 100 °C as far as illitization concerns. Model results also reveal that illitization is stabilized after about 2000 years in bentonite near the waste package, but continues in bentonite near the bentonite-argillite interface, which manifests the strong effect of geochemical interaction between EBS bentonite and host rock on long term illitization in bentonite.

There are two types of analytical solutions of temperature/concentration in and heat/mass transfer through boundaries of regularly shaped 1-D, 2-D, and 3-D blocks. These infinite-series solutions with either error functions or exponentials exhibit highly irregular but complementary convergence at different dimensionless times, td. In this paper, approximate solutions were developed by combining the error-function-series solutions for early times and the exponential-series solutions for late times and by using time partitioning at the switchover time, td0. The combined solutions contain either the leading term of both series for normal-accuracy approximations (with less than 0.003 relative error) or the first two terms for high-accuracy approximations (with less than 10−7 relative error) for 1-D isotropic (spheres, cylinders, slabs) and 2-D/3-D rectangular blocks (squares, cubes, rectangles, and rectangular parallelepipeds). This rapid and uniform convergence for rectangular blocks was achieved by employing the same time partitioning with individual dimensionless times for different directions and the product of their combined 1-D slab solutions. The switchover dimensionless time was determined to minimize the maximum approximation errors. Furthermore, the analytical solutions of first-order heat/mass flux for 2-D/3-D rectangular blocks were derived for normal-accuracy approximations. These flux equations contain the early-time solution with a three-term polynomial in √td and the late-time solution with the limited-term exponentials for rectangular blocks. The heat/mass flux equations and the combined temperature/concentration solutions form the ultimate kernel for fast simulations of multirate and multidimensional heat/mass transfer in porous/fractured media with millions of low-permeability blocks of varying shapes and sizes.