Abstract: One of the challenging problems in mathematical geosciences is the determination of analytical solutions of nonlinear partial differential equations describing transport processes in porous media. We are interested in diffusive transport coupled with precipitation-dissolution reactions. Several numerical computer codes that simulate such systems have been developed. Analytical solutions, if they exist, represent an important tool for verification of numerical solutions. We present a methodology for deriving such analytical solutions that are exact and explicit in space and time variables. They describe transport of several aqueous species coupled to precipitation and dissolution of a single mineral in one, two, and three dimensions. As an application, we consider explicit analytical solutions for systems containing one or two solute species that describe the evolution of solutes and solid concentrations as well as porosity. We use one of the proposed analytical solutions to test numerical solutions obtained from two conceptually different reactive transport codes. Both numerical implementations could be verified with the help of the analytical solutions and show good agreement in terms of spatial and temporal evolution of concentrations and porosities.
Abstract: Diffusion of cations and other contaminants through clays is of central interest, because clays and clay rocks are widely considered as barrier materials for waste disposal sites. An intriguing experimental observation has been made in this context: Often, the diffusive flux of cations at trace concentrations is much larger and the retardation smaller than expected based on their sorption coefficients. So-called surface diffusion of sorbed cations has been invoked to explain the observations but remains a controversial issue. Moreover, the corresponding surface diffusion coefficients are largely unknown. Here we show that, by an appiopriate scaling, published diffusion data covering a broad range of cations, clays, and chemical conditions can all be modeled satisfactorily by a surface diffusion model. The average mobility of sorbed cations seems to be primarily an intrinsic property of each cation that follows inversely its sorption affinity. With these surface mobilities, cation diffusion coefficients can now be estimated from those of water tracers. In pure clays at low salinities, surface diffusion can reduce the cation retardation by a factor of more than 1000.
Abstract: We present exact analytical solutions for a one-dimensional diffusion problem coupled with the precipitation-dissolution reaction A((aq)) vertical bar B((aq)) reversible arrow M((s)) and feedback of porosity change. The solutions are obtained in the form of traveling waves and describe spatial and temporal evolutions of solute concentration, porosity, and mineral distribution for a set of initial and boundary conditions. The form of the solutions limits the choice of admissible boundary conditions, which might be difficult to adapt in natural systems, and thus, the solutions are of limited use for such a system. The main application of the derived solutions is therefore the benchmarking of numerical reactive transport codes for systems with strong porosity change. To test the performance of numerical codes, numerical solutions obtained by using a global implicit finite volume technique are compared to the analytical solutions. Good agreement is obtained between the analytical solutions and the numerical solutions when a sufficient spatial discretization resolves the spatial concentration gradients at any time. In the limit of fast kinetics (local equilibrium), steep concentration fronts cannot be resolved in a numerical discretization schema.
Abstract: The Car-Parrinello molecular dynamics simulation technique was used to predict the structure and dynamics of hydronium solvation in mono-, bi- and trihydrated Na-montmorillonite. In monohydrated montmorillonite, hydronium ions are located within the hexagonal rings of the basal clay plane. Oxygen sites of hydronium ions point towards the clay surface and hydrogen atoms towards the water layer. In bi- and trihydrated montmorillonite, hydronium ions form water-solvated, outer-sphere complexes. Similar to the solvation mechanism in bulk water, hydronium ions donate three hydrogen bonds to interlayer water molecules. In all studied hydration states, hydronium ions do not form hydrogen bonds with the basal oxygen sites. Similar to bulk water, the free energy barrier for a classical proton transfer between interlayer water molecules is of the order of kT and therefore not the limiting factor for the proton diffusion. The diffusivity of hydrogen in the interlayer is controlled by the structural rearrangements of the solvating water molecules.
Abstract: Clays and clay rocks are considered viable geotechnical barriers in radioactive waste disposal. One reason for this is the propensity for cation exchange reactions in clay minerals to retard the migration of radionuclides. Although another retardation mechanism, namely the incorporation of radionuclides into sulfate or carbonate solid solutions, has been known for a long time, only recently has it been examined systematically. In this work, we investigate the competitive effect of both mechanisms on the transport of radium (Ra) in the near-field of a low- and intermediate level nuclear waste repository. In our idealized geochemical model, numerical simulations show that barium (Ba) and strontium (Sr) needed for Ra Sulfate solid solutions also partition to the cation exchange sites of montmorillonite (Mont), which is the major mineral constituent of bentonite that is used for tunnel backfill. At high Mont content, most Ra tends to attach to Mont, while incorporation of Ra in sulfate solid solutions is more important at low Mont content. To explore the effect of the Mont content on the transport of radium, a multi-component reactive transport model was developed and implemented in the scientific software OpenGeoSys-GEM. It was found that a decrease of fixation capacity due to low Mont content is compensated by the formation of solid solutions and that the migration distance of aqueous Ra is similar at different Mont/water ratios.
Abstract: The numerical simulation of reactive mass transport processes in complex geochemical environments is an important tool for the performance assessment of future waste repositories. A new combination of the multi-component mass transport code GeoSys/RockFlow and the Gibbs Energy Minimization (GEM) equilibrium solver GEM-Selektor is used to calculate the accurate equilibrium of multiple non-ideal solid solutions which are important for the immobilization of radionuclides such as Ra. The coupled code is verified by a widely used benchmark of dissolution-precipitation in a calcite-dolomite system. A more complex application shown in this paper is the transport of Ra in the near-field of a nuclear waste repository. Depending on the initial inventories of Sr, Ba and sulfate, non-ideal sulfate and carbonate solid solutions can fix mobile Ra cations. Due to the complex geochemical interactions, the reactive transport simulations can describe the migration of Ra in a much more realistic way than using the traditional linear K(D) approach only. (C) 2009 Elsevier Ltd. All rights reserved.
Abstract: Many applied problems in geoscience require knowledge about complex interactions between multiple physical and chemical processes in the sub-surface. As a direct experimental investigation is often not possible, numerical simulation is a common approach. The numerical analysis of coupled thermo-hydro-mechanical (THM) problems is computationally very expensive, and therefore the applicability of existing codes is still limited to simplified problems. In this paper we present a novel implementation of a parallel finite element method (FEM) for the numerical analysis of coupled THM problems in porous media. The computational task of the FEM is partitioned into sub-tasks by a priori domain decomposition. The sub-tasks are assigned to the CPU nodes concurrently. Parallelization is achieved by simultaneously establishing the sub-domain mesh topology, synchronously assembling linear equation systems in sub-domains and obtaining the overall solution with a sub-domain linear solver (parallel BiCGStab method with Jacobi pre-conditioner). The present parallelization method is implemented in an object-oriented way using MPI for inter-processor communication. The parallel code was successfully tested with a 2-D example from the international DECOVALEX benchmarking project. The achieved speed-up for a 3-D extension of the test example on different computers demonstrates the advantage of the present parallel scheme. (C) 2009 Elsevier Ltd. All rights reserved.
Abstract: Quantification of mass and heat transport in fractured porous rocks is important to areas such as contaminant transport, storage and release in fractured rock aquifers, the migration and sorption of radioactive nuclides from waste depositories, and the characterization of engineered heat exchangers in the context of enhanced geothermal systems. The large difference between flow and transport characteristics in fractures and in the surrounding matrix rock means models of such systems are forced to make a number of simplifications. Analytical approaches assume a homogeneous system, numerical approaches address the scale at which a process is operating, but may lose individual important processes due to averaging considerations. Numerical stability criteria limit the contrasts possible in defining material properties. Here, a hybrid analytical-numerical method for transport modeling in fractured media is presented. This method combines a numerical model for flow and transport in a heterogeneous fracture and an analytical solution for matrix diffusion. By linking the two types of model, the advantages of both methods can be combined. The methodology as well as the mathematical background are developed, verified for simple geometries, and applied to fractures representing experimental field conditions in the Grimsel rock laboratory.
Abstract: The state and dynamics of water and cations in pure and mixed Na-Cs-montmorillonite as a function of the interlayer water content were investigated in the present study, using Monte Carlo and classical, molecular-dynamics methods. While highly idealized, the simulations showed that the swelling behavior of hetero-ionic Na-Cs-montmorillonite is comparable to the swelling of a homo-ionic Na- or Cs-montmorillonite. The mixed Na-Cs-montmorillonite is characterized by intermediate interlayer distances compared to horrio-ionic Na- and Cs-montmorillonites. Dry, hetero-ionic Na-Cs-montmorillonite is characterized by a symmetric sheet configuration, as is homo-ionic Cs-montmorillonite. We found that at low degrees of hydration the absolute diffusion coefficient of Cs+ is less than for Na+, whereas at greater hydration states the diffusion coefficient of Cs+ is greater than for Na+. An analysis of the relative diffusion coefficients (the ratio between the diffusion coefficient of an ion in the interlayer and its diffusion coefficient in bulk water) revealed that water and Na+ are always less retarded than Cs+. With large interlayer water contents, tetralayer or more, Na+ ions preferentially form outer-sphere complexes. The mobility perpendicular to the clay surface is limited and the diffusion is equivalent to two-dimensional diffusion in bulk water. In contrast, Cs+ ions preferentially form âinner-sphere complexesâ at all hydration states and their two-dimensional diffusion coefficient is less than in bulk water. The question remains unanswered as to why experimentally derived relative diffusion coefficients of Cs+ in the interlayer of clays are about 20 times less than those we obtained by classical molecular dynamics studies.
Abstract: Discrete fracture network simulations are computationally intensive and usually time-consuming to construct and configure. This paper presents a case study with techniques for building a 3D finite element model of an inhomogeneous fracture network for modelling flow and tracer transport, combining deterministic and stochastic information on fracture aperture distributions. The complex intersected fractures represent a challenge for geometrical model design, mesh quality requirements and property allocations. For the integrated and holistic modelling approach, including the application of numerical and analytical simulation techniques, new object-oriented concepts in software engineering are implemented to ensure a resourceful and practicable software environment.
Abstract: We present a sequence of purely advective transport models that demonstrate the influence of small-scale geometric inhomogeneities on contaminant transport in fractured crystalline rock. Special weight is placed on the role of statistically generated variable fracture apertures. The fracture network geometry and the aperture distribution are based on information from an in situ radionuclide retardation experiment performed at Grimsel test site (Swiss Alps). The obtained breakthrough curves are fitted with the advection dispersion equation and continuous-time random walks (CTRW). CTRW is found to provide superior fits to the late-arrival tailing and is also found to show a good correlation with the velocity distributions obtained from the hydraulic models. The impact of small-scale heterogeneities, both in fracture geometry and aperture, on transport is shown to be considerable.
Abstract: Transport experiments with colloids and radionuclides in a shear zone were conducted during the Colloid and Radionuclide Retardation experiment (CRR) at Nagraâs Grimsel Test Site. Breakthrough curves of bentonite colloids and uranine, a non-sorbing solute, were measured in an asymmetric dipole flow field. The colloid breakthrough is earlier than that of uranine. Both breakthrough curves show anomalously long late time tails and the slope of the late time tails for the colloids is slightly higher. Anomalous late time tails are commonly associated with matrix diffusion processes; the diffusive interaction of solutes transported in open channels with the adjacent porous rock matrix or zones of stagnant water. The breakthrough curves for different colloid size classes are very similar and show no signs of fractionation due to their (size-dependent) diffusivity. It is proposed that tailing of the colloids is mainly caused by the structure of the flow field and that for the colloid transport, matrix diffusion is of minor importance. This has consequences for the interpretation of the uranine breakthrough. Comparisons of experimental results with numerical studies and with the evaluation of the colloid breakthrough with continuous time random theory imply that the tailing in the conservative solute breakthrough in this shear zone is not only caused by matrix diffusion. Part of the tailing can be attributed to advective transport in fracture networks and advection in low velocity regions. Models based on the advection-dispersion equation and matrix diffusion do not properly describe the temporal and spatial evolution of colloid and solute transport in such systems with a consistent set of parameters. (C) 2003 Elsevier B.V. All rights reserved.
Abstract: Flow and transport patterns in three-dimensional fracture intersections are investigated, with emphasis on the occurrence and effects of local flows around each intersection. Flow and transport simulations indicate that local flow circulations (referred to as âlocal flow cellsâ) arise because of an interplay between intersection geometry and corresponding boundary conditions for flow along the fracture boundaries. Local flow cells are thus unique features in three-dimensional fracture networks. The analysis here suggests that local flow cells have little effect on the overall fluid flow in a fracture network but may constitute an effective mechanism which contributes to the long breakthrough tailing and retardation often observed in transport through discrete fracture networks.
Abstract: [1] A combined experimental and numerical analysis is used to examine the movement and interaction of saltwater and freshwater in unconfined porous media. Particular emphasis is placed on flow phenomena in the partially saturated region above the water table and along the interface between the saturated and partially saturated regions. While some numerical models are apparently capable of simulating these phenomena, there is still a significant lack of experimental data with which to verify the models. Here a series of laboratory-scale experiments is considered to evaluate density-dependent, saltwater-freshwater flow patterns in both the saturated and partially saturated zones. The laboratory experiments demonstrate clearly that significant lateral flows and coupled density-driven flow effects may take place in the partially saturated region above the water table and at the interface between the saturated and partially saturated zones. In parallel, a finite element numerical model is developed. The model reproduces effectively the observed flow behaviors; the quality of the results suggests that the numerical model has the capacity to provide realistic predictions.
Abstract: We examine a set of analytical solutions based on the continuous time random walk (CTRW) approach, which can be evaluated numerically and used to analyze breakthrough data from tracer tests. Practical application of these solutions, with discussion of the physical meaning of the relevant model parameters, is emphasized. The CTRW theory accounts for the often observed non-Fickian (or scale-dependent) dispersion behavior that cannot be properly quantified by using the advection-dispersion equation. The solutions given here, valid for a wide range of dispersive behaviors of conservative tracers, and useful for both characterization and prediction, have been integrated into a library of external functions for use with the GRACE graphical display acid analysis package, Example applications of these solutions are presented. The library and graphics software are freely accessible from a Web site.
Abstract: We analyze a set of observations from a recently published, field-scale tracer test in a fractured till. These observations demonstrate a dominant, underlying non-Fickian behavior, which cannot be quantified using traditional modeling approaches. We use a continuous time random walk (CTRW) approach which thoroughly accounts for the measurements, and which is based on a physical picture of contaminant motion that is consistent with the geometric and hydraulic characterization of the fractured formation. We also incorporate convolution techniques in the CTRW theory, to consider transport between different regions containing distinct heterogeneity patterns. These results enhance the possibility that limitations in predicting non-Fickian modes of contaminant migration can be overcome. (C) 2001 Elsevier Science B.V. All rights: reserved.
Abstract: We use numerical simulations to examine the variability of flow patterns in representative fracture intersection geometries. In contrast to existing studies of perfectly orthogonal intersections, we demonstrate that more realistic geometries lead to a rich spectrum of flow patterns. Moreover, numerical solutions of the Navier-Stokes equations in these fracture intersections indicate that non-linear inertial effects become important for Reynolds numbers as tow as 1-100. Such Reynolds numbers often exist in naturally fractured formations, particularly in karst systems and in the vicinity of wells during pump tests.
Abstract: This study combines experimental work and numerical simulations to reconstruct the thermal history of the Frankenwald Transverse Zone, which was formed by a granitic intrusion into a fault zone. Illite crystallinity, vitrinite reflectance, and geobarometric investigations reveal a metamorphic and paleo-temperature anomaly associated with the granitic intrusion. Results of numerical simulations adequately explain paleo-temperatures in that area. In order to be able to obtain a quantitative comparison between numerical model results and paleo-temperature as observed in the field, we propose an empirical relationship between illite crystallinity and the maximum paleo-temperature based on Literature data of illite crystallinity and a combination of other temperature-dependent parameters like vitrinite reflectance, phase petrology and smectite-to-illite transformation. Application of this strategy to the Frankenwald Transverse Zone yields the following results: (1) The paleo-temperature anomaly can be explained by the cooling of a number of plutons which intruded into the center of the zone. No additional heat sources are required to explain the observed anomaly. (2) The diapiric shape of these plutons could be confirmed because, in contrast, dike-shaped bodies would produce much smaller paleo-thermal anomalies. (3) The resolution of paleo-temperatures obtained from the illite crystallinity data is not good enough to discriminate precisely between advective and conductive modes of heat transfer. According to our preferred model, conductive heat transport is more likely than fluid-driven advective heat transport. (C) 1999 Elsevier Science B.V. All rights reserved.
Notes: 29th General Assembly of the International-Association-of-Seismology-and-Physics-of-the-Earths-Interi or, THESSALONIKI, GREECE, AUG 18-28, 1997
Abstract: For radioactive waste disposal, it is planned to use bentonite in the engineered barrier system to prevent the radionuclides from leaking into the ambient environment. Ion exchange, surface complexation and solid solutions are important chemical mechanisms that contribute to the buffering capacity of bentonite. For safety assessment purposes traditionally only ion exchange and surface complexation process are transferred into a linear sorption concept. In this work, we use the reactive transport code GeoSysGEM to simulate the transport of Ra(2+) in a bentonite column, with both cation-exchange and solid solution formations explicitly considered. A chemical environment based on the FEBEX bentonite was employed, with mineral content and ion selectivity coefficients measured from experiments. Simulation results suggest that both solid solutions and cation-exchange reactions have a strong retardation effect on Ra(2+). Together they provide an enormous fixation capacity. However, most Ra(2+) tends to be included in the sulfate solid solutions formation, which finally controls the amount of Ra(2+) in solution. With this model, we are able to describe the transport of radionuclides in bentonite more precisely. This offers the possibility to investigate the effects of varying buffer material composition and optimize the retardation properties of the engineered barrier system.
Notes: 7th International Conference on Calibration and Reliability in Groundwater Modeling, Wuhan, PEOPLES R CHINA, SEP 20-23, 2009
Abstract: We analyze a set of observations from a recently published, field-scale tracer test in a fractured till. These observations demonstrate a dominant, underlying non-Gaussian behavior, which cannot be quantified using traditional modeling approaches. We use a continuous time random walk (CTRW) approach which thoroughly accounts for the measurements, and which is based on a physical picture of contaminant motion that is consistent with the geometric and hydraulic characterization of the fractured formation.
Notes: International Conference on Groundwater Research, COPENHAGEN, DENMARK, JUN 06-08, 2000
Abstract: In recent years environmental simulations receive increasing attention in concern e.g. with the safety of waste disposal or the utilization of geothermal energy. For this purpose mass and heat transfer must be modelled in hard rocks, which are complex, strongly discontinuous media consisting of numerous fracture systems at various scales. Different approaches exist for the descriptions of fractured media. So called 2 1/2-D models treat only fractures as intersecting planes in 3D space. interaction between fractures and rock matrix is ignored. However, recent studies emphasize the importance of fluid, mass and heat storage within the enclosed rock matrix, especially for long-term analysis. To consider the fracture-matrix system a fully 3-D model is required incorporating plane fractures as well as spatial matrix blocks. Fundamental methods for the finite element discretization into hexahedral elements of these heterogenous media will be presented. Furthermore, applications of 2 1/2-D models and 3-D mass and heat transfer models illustrate the capability of the proposed discretization technique and the finite element code Rockflow which is used for the numerical simulations.
Notes: 2nd ECCOMAS Conference on Numerical Methods in Engineering, PARIS, FRANCE, SEP 09-13, 1996
