<figure><img src="image_1.png" /></figure>

Check for updates

OPEN ACCESS

EDITED BY
Bernd Grambow,
UMR6457 Laboratoire de Physique
Subatomique et des Technologies Associées
(SUBATECH), France

REVIEWED BY
Lara Duro,
Amphos21, Spain
Mavrik Zavarin,
Lawrence Livermore National Laboratory
(DOE), United States

*CORRESPONDENCE

Norbert Maes,
✉ norbert.maes@sckcen.be
Sergey Churakov,
✉ sergey.churakov@psi.ch;

RECEIVED 15 April 2024
ACCEPTED 18 July 2024
PUBLISHED 22 November 2024

CITATION

Maes N, Churakov S, Glaus M, Baeyens B, Dähn R, Grangeon S, Charlet L, Brandt F, Poonoosamy J, Hoving A, Havlova V, Fischer C, Noseck U, Britz S, Siitari-Kauppi M, Li X, Fabritius O and Missana T (2024) EURAD state-of-the-art report on the understanding of radionuclide retention and transport in clay and crystalline rocks.

Front. Nucl. Eng. 3:1417827.
doi: 10.3389/fnuen.2024.1417827

COPYRIGHT

© 2024 Maes, Churakov, Glaus, Baeyens, Dähn, Grangeon, Charlet, Brandt, Poonoosamy, Hoving, Havlova, Fischer, Noseck, Britz, Siitari-Kauppi, Li, Fabritius and Missana. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.

# EURAD state-of-the-art report on the understanding of radionuclide retention and transport in clay and crystalline rocks

Norbert Maes¹*, Sergey Churakov²,³*, Martin Glaus²,
Bart Baeyens², Rainer Dähn², Sylvain Grangeon⁴,
Laurent Charlet⁵, Felix Brandt⁶, Jenna Poonoosamy⁶,
Alwina Hoving⁷, Vaclava Havlova⁸, Cornelius Fischer⁹,
Ulrich Noseck¹⁰, Susan Britz¹⁰, Marja Siitari-Kauppi¹¹,
Xiaodong Li¹⁰,¹¹, Otto Fabritius¹⁰,¹¹ and Tiziana Missana¹²

¹SCK CEN, Mol, Belgium, ²Paul Scherrer Instute (PSI), Villigen, Switzerland, ³Institute of Geological Sciences, University of Bern, Bern, Switzerland, ⁴BRGM, Orléans, France, ⁵Institut des Sciences de la Terre (ISTERRE), Grenoble, France, ⁶Forschungszentrum Jülich GmbH (FZJ), Jülich, Germany, ⁷Dutch Organization for Applied Scientific Research (TNO), Utrecht, Netherlands, ⁸Nuclear Research Institute Rez plc (UJV), Prague, Czech, ⁹Helmholtz-Zentrum Dresden-Rossendorf (HZDR), Rosendorf, Germany, ¹⁰Gesellschaft für Anlagen- und Reaktorsicherheit (GRS), Braunschweig, Germany, ¹¹University of Helsinki (UH), Helsinki, Finland, ¹²Centro de Investigaciones Energéticas, Medioambientales y Tecnológicas (CIEMAT), Madrid, Spain

After isolation of radioactive waste in deep geological formations, radionuclides can enter the biosphere via slow migration through engineered barriers and host rocks. The amount of radionuclides that migrate into the biosphere depends on the distance from a repository, dominant transport mechanism (diffusion vs. advection), and interaction of dissolved radionuclides with minerals present in the host rock and engineered barrier systems. Within the framework of the European Union's Horizon 2020 EURAD project (https://www.ejp-eurad.eu/), a series of state-of-the-art reports, which form the basis of a series of papers, have been drafted. This state-of-the-art paper aims to provide non-specialists with a comprehensive overview of the current understanding of the processes contributing to the radionuclide retention and migration in clay and crystalline host rocks, in a European context. For each process, a brief theoretical background is provided, together with current methodologies used to study these processes as well as references for key data. Owing to innovative research on retention and migration and the extensive knowledge obtained over decades (in the European context), process understanding and insights are continuously improving, prompting the adaptation and refinement of conceptual descriptions regarding safety assessments. Nevertheless, there remains important research questions to be investigated in the future.

KEYWORDS

geological disposal, radionuclide migration, diffusion, retention, sorption, redox chemistry, clay host rock, crystalline host rock
# 1 Introduction to transport controlling processes in clay and crystalline host rocks

After isolation of radioactive waste in deep geological formations and design-based repository closure, radionuclides are expected to enter the biosphere via slow migration. The amount of radionuclides that migrate into the biosphere depends on the distance of a repository from the biosphere, dominant transport mechanism (diffusion vs. advection), and interaction of dissolved radionuclides with minerals present in the host rock and engineered barrier systems (Churakov et al., 2020).

A geological disposal system (GDS) for radioactive waste has to fulfil three main “safety functions”: Containment, Isolation, and Retention (ONDRAF/NIRAS, 2013).

The host rock, as major part of this system, is expected to contribute to these functions by:

* Providing favorable and stable chemical and mechanical conditions to ensure good containment on a geological time scale.

* Isolating waste packages from man and the biosphere.

* Limiting radionuclide release to the biosphere via slow transport and sufficient retention in the host rock.

As deep clay and granitic formations offer the perspective required to fulfil these three safety functions, they are intensively studied by several countries, both in Europe and worldwide. In this paper, however, we will mainly relate them to the European context.

Within the EC EURAD project, Work package 5 (FUTURE: fundamental understanding of radionuclide retention—https://www.ejp-eurad.eu/implementation/fundamental-understanding-radionuclide-retention-future) deals with the fundamental understanding of retention and transport processes in clay and crystalline host rocks, and this manuscript describes the current state of knowledge.

Retention processes, generally called “sorption or uptake,” assure that potential radionuclide transport from a repository to the biosphere is much slower than the transport of groundwater. These processes depend on the pore water chemistry, which controls radionuclide speciation and influence the functional groups of reactive surfaces, present on the minerals. The most extensively studied “sorption” processes are surface complexation and ion exchange. Other processes, such as surface-induced reduction–oxidation (redox) uptake, incorporation, surface precipitation, and solid-solution formation, are less studied. However, they may be important retention mechanisms for a variety of radionuclides and are, therefore, gaining more attention. Section 2 describes these retention processes (clay surfaces will be taken as solid-phase examples to discuss the processes).

The two main transport mechanisms of radionuclides in host rock are diffusion and advection, driven by chemical concentration and hydraulic pressure gradients, respectively. Which mechanism predominates depends on the host rock matrix and the environmental conditions. This is described in detail in Section 3 for clay and crystalline host rocks.

In clay host rocks, diffusion is the dominant transport process. The negatively charged clay mineral surfaces not only result in a

high retention capacity for cationic radionuclides but also in a distinct diffusion behavior of anionic, cationic, and neutral species. The most important retention processes, influencing the transport rate, are ion exchange and surface complexation. As clay materials may exhibit strong reducing properties (i.e., a negative redox potential), due to the presence of reducing elements, Fe(II), in their crystal structure, surface-induced reduction–adsorption and reduction–precipitation in clays has become a field of intensive research.

In crystalline rocks, solute transport is mainly controlled by advection in connected fractures forming a preferential flow path and by matrix diffusion. The latter may be enhanced in weathered/ altered zones along the fractures. Moreover, the presence of fracture- filling materials, which are composed of clay minerals, Fe-oxy- hydroxides, and calcite, among others, may cause retention via different mechanisms such as ion exchange, surface complexation, surface-induced precipitation, and solid- solution formation.

This state-of-the-art report aims to provide a comprehensive overview of our current understanding of the processes that contribute to radionuclide retention and migration in clay and crystalline rocks within the European context. For each process, a brief theoretical background is provided, together with current methods used to study these processes as well as references to key publications.

Owing to decades of innovative research on the retention and migration of radionuclides, extensive knowledge and process understanding have been obtained, prompting the adaptation and refinement of conceptual descriptions regarding safety assessments. Uncertainties remain, and the key uncertainties that presently need to be resolved are listed in Section 4.

# 2 Retention processes

This chapter, on retention processes, is divided into three sections. In the first section, we briefly discuss aqueous speciation and the importance of thermodynamic data. This is a cornerstone for acquiring a good understanding and description of retention processes. In the two remaining sections, so-called “classical” retention processes, such as ion exchange and surface complexation, are discussed. These two processes are the dominant sorption phenomena and are the most extensively studied. These processes are described by taking clay minerals as an example but are also valid for other surfaces (e.g., oxides). In the second section, other relevant but less studied retention processes, such as surface-induced redox uptake, incorporation, surface precipitation, and solid-solution formation, are discussed.

## 2.1 Aqueous speciation, thermodynamic data, and thermodynamic databases

The study of RN aqueous speciation (i.e., its distribution among the different occurring chemical forms or species) determines the conditions under which contaminants can be mobile or can be retarded via (co) precipitation or sorption, making it fundamental
for understanding the environmental fate and transport of contaminants.

The most important factors controlling aqueous speciation are: pH, redox ($E_h$), salinity, and temperature; however, the presence of solid phases, colloids, and organic and inorganic ligands should also be accounted for. Many experimental methods can be used for determining radionuclide speciation. These include chemical analytical methods, such as potentiometry; calorimetry; spectroscopic techniques, including Ultraviolet-Visual-Near Infrared (UV-Vis-NIR); Raman; attenuated total reflectance Fourier-transform infrared (ATR-FT-IR); fluorescence or optical spectroscopy; time-resolved laser fluorescence spectroscopy (TRLFS); and nuclear magnetic resonance (NMR). Recent techniques, such as extended X-ray absorption spectroscopy (EXAFS), small-angle X-ray scattering (SAXS), and high-energy X-ray scattering (HEXS), have been successfully applied in this field, even in the case of actinides, and have benefited from recent developments on beamlines (Schacherl et al., 2022a; Schacherl et al., 2022b).

Some of the above-mentioned techniques provide indirect information on the formation of complexes in solution; others may give additional information on metal-ligand coordination modes, coordination numbers, and ion distances. Additional detailed information on the techniques used for actinide speciation can be found in a recent review (Patrice et al., 2016).

It is well accepted that speciation and solubility data should be obtained in well-defined solutions, with the knowledge of the nature (oxidation state and coordination environment) of the species involved in equilibrium reactions under consideration, and of the solid phases limiting the solubility. However, experimental studies cannot be performed under all the potential chemical conditions. Geochemical and thermodynamic modeling represents an important tool for the interpretation of empirical radionuclide speciation and solubility data and for predicting their behavior under different environments.

Thus, the main objective of experimental studies is to derive stability constants and all the thermodynamic parameters needed to feed geochemical models. Laboratory analyses provide a direct measure of different forms of an element, but chemical modeling is needed to apply all known thermodynamic relationships among chemical forms and to predict the overall equilibrium distribution. The model-based description of uptake, thus, includes several coupled physico-chemical processes that are described by thermodynamic and kinetic laws.

The mathematical approach to the analysis of aqueous speciation, above all in natural systems, is quite complex, and the development of dedicated software has been ongoing since the '80s. Far from being exhaustive, among the most known speciation software packages (and their updates) are: MINEQL (Westall et al., 1976); PHREEQE (Parkhurst et al., 1980) and PHREEQC (Parkhurst and Appelo, 2013); MINTEQA2 (Allison et al., 1991) and MINEQL+ (Schecher and McAvoy, 1992); EQ3/6 (Wolery and Daveler, 1992); CHESS and JCHESS (Van der Lee and de Windt, 2002); Geochemists' Workbench (Bethke and Yeakel, 2016); and GEMSs (http://gems.web.psi.ch/) (Kulik et al., 2013).

The quality of calculations strongly depends not only on the underlying assumptions, i.e., on the selection of chemical processes by a user, but also on the completeness and reliability of

thermodynamic data used as an input, and kinetics may also be important. Knowledge of fundamental thermodynamic properties, such as solubility products and complexation constants, is required. In addition, depending on the geochemical code used, there are code-specific dependent limitations such as temperature (only 25°C), pressure (only at 1 bar), ionic strength, and correction models. The combination of thermodynamic data completeness and the capabilities inherent to a geochemical code limit the modeling outcome.

The development of sound thermodynamic databases (TDBs) for radionuclide aqueous speciation and solid formation is, therefore, a prerequisite for geochemical calculations. In the last years, extensive research on the chemical thermodynamics and speciation of radionuclides in dilute, low-temperature groundwater systems has been carried out within the OECD-NEA TDB Project (www.oecd-nea.org/dbtdb/guidelines/tdb2.pdf). The nuclear Energy Agency compiles and evaluates existing thermodynamic data “to make available a comprehensive, internally consistent, quality-assured and internationally recognized chemical thermodynamic database of selected chemical elements in order to meet the specialized modeling requirements for safety assessments of radioactive waste disposal systems” (https://www.oecd-nea.org/dbtdb/). The Chemical Thermodynamics Series comprises review reports with selected data for elements of interest in radioactive waste management (https://www.oecd-nea.org/dbtdb/info/publications/). This work reflects present quantitative knowledge based on available experimental information and intends to be self-consistent. Despite the high quality of the databases, they are limited and insufficient to use alone, as many species are missing because data have not been evaluated according to the high standards of data selection.

Following the principles of NEA TDBs, some European agencies for nuclear waste management and some research institutions develop and continuously update their own, more exhaustive (and more widely applicable) databases. This is the case for the “ThermoChimie” TDB (https://www.thermochimie-tdb.com/), initially created and developed by ANDRA (French National Radioactive Waste Management Agency (Giffaut et al., 2014)) and later joined by “Radioactive Waste Management Limited” (NDA, UK) and ONDRAF/NIRAS (National Agency for Radioactive Waste Management, Belgium). It should be mentioned that SCK CEN developed first the MOLDATA database for ONDRAF/NIRAS, complemented with a software called Geochemical Data Processor (GDP) enabling internal consistency checking of data for the compilation of the in-house MOLDATA database (Wang et al., 2020a; Wang et al., 2020b). In Switzerland, the PSI/Nagra TDB (Hummel and Thoenen, 2021) is constantly updated and is provided in various formats compatible with most common geochemical tools (https://nagra.ch/wp-content/uploads/2023/05/NTB-21-03.pdf; https://www.psi.ch/en/les/database). Other examples are THEREDA (thermodynamic reference database for nuclear waste disposal in Germany, https://www.thereda.de/de/) (Voigt et al., 2007), the WEIMAR.dat (https://www.grs.de/publikationen/grs-500) database, which includes relevant mineral phases of the Gorleben site (Noseck et al., 2012), and EQ3/6 from the Lawrence Livermore National Laboratory (https://seaborg.llnl.gov/resources/geochemical-databases-modeling-codes).
Despite these important efforts, some gaps in data still exist, as well as uncertainties associated with the use of the models under certain conditions. Work is ongoing internationally to further develop the TDB and to extend the ability to model chemical speciation. It must be mentioned that, when dealing with natural systems, some of the necessary thermodynamic data are missing or not completely reliable (many experimental data likely exist but have not been evaluated entirely). Revisions and updates of these databases are still required, considering up-to-date experimental and theoretical studies published in peer-reviewed journals. For example, the lack of information on the kinetics of environmental processes obliges calculations assuming equilibrium. New studies should focus on specific chemical conditions, such as higher pH or high ionic strengths, or on the influence of other agents, such as organics or colloids, that may affect solubility, speciation, and sorption reactions and must be further analyzed.

Laboratory experiments, chemical analysis, modeling, and database maintenance and improvement are all part of an important iterative work that must be continued and even reinforced.

## 2.2 Classical retention processes

To assess the suitability of geological formations as potential host rocks for the deep disposal of high-level radioactive waste worldwide, detailed long-term safety studies are carried out, such as ONDRAF/NIRAS (2001), NAGRA (2002), SKB (2011), and POSIVA (2012b). Classical retention processes (adsorption) of radionuclides on rock substrates (vs. colloids) in the near and far field of a repository are an important component in such safety studies. It is a common practice to treat adsorption in terms of a distribution ratio, $R_d$, or distribution coefficient, defined as

$$R_d = \frac{C_{\text{sorbed}}}{C_{\text{eql}}} \qquad (1)$$

where $C_{\text{sorbed}}$ is the radionuclide concentration retained in the solid phase (mol/kg) and $C_{\text{eql}}$ is the equilibrium radionuclide concentration in the liquid phase [mol/L]. In the literature, most often, the term $K_d$ is used (same relationship), but this term entails that the sorption on a solid has linear dependency on the concentration in the solution and is reversible (when sorption is not linear, this term is often referred to as $R_d$). Sorption databases consisting of the selected $K_d$ values for safety relevant radionuclides are crucial for the performance assessments of a GDS but are important as well for surface disposal sites.

Clay minerals such as montmorillonite and illite (including illite-smectite mixed layers) play an important role in the capacity of engineered and geological barriers to retain radionuclides in many disposal concepts for radioactive waste, either as part of a host formation (argillaceous rocks), an engineered barrier system (bentonite), or natural fracture-filling materials in crystalline host rocks (see Section 3.3.2.3.1). In addition to clay minerals, other minerals exhibit sorption properties. For crystalline rocks, in particular, rock-forming minerals, such as feldspars, quartz, mica, and (hydr)oxides, are important contributors to radionuclide sorption. The two most important sorption mechanisms are ion exchange (especially for clays) and surface complexation.

In an ion-exchange process, the adsorption of cations occurs via Coulomb attractions at permanent, negatively charged surface sites. Ion-exchange sites are a typical feature of clay minerals (also see the next section where the process of ion exchange is explained based on the properties of clays). Moreover, 2:1 clay minerals, such as vermiculite, montmorillonite, illite, and illite-smectite mixed layers, exhibit large cation exchange capacities (CECs), whereas 1: 1 clay minerals, such as kaolinite, have low CECs.

Interactions are dominated by electrostatic forces and mainly result in so-called “outer sphere” complexes. The maximum adsorption capacity of a mineral or rock by ion exchange is denoted as the CEC. Cation exchange is strongly dependent on the ionic strength through competitive adsorption of a background electrolyte.

Surface complexation takes place on amphoteric surface groups (>SOH) present in minerals and is strongly dependent on pH because of the protolysis behavior of amphoteric groups (Dzombak and Morel, 1990). In clays, for example, functional groups are either related to Si atoms in tetrahedral sheets or Al atoms in octahedral sheets (>SiOH, >AlOH, >Al₂OH, and >AlSiOH). These groups are common to the edge surfaces of clay minerals (broken bonds). Generally, surface complexation leads to the formation of “inner sphere” surface complexes. In contrast, in an “outer sphere” complex, an ion retains its hydration sphere and attaches to the surface only via electrostatic forces.

In the next section, ion exchange and surface complexation are discussed using clay minerals as an example substrate because clays are an essential part of European geological disposal concepts in clays and crystalline rocks (bentonite buffer, and fracture fillings). The sorption by other rock forming minerals is discussed in Section 3.3.2.3.1. Various sorption-modeling approaches are discussed, and information on the atomistic nature of these classical adsorption mechanisms is addressed mainly for clays.

### 2.2.1 Ion exchange

#### 2.2.1.1 Theory and current understanding

Mineral surfaces may carry a permanent charge, and charge neutrality is maintained by the presence of an excess of counterions in solution held electrostatically in close proximity of the charged surface. Electrostatically bound ions, which constitute an electrical double layer, can undergo stoichiometric exchange with the ions in solution. This process is called ion exchange.

The surfaces of clay mineral platelets (also see Section 3.2.1 for more background on clay structural properties) carry a permanent negative charge arising from isomorphous substitution of the main lattice cations by other cations of lower valency, such as Si(IV) by Al(III) in a tetrahedral layer, and Al(III) by Mg(II) or Fe(II) in an octahedral layer (Liu et al., 2012b; Liu et al., 2013; Liu et al., 2014; Liu et al., 2015b). This permanent negative charge of a clay mineral gives rise to the so-called CEC (e.g., Grimm, 1953; Van Olphen, 1963). The CEC of montmorillonites (smectitic or swelling clay, with full access to the interlayer) is approximately 1 equiv/kg, and that of illites (non-swelling clay with only access to the outer surface) is approximately 0.2 equiv/kg. Notably, in the latter case, a large part of the structural charge is neutralized by potassium entrapped in the collapsed interlayer and is, thus, not readily available for cation exchange. By expressing the CEC in equivalent per unit mass, the value is independent of the charge of the exchangeable cations
compensating for the negative surface charge arising from isomorphous substitution in a clay mineral matrix. On this convention, the Gaines and Thomas formalism (Gaines and Thomas, 1953) is based. The description of equations for the law of mass action is expressed by a selectivity coefficient ($K_c$). Another well-known expression for the description of the cation exchange processes is given by the Vanselow Vanselow (1932) selectivity coefficient ($K_v$). The main difference between these two approaches is that $K_c$ is based on the equivalent fraction scale, whereas $K_v$ uses the molar fraction scale for the cations on the exchange complex. The current understanding of ion exchange in clay minerals is advanced because this topic has been a subject of intensive research for more than 40 years.

### 2.2.1.2 Ion-exchange models

**Law of mass-action models:** Cation exchange reactions are often expressed in terms of a selectivity coefficient obtained by the application of the mass action law. The cation exchange reaction of metal B, of valence $z_B$, exchanging with a metal A, of valence $z_A$, on a clay mineral in the A-form can be written as:

$$z_B A^{z_A} - clay + z_A B^{z_B} \rightleftharpoons z_A B^{z_B} - clay + z_B A^{z_A} \quad (2)$$

Following the convention given by Gaines and Thomas (1953), a selectivity coefficient, $B_A K_c$ [-], for reactions Equation 2 can be defined as:

$$B_A K_c = \frac{(N_B)^{z_A} [A]^{z_B} (y_A)^{z_B}}{(N_A)^{z_A} [B]^{z_B} (y_B)^{z_A}} \quad (3)$$

where $N_A$ and $N_B$ are the equivalent fractional occupancies, defined as the equivalents of A (or B) sorbed per kg of clay, divided by the CEC [equiv/kg]. [A] and [B] are the aqueous concentrations [mol/L], and $y_A$ and $y_B$ represent the aqueous-phase-activity coefficients [-].

A selectivity coefficient can be derived from experimental data at trace B concentrations ($N_A \sim 1$), using:

$$B_A K_c = (B R_d)^{z_A} \frac{(z_B)^{z_A}}{(CEC)^{z_A}} \cdot [A]^{z_B} \frac{(y_A)^{z_B}}{(y_B)^{z_A}} \quad (4)$$

where $B R_d$ [m³/kg] represents the sorption of cation B by cation exchange, defined as the number of moles of B sorbed per unit mass of the solid phase, divided by the number of moles of B in the aqueous solution per unit volume.

An extensive literature review of cation exchange on clay minerals has been compiled by Bruggenwert and Kamphorst (1982). In their survey, cation exchange constant $K_c$ (Gaines and Thomas, 1953), or $K_v$ (Vanselow, 1932), are listed and can be easily incorporated in geochemical software to model cation exchange. A tabulation and an evaluation of ion-exchange data on smectites have been carried out by Benson (1982). Data are presented as log(K) values and reflect very well the ranges found for some common ion equilibria on smectite. More recent data on cation exchange are available from Charlet and Tournassat (2005), Tournassat et al. (2007), Tournassat et al. (2008), Tournassat et al. (2009), and Tournassat et al. (2011).

Illite possesses a very specific type of exchange site, often called frayed edge sites (FESs), which have a high affinity for (and are only accessible by) large alkali ions because of their specific ionic size

matching the surface structure. Brouwer et al. (1983) have developed a three-site (FES, type II, and planar sites) cation exchange model for Cs and Rb sorption on illite. This model has been extended to a generalized Cs sorption model (Bradbury and Baeyens, 2000) to predict Cs sorption isotherms for argillaceous rocks. Recently, this three-site cation-exchange model could be implemented for $Tl^+$ adsorption on illite (Wick et al., 2018). Non-linear sorption of Cs has been also observed in some smectites, attributable to the presence of randomly interstratified illite-smectite mixed layers, providing highly selective sorption sites for alkaline metals. This behavior could be modeled considering a two-site exchange model accounting for sorption on planar and FES-like sites (Missana et al., 2014). An overview of cation-exchange systems relevant to geological disposal is given in Table 1, together with information on surface complexation (see Section 2.2.2).

### 2.2.1.3 Atomistic cation exchange models

Atomistic simulations are widely used to elucidate the hydration and adsorption mechanisms of ions in the interlayer of clay particles and the structure of the electric double layer at the basal plane of clay minerals. Such simulations are performed using Molecular Dynamics (MD) and Monte-Carlo approaches. In the Monte-Carlo technique, the equilibrium distribution of ions is obtained by stochastic modification of atomic positions following the Markov Chain stochastic process. In MD, simulation trajectories of ions and molecules are obtained by solving Newton's equation of motion for all atoms in a system. Thermodynamic properties and structural parameters of the system are obtained by time and ensemble averaging of the corresponding physical quantities depending on the atomic position over the entire simulated trajectory. The most important factors controlling the accuracy of the simulations are the underlying model for interatomic interactions and surface structural models used in the simulations. Most of the simulations are performed using so-called empirical pair interaction potentials. Such an approach provides a compromise between accuracy and computational efficiency. Only a few studies are known to have used a quantum mechanical system description (see Churakov and Liu (2018) for an exhaustive review). The most popular and widely used force fields for the simulations of clay minerals are CLAYFF (Cygan et al., 2004), INTERFACE (Heinz et al., 2013), and OPLS (Jorgensen and Tirado-Rives, 1988). These empirical force fields can, in principle, be combined with a variety of models for water (Guillot, 2002) and ions (Aqvist, 1990). More recently, several efforts were undertaken to develop a so-called polarizable force field, which takes into account the electronic polarization of electron shells depending on the instantaneous atomic configuration (Tesson et al., 2016).

The investigation of cation exchange in molecular simulations is based on calculations of the equilibrium distributions of ions between an interlayer space and a bulk electrolyte, in direct simulations with interfaces (Lammers et al., 2017) or by thermodynamic integration without interfaces (Rotenberg et al., 2009).

## 2.2.2 Surface complexation

### 2.2.2.1 Theory and current understanding

Different thermodynamic sorption models (TSM) have been introduced to describe the solid-solution interface during
<table><caption>TABLE 1 Overview of cation exchange and surface complexation models on 2:1 clay minerals.</caption><thead><tr><th>Clay</th><th>Model</th><th>Proto-lysis</th><th>SC</th><th>CE</th><th>Cations</th><th>Reference</th></tr></thead><tbody><tr><td>bentonite, zeolites</td><td>CE</td><td>-</td><td>-</td><td>K<sub>v</sub></td><td>Na<sup>+</sup>, K<sup>+</sup>, NH<sub>4</sub><sup>+</sup>, Ca<sup>2+</sup>, Ba<sup>2+</sup></td><td>Vanselow (1932)</td></tr><tr><td>montmorillonite</td><td>CE</td><td>-</td><td>-</td><td>K<sub>c</sub></td><td>Cs<sup>+</sup>, K<sup>+</sup>, Sr<sup>2+</sup></td><td>Faucher and Thomas (1954); Gaines and Thomas (1955)</td></tr><tr><td>montmorillonite</td><td>CE</td><td>-</td><td>-</td><td>K<sub>c</sub></td><td>Na<sup>+</sup>, Co<sup>2+</sup>, Ni<sup>2+</sup>, Cu<sup>2+</sup>, Zn<sup>2+</sup>, Cd<sup>2+</sup></td><td>Maes et al. (1975)</td></tr><tr><td>Illite</td><td>CE</td><td>-</td><td>-</td><td>K<sub>c</sub></td><td>Na<sup>+</sup>, Rb<sup>+</sup>, Cs<sup>+</sup>, Tl<sup>+</sup>, Ca<sup>2+</sup>, Sr<sup>2+</sup>, Ba<sup>2+</sup></td><td>Brouwer et al. (1983); Poinssot et al. (1999); Wick et al. (2018)</td></tr><tr><td>montmorillonite</td><td>CE</td><td>-</td><td>-</td><td>K<sub>c</sub></td><td>Na<sup>+</sup>, Fe<sup>2+</sup></td><td>Singhal and Kumar (1975)</td></tr><tr><td>montmorillonite</td><td>CCM</td><td></td><td>1 site</td><td>-</td><td>H<sup>+</sup></td><td>Wanner et al. (1994)</td></tr><tr><td>montmorillonite</td><td>DDLM</td><td></td><td>1 site</td><td>-</td><td>H<sup>+</sup></td><td>Kraepiel et al. (1998)</td></tr><tr><td>montmorillonite</td><td>MUSIC</td><td></td><td>4 sites</td><td>-</td><td>H<sup>+</sup></td><td>Tournassat et al. (2004b)</td></tr><tr><td>Illite</td><td>GEMS</td><td></td><td>2 sites</td><td></td><td>H<sup>+</sup>, Ln</td><td>Sinitsyn et al. (2000) Kulik et al. (2000)</td></tr><tr><td>montmorillonite<br>illite</td><td>2SPNESC/<br>CE</td><td>2 pK<sub>a</sub></td><td>strong/weak</td><td>K<sub>c</sub></td><td>Ni(II), Zn(II), Eu(III), Sn(IV), Am(III), Np(V), U(VI)</td><td>Baeyens and Bradbury (1997), Bradbury and Baeyens (1997b); Bradbury and Baeyens (2009a), Bradbury and Baeyens (2009b)</td></tr><tr><td>montmorillonite</td><td>EM/CE SC</td><td>2 pKa</td><td>strong/weak</td><td>K<sub>ex</sub></td><td>Ni(II), Zn(II)</td><td>Kraepiel et al. (1999)</td></tr><tr><td>montmorillonite</td><td></td><td>4 pKa</td><td>≡AlOH<br>≡SiOH</td><td>K<sub>ex</sub></td><td>UO<sub>2</sub><sup>2+</sup></td><td>McKinley et al. (1995)</td></tr><tr><td>montmorillonite</td><td>spillover</td><td></td><td></td><td></td><td>UO<sub>2</sub><sup>2+</sup></td><td>Tournassat et al. (2018)</td></tr><tr><td>montmorillonite</td><td>EM/CE<br>SC</td><td>2 pKa</td><td>1 site<br>KSOM</td><td>K<sub>ex</sub></td><td>Pb(II), Cd(II)</td><td>Barbier et al. (2000)</td></tr><tr><td>montmorillonite</td><td>DLM</td><td>3 pKa (two<br>≡AlOH, ≡SiOH)</td><td>≡AlOH<br>≡SiOH</td><td>neglected<br>I = 0.1</td><td>Np(V)</td><td>Turner et al. (1998)</td></tr><tr><td>montmorillonite<br>silica<br>alumina</td><td>EM<br>CCM</td><td>3 pKa (KA1,<br>KSi1 and KSi2</td><td>≡AlOH<br>≡SiOH</td><td>K<sub>ex</sub></td><td>Eu<sup>3+</sup>; UO<sub>2</sub><sup>2+</sup></td><td>Kowal-Fouchard et al. (2004b), Kowal-Fouchard et al. (2004a)</td></tr><tr><td>bentonite</td><td>DLM</td><td>3 pKa</td><td>1 site</td><td>n.i</td><td>U(VI), Se(IV)</td><td>Boult et al. (1998)</td></tr><tr><td>beidellite</td><td>TLM</td><td></td><td>2 sites<br>≡AlOH<br>≡SiOH</td><td>K<sub>ex</sub></td><td>UO<sub>2</sub><sup>2+</sup></td><td>Pabalan and Turner (1996)</td></tr></tbody></table>

SC: surface complexation; CE: cation exchange; K<sub>v</sub>: Vanselow constant; K<sub>c</sub>: Gaines & Thomas selectivity constant; K<sub>ex</sub>: cation exchange constant; CCM: constant capacitance modelling; EM: electrostatic modelling; DLM: Diffuse layer model (Dzombak and Morel, 1990); ECCM: (Extended) constant capacitance model (Nilsson et al., 1996); TLM: Triple-layer model (Hayes et al., 1991); MUSIC: Multi-site complexation model (Hiemstra et al., 1996); GEMS: Gibbs energy minimization selector (http://gems.web.psi.ch/).

sorption processes (Ochs et al., 2012). TSMs used to calculate defined surface complexation processes of metal cations (Me), for example, are called surface complexation models (SCMs). Each SCM is characterized by an individual set of surface complexation parameters (SCPs) defining the mechanistic approach to describe retardation processes at the solid-solution interface. These parameters comprise surface complexation constants (logK<sub>Me</sub>), specific surface area (SSA), surface site density (SSD), and protolysis constants (logK<sub>prot</sub>), which describe the acid-base characteristics of each mineral surface (Britz, 2018). TSMs extend the ion-association theory of aqueous species to surface species through thermodynamic equilibrium reactions using measurements of SSAs and estimated SSDs. In addition to the commonly applied diffuse double-layer model (DDLM) and non-electrostatic models (non-electrical double layers (EDL) models, noEDLM), more sophisticated

model approaches exist, such as the basic Stern model (BSM), triple-layer models, and three-plane models (TPMs), which are widely used TSMs. All TSMs are based on mass law equations and mole balance equations (Davis and Kent, 1990b), but they involve different descriptions of the electric double layer, i.e., an electrostatic interaction term (Westall and Hohl, 1980; Davis and Kent, 1990a). In the noEDLM, electrostatic interaction is not considered; thus, the DDL is neglected. Each TSM (except noEDLMs) considers different electrostatic planes that contribute to sorption processes.

### 2.2.2.2 Surface complexation models

**Acid-base modeling of minerals.** The protolysis behavior of amphoteric surface hydroxyl sites is usually described by following protonation and deprotonation reactions, with their corresponding mass action relationship:
<table><caption>TABLE 2 pKa's of edge sites for individual surfaces of phyllosilicates. Modified and augmented after Churakov and Liu (2018), with permission from Elsevier.</caption><thead><tr><th>Surfaces</th><th>Sites</th><th>pKa's</th><th>Reference</th></tr></thead><tbody><tr><td colspan="4"><strong>Pyrophyllite/Montmorillonite</strong></td></tr><tr><td rowspan="3">(110)<br>Neutral</td><td>≡Si(OH)<sub>2</sub>/≡Si(OH)</td><td>-6.8/8.3</td><td rowspan="3">Liu et al. (2014)</td></tr><tr><td>≡AlSiOH</td><td>1.7</td></tr><tr><td>≡Al(OH)<sub>2</sub></td><td>5.5</td></tr><tr><td rowspan="2">(110) Mg-sub</td><td>≡Si(OH)<sub>2</sub>/≡Si(OH) (connecting Mg via O)</td><td>-2.5/11.0</td><td rowspan="2"></td></tr><tr><td>≡MgSiOH</td><td>4.2</td></tr><tr><td rowspan="3">(110) Al-sub</td><td>≡Al(OH)<sub>2</sub><sub>Tetra</sub>/≡Al(OH)<sub>2</sub><sub>Tetra</sub></td><td>-2.4/17.5</td><td rowspan="3"></td></tr><tr><td>≡AlAlOH</td><td>10.2</td></tr><tr><td>≡Al(OH)<sub>2</sub><sub>Octa</sub></td><td>12.7</td></tr><tr><td rowspan="2">(010) Neutral</td><td>≡Si(OH)<sub>2</sub>/≡Si(OH)</td><td>-14.3/7.0<br>6.8</td><td rowspan="2">Liu et al. (2013), Tazi et al. (2012)</td></tr><tr><td>≡Al(OH)<sub>2</sub>(OH<sub>2</sub>)<sub>Octa</sub>/≡Al(OH)<sub>2</sub>(OH)<sub>Octa</sub></td><td>3.1/8.3<br>7.6/22.01</td></tr><tr><td rowspan="2">(010) Mg-sub</td><td>≡Si(OH)<sub>2</sub>/≡Si(OH) (connecting Mg via O)</td><td>-10.9/10.8</td><td rowspan="2">Tazi et al. (2012), Liu et al. (2014)</td></tr><tr><td>≡Mg(OH)<sub>2</sub>/2</td><td>13.2</td></tr><tr><td rowspan="2">(010) Al-sub</td><td>≡Al(OH)<sub>2</sub><sub>Tetra</sub>/≡Al(OH)<sub>2</sub><sub>Tetra</sub></td><td>-2.4/15.1</td><td rowspan="2">Liu et al. (2014)</td></tr><tr><td>≡Al(OH)<sub>2</sub>(OH<sub>2</sub>)<sub>Octa</sub>/≡Al(OH)<sub>2</sub>(OH)<sub>Octa</sub></td><td>4.9/8.5</td></tr><tr><td rowspan="2">(010) Fe(III)</td><td>≡Si(OH)</td><td>8.6</td><td rowspan="2">Liu et al. (2015a)</td></tr><tr><td>≡Fe(OH)<sub>2</sub>(OH<sub>2</sub>)/≡Fe(OH)<sub>2</sub>(OH)</td><td>1.2/5.1</td></tr><tr><td rowspan="2">(010) Fe(II)</td><td>≡Si(OH)</td><td>11.2</td><td rowspan="2"></td></tr><tr><td>≡Fe(OH)<sub>2</sub> [5-fold]/≡Fe(OH)<sub>2</sub>(OH<sub>2</sub>) [6-fold]</td><td>6.6/10.2</td></tr><tr><td colspan="4"><strong>Nontronite</strong></td></tr><tr><td rowspan="2">(010) Fe(III)</td><td>≡Si(OH)</td><td>8.7</td><td rowspan="2">Liu et al. (2015a)</td></tr><tr><td>≡Fe(OH)<sub>2</sub>(OH<sub>2</sub>)/≡Fe(OH)<sub>2</sub>(OH)</td><td>1.2/4.4</td></tr><tr><td rowspan="2">(010) Fe(II)</td><td>≡Si(OH)</td><td>9.9</td><td rowspan="2"></td></tr><tr><td>≡Fe(OH)<sub>2</sub> [5-fold]</td><td>8.4</td></tr></tbody></table>

$$ \equiv SOH + H^+ \Leftrightarrow \equiv SOH_2^+ K^+ = \frac{[\equiv SOH_2^+]}{[\equiv SOH]\{H^+\}} \quad (5) $$

$$ \equiv SOH \Leftrightarrow \equiv SOH^- + H^+ K^- = \frac{[\equiv SOH^-]\{H^+\}}{[\equiv SOH]} \quad (6) $$

where [ ] represents molar concentration and { } denotes activities.

Of note, the above mass action equations are written without an electrostatic term. In electrostatic models (see Table 2), an additional electrostatic term is included in the equations, which accounts for the distribution of the ions in a double layer according to the Maxwell-Boltzmann equation (e.g., Dzombak and Morel, 1990).

We would like to refer to two comprehensive review papers on this subject, related to clay minerals. The first, by Duc et al. (2005), details the experimental difficulties/constraints regarding titration measurements and their influence on modeling. The second review is by Bourg et al. (2007), giving an overview of available protolysis

models. This review contains data on pK, site densities, and specific edge surface areas from 12 different studies on acid-base titration models of Na-montmorillonite.

The site types, their densities, and the protolysis constants can vary in the various model approaches, from one single site to 27 sites (Tournassat et al., 2004b; Tournassat et al., 2004a). In general, two sites are considered and are associated with tetrahedral silanol sites and octahedral aluminol sites (Bradbury and Baeyens, 2000).

Table 1 lists a number of model approaches for the acid-base behavior of clay minerals, which are isolated from metal sorption data and can only be used for describing the proton buffer capacity of clay minerals.

**Sorption of anions.** The sorption of anionic species generally takes place via a ligand-exchange mechanism. For example, the sorption of boron on the edge sites of montmorillonite, illite, and kaolinite was modeled by Goldberg and Glaubig (1986) by means of the following reaction:
$$ \equiv SOH + H_3BO_3 \Leftrightarrow \equiv SH_2BO_3 + H_2O \quad (7) $$

with the associated mass action equation:

$$ K = \frac{[\equiv SH_2BO_3]}{[\equiv SOH][H_3BO_3]} \quad (8) $$

In this study, the boron sorption on the three clay minerals was modeled either by fitting protolysis constants to tetrahedral silanol sites and octahedral aluminol sites or by using protolysis constants from Al and Fe oxides (Goldberg and Sposito, 1984).

Modeled sorption data of anions on clay minerals are very sparse in the open literature. Most studies are limited to B(III), Mo(VI), Se(IV), and As(V) uptake on montmorillonite and illite (Goldberg and Glaubig, 1986; Goldberg and Glaubig, 1988a; Goldberg and Glaubig, 1988b; Motta and Miranda, 1989; Manning and Goldberg, 1996). More recently, Bruggeman (2006) investigated the uptake of Se(IV/VI) on illite and was capable of applying a non-electrostatic model (NEM) to quantitatively describe the pH-dependent adsorption of anionic species, using both wet chemistry and spectroscopic methods. Missana et al. (2009) investigated selenite sorption on smectite and illite and modeled the sorption behavior on smectite–illite mixtures, using a two-site NEM.

**Sorption of cations.** The general form of the surface complexation equations used for representing the surface binding of cations can be written as:

$$ \equiv SOH + Me^{z+} + yH_2O \Leftrightarrow \equiv SOMe(OH)_y^{z-}(y+1) + (y+1)H^+ \quad (9) $$

where Me is a metal with valency z, y is an integer, and the corresponding surface stability constant is expressed as $K_y$. For y = 0 the surface complex is $\equiv SOMe^{(z-1)}$.

In a NEM, the corresponding surface complexation constant, $K_y$, can be expressed as:

$$ K_y = \frac{[SOMe(OH)_y^{z-}(y+1)]}{[\equiv SOH]} \cdot \frac{f_{\equiv SOMe(OH)_y^{z-}(y+1)}}{f_{\equiv SOH}} \frac{\{H^+\}^{y+1}}{\{Me^{z+}\}} \quad (10) $$

where [ ] terms are surface concentrations, f terms are surface activity coefficients, and { } terms are aqueous activities. Following discussions given in Dzombak and Morel (1990), an assumption concerning surface activity coefficients is made, i.e., the ratios of the surface activity coefficients in the mass action equations describing surface complexation reactions are taken to be unity.

A model developed for montmorillonite and illite, the two-site protolysis non-electrostatic surface complexation and cation exchange, and 2SPNE SC/CE (Bradbury and Baeyens, 1997b; Bradbury and Baeyens, 2009a), can be used to illustrate the modeling approach using the following three datasets:

* clay acid-base titration data,

* sorption edges (trace sorption as a function of pH at a fixed ionic strength),

* sorption isotherms (concentration-dependent sorption at a fixed pH and fixed ionic strength).

A consequent iterative process was applied to model all three datasets until all experimental measurements were satisfactorily described. The model results were described with site capacities

(two weak protolysis sites ($\equiv S^{W1}OH$ and $\equiv S^{W2}OH$) and one strong site ($\equiv S^{S}OH$) for metal sorption at trace concentration), protolysis constants, cation-exchange capacities, cation-exchange selectivities, and surface complexation constants for strong and, in some cases, weak sites. Furthermore, there were enough elements investigated so that linear free energy relations could be derived for cation sorption (Bradbury and Baeyens, 2005b; Bradbury and Baeyens, 2005a).

An overview of the cation exchange and SCMs of anionic and cationic species for relevant geological disposal elements on 2:1 clay minerals is given in Table 1.

**Atomistic acid-base and surface complexation models.** The intrinsic acidity of surface OH groups is one of most important parameters controlling the development of surface charge on the edge sites of clay particles as a function of pH. Intrinsic acidity depends on the structural position of OH groups and the isomorphic substitutions in TOT sheets. Significant progress has been made recently in understanding the stability and structure of edge sites in clay minerals (Keri et al., 2020; Churakov, 2007; Churakov, 2006; Liu et al., 2012b; Okumura et al., 2017). The obtained data were used to investigate the acidity of surface groups using ab initio simulations (Liu et al., 2015b; Liu et al., 2014; Liu et al., 2013; Liu et al., 2011). The essence of the method is the calculation of free energy for the proton transfer reaction in the presence of explicit solvent molecules:

$$ \equiv SOH_n^p + H_2O_{aq} \Leftrightarrow \equiv SOH_{n-1}^{p-1} + H_3O_{aq}^+ \quad (11) $$

using thermodynamic integration techniques. The site acidities obtained in such a way represent so-called intrinsic pK's, corresponding to the deprotonation reaction of an isolated > SiOH site in pure water on an “ideal” charge neutral surface. To describe experimental data, these constants have to be combined with an electrostatic “analytical” model or one based on molecular simulations, which accounts for interactions between the surface charge and the ions in an electrolyte solution, using thermodynamic integration techniques (Churakov et al., 2014b). Table 2 summarizes currently available data obtained by ab initio simulations for different clay minerals.

Recent work by Orucoglu et al. (2022), Liu et al. (2022) and Gao et al. (2023b) allowed for obtaining a multi-scale understanding of structures, the reactivity of pyrophyllite (Schliemann and Churakov, 2021b; Schliemann and Churakov, 2021a), the and sorption capacity of montmorillonite (Kunipia) edges. It was demonstrated that Kunipia clay has crystalline edges and a cis-vacant layer structure. This atomic-scale study allowed for building and predicting edge acid-base properties and sorption properties, including sorption reversibility.

Ab initio simulations have been applied to reveal energetically favorable configurations of surface complexation ions on the edge surfaces of phyllosilicate minerals (Table 3). Most of the relevant studies were performed for the simple structural prototype of pyrophyllite (pyrophyllite has no structural charge). The general conclusions about the sorption mechanism should be transferable to other clay minerals, however. The extent and the distribution of surface charge has a strong influence on the strength of the ion-surface interactions.

TSMs described in previous sections rely on the existence of at least two sorption sites (Bradbury and Baeyens, 1997b; Baeyens and Bradbury, 1997). These sites are referred to as strong and weak sites,
<table><caption>TABLE 3 Ab initio surface studies of sorption complexation mechanism of cation on edges of clay minerals.</caption><thead><tr><td>Mineral</td><td>Cation</td><td>References</td></tr></thead><tbody><tr><td>Montmorillonite</td><td>Zn<sup>2+</sup></td><td>Churakov and Daehn (2012)</td></tr><tr><td>Pyrophyllite</td><td>Fe<sup>2+</sup>/Fe<sup>3+</sup></td><td>Keri et al. (2020)</td></tr><tr><td>Pyrophyllite</td><td>Fe<sup>2+</sup></td><td>Liu et al. (2012a), Liu et al. (2013), Liu et al. (2014)</td></tr><tr><td>Pyrophyllite</td><td>Cd<sup>2+</sup>/Ni<sup>2+</sup></td><td>Zhang et al. (2016), Zhang et al. (2017)</td></tr><tr><td>Kaolinite/pyrophyllite/illite</td><td>UO<sub>2</sub><sup>2+</sup></td><td>Kremleva et al. (2015a), Kremleva et al. (2015b), Kremleva et al. (2012), Kremleva et al. (2011)</td></tr></tbody></table>

which can be related to the higher and lower affinity constants of the surface complexation reaction. Ab initio simulations combined with X-ray absorption spectroscopy (XAS) studies specifically address the molecular structure of such sorption sites. It could be shown that, for divalent and trivalent transition metals, high-affinity sites are represented by, etch pits (single-atom defects) on the edge of clay minerals. The weak sites in the TSM are structurally related to the inner-sphere sorption complexes on the edge surfaces (Churakov and Daehn, 2012; Keri et al., 2020). The structure and site density of the cationic surface complexes obtained by ab initio simulations are consistent with spectroscopic data and thermodynamic simulations, respectively.

### 2.2.3 Reversibility

An analysis of radionuclide desorption behavior is essential for evaluating the (ir)reversibility of sorption processes in order to predict the long-term migration behavior of radionuclides.

Desorption tests are typically performed following sorption tests. In the sorption step, after the selected contact time, the solid-to-liquid distribution coefficient, $R_{d,sorb}$, is determined. The solid in which the adsorbate is retained and can be suspended again in a fresh solution to determine the distribution coefficient, $R_{d,desorb}$, for the desorption step. From a thermodynamic perspective of reaction chemistry, when $R_{d,sorb}$ and $R_{d,desorb}$ are equal for reaction times within the same time span of hours to a few days, this process is considered “reversible”; if not, it is considered “irreversible.” For longer reaction times, other processes, such as recrystallisation, the restructuring of surface complexes, or the neoformation of clay minerals, might occur, changing an initial reversible surface complex to an irreversible one over years.

The reasons why retention is (or is not) reversible are numerous and must be evaluated in detail. In some of the cases, irreversibility is only apparent because the sorption and desorption processes were not at equilibrium. However, hysteretic behavior might indicate that reactions other than sorption are taking place, which contribute to the overall retention in a system.

The main mechanism leading to contaminant retention in a mineral surface is adsorption, which is generally rapid and reversible, when no other process such as the recrystallization restructuring of surface complexes or neoformation of clay minerals occur (Dähn et al., 2001; Dähn et al., 2002; Dähn et al., 2003; Dähn et al., 2011; Begg et al., 2017). However, several kinetically controlled processes contributing to slower uptake exist, such as diffusion processes inside a solid where the formation of stronger bonds and the crystallization of new phases can occur (Sparks, 2003). In addition, aging can have a

significant effect on contaminant release. These processes are usually not included in surface complexation modeling.

The importance of the “reaction time” and kinetically controlled processes increases when dynamic phenomena (diffusion, advective-dispersive transport) must be described. An apparently “irreversible” process for short residence times can be completely reversible at the long-time scales considered in safety assessment for radioactive waste repositories.

Among the most important processes contributing to irreversible sorption on clays is the formation of surface precipitates incorporating and sequestrating a contaminant (also see Section 2.2, other uptake processes). The identification and characterization of surface precipitates is an important task that has to be accomplished using microscopic and spectroscopic techniques (Scheidegger et al., 1997; Ford et al., 2001).

Structural modifications within the layered structure of phyllosilicates (clays) can lead to irreversible retention of contaminants. Many recent studies have been devoted to the analysis of the irreversible sorption of Cs into illite, focusing on processes favoring its fixation caused by the collapse of the FESs, and their de-lamination favored in the presence of divalent ions (Benedicto et al., 2014c; Fuller et al., 2015; Durrant et al., 2018).

One of the key factors determining the degree of reversibility of the surface complexation on the edges of clay minerals at an atomistic scale is the structural evolution of the mineral surfaces as part of dissolution and precipitation processes. These processes determine the availability of high-affinity sorption sites and the structural (irreversible) incorporation of cations. Ab initio simulations can provide direct information on the dissolution process for individual surface sites (Schliemann and Churakov, 2021a). The macroscopic description of these data needs a global view of the temporal evolution of site density. Such information can be obtained using coarse-grained models for clay particle dissolution and growth as well as mineral-fluid interactions (Kurganskaya and Luttge, 2013; Delhorme et al., 2010; Yang et al., 2020).

### 2.2.4 Research methods

Different approaches exist to characterize mineral-specific and sorption-relevant processes, such as the acid-base behavior of surfaces and the retardation of cations (e.g., cation-exchange and surface complexation processes), in general. To examine sorption reactions at the solid-liquid interface, commonly applied methods include classical potentiometric titration experiments, so-called batch experiments to elucidate adsorption, and spectroscopic investigations.
#### 2.2.4.1 Wet chemistry

**Determination of CEC.** A large range of methods exists to determine the CEC. Those using the highly selective metal dyes (e.g., Cu(trien) (Meier and Kahr, 1999), Ni(en)₃ (Baeyens and Bradbury, 1991), or CoHex (Rémy and Orsini, 1976)) are the most popular because of their simplicity (single extraction, possibility of colorimetric measurements of selective metal-ion consumption, etc.).

Within the Mont Terri Project, an extensive intercomparison study was set up on different methods in order to find their advantages, disadvantages, and potential drawbacks (Hadi et al., 2019).

**Acid-base titrations.** Wet chemistry investigations on dispersed systems are aimed at elucidating the nature of the radionuclide retention at the solid-liquid interface, i.e., to characterize the different adsorption site types and capacities. Detailed knowledge of the protonation and deprotonation processes (acid-base behavior) of mineral surfaces are of great importance for realistically describing retardation processes at mineral-water interfaces, since they are—for most radionuclides/ contaminants—pH dependent.

Potentiometric acid-base titration experiments of mineral surfaces have been widely applied (e.g., Arnold et al., 2001; Lützenkirchen et al., 2012a; Lützenkirchen et al., 2012b; Schwarz et al., 1984) to characterize mineral electrostatic behavior in terms of points of zero charge, surface charge densities, and net surface proton excess. A combination of discontinuous electrolyte and mass titration experiments with continuous potentiometric acid-base titrations offer even deeper insight into mineral surface protonation and deprotonation states as a function of, for example, pH, ionic strength, and solid-liquid ratio (e.g., Lützenkirchen et al., 2012a; Lützenkirchen et al., 2012b; Preocanin and Kallay, 1998; Noh and Schwarz, 1989). An alternative experimental approach (Baeyens and Bradbury, 1995; Baeyens and Bradbury, 1997; Bradbury and Baeyens, 2009a; Tournassat et al., 2004b) is a batch titration technique with back-titration on aluminum oxides, as originally proposed by Schulthess and Sparks (1986).

**Adsorption experiments.** The research methods for ion exchange and surface complexation mainly target the elucidation of the nature of the uptake mechanisms (inner sphere/outer sphere complexes) and the type and capacities of the adsorption sites of sorbents investigated. Wet-chemistry investigation on dispersed systems is the most common experimental method used.

Sorption experiments are carried out by equilibrating an adsorbate solution of known composition and volume with a known amount of adsorbent. The suspension is gently shaken/stirred for a certain period of time until sorption equilibrium is reached. After equilibrium, the solution is separated from the solid phase via centrifugation and/or filtration and analyzed. From these experiments, the distribution coefficient $R_d$ is calculated (see Equation 1). When a series of experiments is conducted at a well-defined ionic strength and similar initial concentrations but at varying pH values, so-called "sorption edges" are obtained. When a series of experiments is conducted at a well-defined ionic strength and constant pH but by varying the initial concentration, so-called "adsorption isotherms" are determined.

**Spectroscopic and atomistic techniques.** Surface complexation studies of ionic species are generally carried out using batch sorption

techniques similar to those used in ion-exchange studies. These methods are frequently complemented by a number of alternative techniques such as micro-calorimetry and infrared spectroscopy, to elaborate on the nature of the surface complexes. X-ray Adsorption spectroscopy (XAS) techniques (Extended X-ray Adsorption Fine Structures (EXAFS)/X-ray Adsorption Near Edge Structure (XANES), see Section 2.3.2.3.1 for more details) have been applied to demonstrate the outer-sphere character of ions in an interlayer.

EXAFS has proven successful in providing a molecular understanding of sorption on edge sites (Dähn et al., 2003; Schlegel et al., 2001b; Dähn et al., 2011), uptake on interlayer sites via cation exchange (Sposito, 1984; Muller et al., 1997; Chisholm-Brause et al., 1994), and uptake in newly formed lamellar phases such as layer silicates (Schlegel et al., 2001a; Dähn et al., 2002) and mixed layered double hydroxides (LDHs) (Scheidegger et al., 1998; Thompson et al., 1999a; Towle et al., 1997).

MD and Monte-Carlo simulations are used to reveal the equilibrium distribution of ions in the interlayers of smectites minerals, as well as cation-exchange processes and the structure of the electrical double layer at solid-liquid interfaces. Surface complexation processes involve protonation and deprotonation reactions at the edges of clay minerals. The surface complexation reactions result in the direct interaction of ions with deprotonated oxygen sites. These processes lead to a significant change in the electron density between the interacting atoms and, thus, require a quantum mechanical description of systems. Recent reviews of quantum mechanical studies of adsorption processes on clay minerals are available in Churakov and Liu (2018) and (Liu et al., 2022).

## 2.3 Other uptake processes

### 2.3.1 Neoformation/surface precipitation

#### 2.3.1.1 Theory and current understanding

Most of the studies on the uptake of contaminants in clay systems have been based on wet chemistry experiments and modeling at the macroscopic level. Batch studies (also see Section 2.2.4) provide an efficient tool for determining the distribution coefficients of metal ions and, under certain circumstances, allow for differentiation between ion exchange and surface complexation processes occurring at the solid-water interface (Sposito, 1984). However, the ability of a particular SCM to fit macroscopic data does not certify that the assumptions underlying a model are correct. For example, macroscopic data do not allow unambiguous discrimination between surface complexation and nucleation processes. A way to gain mechanistic information about metal interactions on solid surfaces is to combine batch studies with synchrotron-based spectroscopic investigations.

With increasing reaction time and metal concentrations, transition-metal surface precipitates (neo-formed phyllosilicate or LDHs) can form on clay minerals. This process can irreversibly bind contaminants to mineral surfaces, whereas (pure) sorption models assume that the adsorption of trace metals is reversible. Using EXAFS, it has been demonstrated that Ni-, Co- and Zn-containing precipitates can form when clay minerals and Al- and Si-(hydr)oxides are treated with Ni, Co, and Zn, and this occurs even when the initial metal concentration in solution is undersaturated
relative to that of the pure (oxyhydr)oxide forms of the metal (Towle et al., 1997; Scheidegger et al., 1998; Manceau et al., 1999; Thompson et al., 1999a; Thompson et al., 1999b; Morton et al., 2001; Schlegel et al., 2001a; Dähn et al., 2001; Dähn et al., 2002; Lee et al., 2004). For example Scheidegger et al. (1998) observed the formation of a LDH phase when pyrophyllite (a 2:1 clay that lacks isomorphic substitution) was treated with Ni. A Ni-Al LDH phase formed after a contact time of only a few minutes between pyrophyllite and Ni, suggesting that nucleation can occur rapidly in metal clay sorption systems. LDH formation has been further reported for Co(II) and Zn(II) sorption on kaolinite and Al-(hydr)oxides (Towle et al., 1997; Thompson et al., 1999a; Thompson et al., 1999b). Similarly, the uptake of Co on quartz was shown to result in the neoformation of a trioctahedral clay-like structure (Manceau et al., 1999).

Using polarized-EXAFS (P-EXAFS) Dähn et al. (2002) and Schlegel et al. (2001a) demonstrated the presence of neo-formed phyllosilicates after treating di- and tri-octahedral 2:1 clay minerals (montmorillonite and hectorite) with Ni and Zn solutions at elevated pH and metal concentrations. The neo-formed phases were oriented parallel to the octahedral smectite sheets. In a study by Schlegel et al. (2001a), a structural link of neo-formed phases to hectorite particles could be observed. A study of Dähn et al. (2002) indicated that neoformation processes of clay minerals had already started after 1 day of reacting and continued for as long as metal ions were present in solution.

### 2.3.1.2 Research methods

As mentioned in Section 2.2.4.2, EXAFS has proven successful in gaining a molecular understanding of sorption on edge sites, uptake on interlayer sites via cation exchange, and uptake in newly formed lamellar phases. In the relevant studies, macroscopic wet chemistry investigations were combined with synchrotron-based spectroscopic analysis. Studies on the uptake of metals on mineral surfaces observed the formation of phyllosilicates upon the uptake of Co on quartz (Manceau et al., 1999), Zn on hectorite (Schlegel et al., 2001a), and Ni on montmorillonite (Dähn et al., 2001; Dähn et al., 2002), as well as the formation of LDH phases for Co uptake on Al-(hydr)oxides (Towle et al., 1997), Zn on kaolinite and Al-(hydr)oxides (Thompson et al., 1999b; Thompson et al., 1999a), and Ni on gibbsite, pyrophyllite, and montmorillonite (Scheidegger et al., 1997; Scheidegger et al., 1998). The sensitivity of EXAFS measurements on clay minerals can be significantly improved by performing P-EXAFS. In P-EXAFS, neighboring atoms along the polarization direction of an X-ray beam are preferentially probed, and atoms located perpendicular to this direction are attenuated. Applying P-EXAFS to self-supporting clay films has the advantage of minimizing the contributions from out-of-plane Si atoms from a tetrahedral sheet, when an X-ray polarization vector is in the *ab* plane of a montmorillonite self-supporting film (Manceau et al., 1998; Dähn et al., 2002; Dähn et al., 2003; Dähn et al., 2011). Conversely, when the polarization vector is aligned normal to the film plane, the contribution from octahedral-layer cations vanishes. This method allows for determining the speciation of chemical entities at the atomistic scale, even when the element of interest is present at low concentrations (concentrations of X-ray absorber down to approximately 100 ppm in clay matrixes).

### 2.3.1.3 Remaining uncertainties

Research on the neoformation of minerals raises many questions, such as how can stable mineral clay phases generate sufficient soluble Al or Si to react with sorbed metals? None of the above-mentioned studies could provide deep insights into what triggers the precipitation of a phyllosilicate or a LDH phase. In both cases, it was observed that high pH, high metal concentrations, and prolonged reaction times favored the formation of phyllosilicates or LDH precipitates. In the past, most efforts in radioactive waste management studies focused on the understanding of sorption process at low metal concentrations in solution, as it was assumed that precipitation processes are not relevant for the safety analyses, considering a conservative approach, of future radioactive waste disposal sites. However, nowadays it has been accepted that, for example, high Fe concentrations arising from waste forms can occur, deeming the precipitation processes relevant (Soltermann et al., 2014b; Soltermann et al., 2013b; Soltermann et al., 2014c). In order to understand the neoformation of phyllosilicates or LDH phases from a thermodynamic perspective, solubility products need to be determined in the future (Dähn et al., 2003). Presently, such information is lacking in TDBs. The heterogeneous formation of phyllosilicates has important geochemical implications because layer silicates are stable minerals under mildly-acidic-to-basic pH conditions and can irreversibly bind metals in waste and soil matrices, as demonstrated in Zn-smelter impacted soils (Manceau et al., 2000; Vespa et al., 2010).

For long-duration sorption experiments, the development of microbial activity with time could interfere with sorption processes and neoformation/surface precipitation. There is, however, not much documentation on the effect of microbial perturbation on the results of batch sorption tests. It is suggested to thoroughly review the literature and to consider this in future studies.

## 2.3.2 Surface-induced redox uptake

### 2.3.2.1 Theory and current understanding

#### 2.3.2.1.1 General definition of redox processes and selected examples.

Redox reactions require the transfer of electrons between aqueous species (homogenous) or between aqueous species and mineral surfaces (heterogeneous) (Ahmed and Hudson-Edwards, 2017). Redox processes are key geochemical processes because they strongly influence the mobility of redox-sensitive contaminants (e.g., Se, Tc, U, Np, and Pu), which are abiotically reduced to often less toxic or less mobile forms (Ma et al., 2019; Guillaumont et al., 2003; Olin et al., 2005; Winkel et al., 2012). These include the tetravalent oxides of U, Np, Pu, and Tc; Se(0); FeSe; and sulfide phases such as TcS, while notable exceptions to this rule include the aqueous $Pu^{3+}$ complexes and gaseous Se species (Kirsch et al., 2011; Winkel et al., 2012).

Minerals containing redox-active species such as Fe (II, III) or Mn (II, III, IV) can participate in a wide range of surface-induced electron-transfer reactions. For instance, Fe-bearing minerals cover a large range of redox potentials depending on the amount and location of Fe in a crystal lattice (Gorski et al., 2012a; Gorski et al., 2012b; Gorski et al., 2013; Gorski et al., 2016; Sander et al., 2015; Li et al., 2019a). The rates and extents of redox reactions are controlled by the physico-chemical properties (e.g., redox and sorption potential, or electron conductivity, depending on factors such as composition or particle size) of the minerals (Ilgen et al., 2019; Scheinost et al., 2008).
In the last two decades, many studies have addressed the influence of Fe-bearing minerals on the immobilization of redox-sensitive nuclides such Fe, Se, Tc, U, Np, and Pu. Studies on relevant minerals are summarized below:

* Clays: Alexandrov and Rosso, 2013, Boland et al., 2011, Chakraborty et al., 2010, Fox et al., 2013, Hoving et al., 2017, Ilton et al., 2010, Jaisi et al., 2009, Jeon et al., 2005, Jones et al., 2017, Latta et al., 2012, Latta et al., 2017, Liger et al., 1999, O'Loughlin et al., 2003, Pearce et al., 2017, Peretyazhko et al., 2008, Peretyazhko et al., 2012, Roberts et al., 2019, Schaefer et al., 2011, Soltermann et al., 2014a, Soltermann et al., 2013a, Stumm and Sulzberger, 1992, Tsarev et al., 2016, Yang et al., 2012, Zhu and Elzinga, 2014, Joe-Wong et al., 2017, Charlet et al., 2007, Frohlich et al., 2012, Brookshaw et al., 2015, Qian et al., 2023.

* Iron oxides (including oxyhydroxides and hydroxides): Bender and Becker, 2019, Dumas et al., 2019, Emmanuel and Ague, 2009, Emmanuel et al., 2010, Kirsch et al., 2011, Kobayashi et al., 2013, Um et al., 2011, Yalçıntaş et al., 2016, Scheinost and Charlet, 2008, Scheinost et al., 2008, Christiansen et al., 2011, Wylie et al., 2016, O'Loughlin et al., 2003, O'Loughlin et al., 2010, Huber et al., 2012, Pidchenko et al., 2017, Roberts et al., 2017, Börsig et al., 2018, Peretyazhko et al., 2008, Poulain et al., 2022.

* Fe(II) (hydroxo-) carbonates: Scheinost et al., 2016; Scheinost and Charlet, 2008; Kirsch et al., 2011; Ithurbide et al., 2009; Badaut et al., 2012; Llorens et al., 2007)

* Fe(II) sulfides: (Moyes et al., 2002, Scheinost and Charlet, 2008, Scheinost et al., 2008, Yalçıntaş et al., 2016, Finck et al., 2012, Han et al., 2011, Kirsch et al., 2011, Wharton et al., 2000, Rodriguez et al., 2020, Breynaert et al., 2008, Liu et al., 2008, Bruggeman et al., 2007, Pearce et al., 2018, Livens et al., 2004.

Element-specific reviews, for example those focusing on Tc (Meena and Arai, 2017; Pearce et al., 2020), can also be used as a state-of-the art introduction on the nature of the redox interactions between Fe minerals and redox-sensitive elements.

### 2.3.2.1.2 Redox-relevant minerals in disposal systems.
The following section is based on the mineralogical data available for clay formations, namely Boom Clay, Opalinus Clay, and the clay-rich horizon of the Callovian-Oxfordian Formation in Bure and for the granitic formation of Aspö.

As discussed above, for a redox reaction to occur at the mineral-aqueous solution interface, the solid phase must contain a chemical element that can have different oxidation states in the considered chemical system, and it must be in a reduced state if the element with which it interacts is in the oxidized state, and vice versa. In the environment, two major metallic elements reported to be capable of participating in redox reactions are Fe and Mn. In the following section, Mn will not be discussed specifically. Although Mn oxyhydroxides are extremely reactive (Manceau et al., 1997; Means et al., 1978), the deep geological formations that are reviewed here contain little amounts of Mn and no discrete Mn oxyhydroxides (Gaucher et al., 2004; Landstroem and Tullborg, 1995; Lerouge et al., 2014; Zeelmaekers et al., 2015). It is worth mentioning that Mn oxyhydroxides may play a role in the case of the subsurface storage of low-activity wastes (Debure et al., 2018).

Carbonate minerals, sulfide minerals, and phyllosilicates in pristine rock. Fe can be found in a variety of minerals from different groups, as Fe(II), Fe(III), or in mixed Fe(II)/Fe(III) oxidation states. In most deep clay-rich formations, three main mineralogical “groups” contain Fe: sulfide minerals, carbonate minerals, and phyllosilicates. The most frequent sulfide mineral bearing Fe is pyrite (FeS₂), but in spite of a few reports possibly related to the difficulties associated with identifying this poorly crystalline mineral, its precursor mackinawite may play a significant role owing to its high redox reactivity (Grambow, 2016; Breynaert et al., 2010; Frohlich et al., 2012; Zeelmaekers et al., 2015; Bruggeman et al., 2007). Concerning carbonate minerals, siderite (FeCO₃), a pure Fe endmember, is regularly observed in clay-rich formations, but Fe-containing calcite is also common, and ankerite is regularly observed (Kars et al., 2015; Bossart et al., 2017; Zeelmaekers et al., 2015). In all carbonate minerals, Fe is present in the Fe(II) oxidation state, with a very variable concentration (Lerouge et al., 2014). Finally, several phyllosilicates are known to host iron. Among them, the most common are likely chlorite, biotite, illite (including illite layers in mixed-layered minerals), and micas (see data compilation in Lerouge et al., 2017; Finck, 2020), with distinguishing between illite and mica being almost difficult because of the close similarity in their crystal structure and because they are frequently mixed in a given sample (Grangeon et al., 2015). Thus, it is likely that the proportions of illite, illite-smectite mixed-layered minerals, and micas presented in literature studies exhibit uncertainties. In the Callovian-Oxfordian level of the deep geological laboratory of Bure, di-octahedral interlayer-deficient micas (mainly illite and illite-rich illite-smectite mixed-layered minerals) represent 36%-71% of Fe-bearing minerals, and their FeO content ranges from 1.7 to 5.7 wt%. Chlorites typically represent 2% of Fe-bearing minerals, and their FeO content is approximately 20 wt%. Finally, other phyllosilicates are present at trace concentration, but their contributions to the total Fe abundance may not be negligible. For example, although biotite is present at trace concentration, its high FeO content of 13–22 wt% means that it could contribute significantly to the total Fe content in a given sample where it would be enriched. Most of the Fe associated with phyllosilicates occurs in the Fe(II) oxidation state and is present in a phyllosilicate structure, but a minor fraction is in an exchangeable position, making the amount of exchangeable Fe in the Callovian-Oxfordian formation, at the depth of the Bure laboratory, able to reach approximately 1.2 meq/kg (Tournassat et al., 2008).

In granitic geological environments, phyllosilicates are rock-forming minerals. For example, at Aspö, biotite and chlorite are observed, with the latter possibly resulting from hydrothermal alteration (Morad et al., 2018). The presence of illite-rich illite-smectite mixed-layered minerals can also be observed as fracture-filling material, but their abundance as fracture-filling material is relatively low, typically less than 10 wt% (Stanfors et al., 1999).

(Oxyhydr)oxides in pristine rock. In addition to the above-mentioned minerals, Fe (oxyhydr)oxides, hereafter referred to as “Fe oxides,” can also be found. In most clay-rich formations, their abundance is low; therefore, they probably have a negligible role in contributing to the whole-rock sorption capacity. However, they are of great importance for the understanding of rock geochemistry,
since they control the Fe in porewater, even if present at trace concentration. For example, in the Callovian-Oxfordian level of Bure, goethite (FeO(OH)) abundance is typically less than 0.1 wt%, a trace concentration that is, however, sufficient to control Fe in the pore water (Gailhanou et al., 2017; Kars et al., 2015).

In granitic rocks that are not subject to hydrothermal reactions, Fe oxides are commonly observed as part of the mineralogical assemblage that form fracture-filling material near the surface, but not in deep regions (Mathurin et al., 2014), coherent with their formation by the alteration of Fe-rich phyllosilicates (e.g., biotite). The typical small size of Fe oxides makes them prone to being mobilized by an advective water flow, which would be coherent with the observation of Fe-rich colloids in open fractures (Degueldre et al., 1989).

Potential influence of construction and waste package materials on the neoformation of Fe oxides. In many disposal concepts, construction materials remain in a repository, either because they will be used to build access structures and galleries or because they will be part of waste packages (metallic canisters and steel overpacks). Among these materials, steel in contact with either clay or cement will potentially induce the formation of Fe(II)-bearing minerals (Bildstein and Claret, 2015). In the case of the corrosion of steel in contact with clays, oxides such as magnetite, (hydroxy-)carbonates such siderite and chukanovite, Fe sulfides, and Fe-rich phyllosilicates have been observed, depending on the geochemical environment and redox conditions (Schlegel et al., 2010; El Mendili et al., 2014; Schlegel et al., 2014; Necib et al., 2016; El Hajj et al., 2013). In both freshwater and saline conditions, mixed Fe(II)/Fe(III) LDHs, termed green rusts, may form (Bach et al., 2014; Christiansen et al., 2011; Grambow et al., 1996). Contrastingly, when steel is corroded under cement pore-water conditions, Fe(OH)₂ and Fe(OH)₃ can be stabilized (Ma et al., 2019).

In addition, cement in contact with a clay-rich formation can destabilize clay minerals and lead to the precipitation of new phases that are enriched in Mg, as observed in a Cement-clay Interaction (CI) experiment in Opalinus Clay conducted at Mont Terri (Dauzères et al., 2016; Jenni et al., 2014; Lerouge et al., 2017; Mäder et al., 2018). These phases, often termed “M-S-H,” are nanocrystalline defective Mg phyllosilicates that have a layer structure close to that of stevensite (Roosz et al., 2015) and that contain appreciable amounts of Fe(III) (Bonen, 1992; Lerouge et al., 2017).

The mixed Fe(II)/Fe(III) oxide magnetite is a common corrosion product under anoxic repository conditions, experimentally evidenced at the steel-clay rock interfaces of waste canisters embedded in clay backfill (Schlegel et al., 2010; El Mendili et al., 2014; Schlegel et al., 2014; Necib et al., 2016; El Hajj et al., 2013). Magnetite, especially as a pristine nanoparticle formed under strictly anoxic conditions with near-ideal stoichiometry, rapidly reduces a suite of relevant radionuclides, including Se, Tc, U, Np, and Pu (Bender and Becker, 2019; Dumas et al., 2019; Emmanuel and Ague, 2009; Emmanuel et al., 2010; Kirsch et al., 2011; Kobayashi et al., 2013; Um et al., 2011; Yalçıntaş et al., 2016; Scheinost and Charlet, 2008; Scheinost et al., 2008; Christiansen et al., 2011; Wylie et al., 2016; O'Loughlin et al., 2003; O'Loughlin et al., 2010; Huber et al., 2012; Pidchenko et al., 2017; Roberts et al., 2017; Börsig et al., 2018). The mechanisms of reduction follow a complex multi-step process, which is dependent on the chemical conditions of a solution (e.g.,

pH). For example, in the case of Se(VI), reduction is a two-steps process, starting from reduction to Se(IV) and then to trigonal Se(0). During reduction, pH oscillations may occur, leading to the formation of Fe(II) that re-adsorbs to magnetite and renews its reduction capacity (Poulain et al., 2022). This latter effect is observed at pH five and not at pH 7, highlighting the need to have a complete understanding of the chemical conditions, to be able to predict the fate of redox-sensitive elements upon interaction with magnetite.

In the case of Tc, Tc(VII) is reduced by magnetite to Tc(IV), which either forms small polymeric chains with a TcO₂·xH₂O-like structure sorbed onto the magnetite surface or is structurally incorporated by substituting for Fe in octahedral positions in the magnetite structure (Zachara et al., 2007; Peretyazhko et al., 2012; Kobayashi et al., 2013; McBeth et al., 2011). Incorporated Tc is better protected against reoxidation (Marshall et al., 2014; Um et al., 2017) and migration; hence, understanding the control between both processes is essential for the safety of repositories. From this perspective, numerical methods such as Density Functional Theory (DFT) can help in providing meaningful insights at the molecular-scale level by allowing elucidating the structure of nanocrystalline and cryptocrystalline phases, as exemplified with TcO₂·nH₂O (Oliveira et al., 2022). Earlier work suggests that the structural incorporation of Tc might be favored under pH conditions, where magnetite is less stable, offering a more dynamic surface with frequent Fe dissolution/reprecipitation events, through which Tc(IV) with a similar ionic radius and coordination could be trapped into the magnetite structure (Yalçıntaş et al., 2016).

In the case of Pu, Pu(V) is reduced to Pu(III) and forms highly symmetric tridentate sorption complexes at the {111} faces of magnetite (Kirsch et al., 2011). In the case of neo-forming magnetite from Pu-containing solutions, Pu(III) also becomes partly trapped by magnetite but is increasingly released by Fe-enforced aging (Dumas et al., 2019). The mechanism behind this seems to be the large size of Pu(III) ions, which is largely incompatible with the magnetite structure; incorporation proceeds through the formation of pyrochlore-like islands in 2-nm magnetite nanoparticles, which are energetically less favorable than Pu(III) sorption complexes and, hence, form via kinetic entrapment rather than by being thermodynamically stable. Analogous to the Tc case, control between sorption and incorporation might be highly relevant for Pu mobility and, hence, safety.

In general, the geometry, speciation, and oxidation state of the outermost surface layer of magnetite is a function of the redox potential and pH (Katheras et al., 2024). These aspects of the magnetite surface chemistry should be considered when interpreting the sorption of radionuclides.

### 2.3.2.2 Experimental methods to distinguish between total Fe and redox-active Fe

Fe can play an important role in the reduction and oxidation of redox-sensitive radionuclides and can, thus, influence their migration and retention. Whether electron transfer between Fe and radionuclides occurs depends on various factors: i) the form in which Fe is present (e.g., redox state, aqueous, adsorbed, structural, coordination, and mineral properties), ii) the type and redox state of the radionuclide, and iii) chemical conditions in the
solution (pH, ionic strength, etc.). Measuring total Fe, or total Fe(II) and total Fe(III), does not necessarily lead to correct estimates of the actual redox interaction between Fe and the radionuclide. Examples are the incomplete reduction of Tc(VII) by structural Fe in various clay minerals (Bishop et al., 2011), the incomplete reduction of Se(IV) by Fe(II) sorbed to clay (Charlet et al., 2007), and the absence of Se(VI) reduction by pyrite (e.g., Bruggeman and Maes, 2010; Charlet et al., 2012).

The redox properties of both sorbed Fe and structural Fe have been investigated in a variety of Fe minerals under different conditions (pH, $E_h$). Fe minerals subjected to study include Fe sulfides, siderite, magnetite, Fe associated with clay minerals, and iron (oxyhydr-oxides. Several methods have been used to distinguish between total Fe and Fe accessible to redox reactions, and are explained in the following paragraphs.

A first method is radionuclide-mineral batch experiments. The extent to which the Fe in a mineral is redox-active is assessed by comparing the composition of a solution and solid state before and after the batch experiment using spectroscopic measurements such as XAS and Mössbauer spectroscopy or using (sequential) mineral dissolution methods. In these batch experiments, often small amounts of radionuclide are used; therefore, the method does not necessarily assess the total amount of redox-active Fe. However, it provides the most reliable information on whether the Fe mineral is reactive towards a specific radionuclide.

Examples of studies investigating the impact of redox-active Fe on the reduction/oxidation of radionuclides are:

* Fe associated with clay minerals and micas:

* Tc(VII) (e.g., Bishop et al., 2011; Jaisi et al., 2009; Jaisi et al., 2005; Peretyazhko et al., 2008; Yang et al., 2012).

* Np(V) (Schacherl et al., 2023)

* U(VI) (e.g., Burgos, 2016; Chakraborty et al., 2010; Ilton et al., 2010; Luan et al., 2014; Tsarev et al., 2016)

* Se(IV) (e.g., Charlet et al., 2012; Scheinost et al., 2008).

* Fe sulfides:

* Tc(VII) (e.g., Huo et al., 2017; Liu et al., 2008)

* U(VI) (e.g., Bruggeman and Maes, 2010; Livens et al., 2004).

* Se (IV, VI, -II) (e.g., Bruggeman et al., 2012; Charlet et al., 2012; Curti et al., 2013; Han et al., 2011; Han et al., 2012; Kang et al., 2011; Kang et al., 2013; Naveau et al., 2007; Scheinost and Charlet, 2008)

* Pu(VI) (e.g., Hixon et al., 2010)

* Magnetite:

* Tc(VII) and Np(V) (e.g., Cui and Eriksen, 1996)

* U(VI) (e.g., Grambow et al., 1996; Huber et al., 2012; Missana et al., 2003; Rovira et al., 2003)

* Se(IV) (e.g., Scheinost and Charlet, 2008)

* Se(VI) (Poulain et al., 2022)

* Pu(V) (e.g., Kirsch et al., 2011; Powell et al., 2005).

* Carbonates:

* U(VI) (e.g., Ithurbide et al., 2009; Ithurbide et al., 2010)

* Se(IV) (e.g., Badaut et al., 2012; Chakraborty et al., 2010; Scheinost and Charlet, 2008)

* Np(V) (e.g., Scheinost et al., 2016)

* Fe(oxyhydr)oxides:

* Tc(VII) (Peretyazhko et al., 2012; Um et al., 2011; Yalçintaş et al., 2016; Zachara et al., 2007)

o Pu(V) (e.g., Powell et al., 2005)

A second method is using probe compounds, where the redox activity of an Fe-containing mineral is assessed using ‘probe compounds.’ In batch experiments, a probe compound (oxidant or reductant) reacts with an Fe mineral, and the redox activity is assessed based on both the loss of the probing compound and the quantification of Fe(II) and Fe(III) in the mineral before and after reactions. Probe compounds used to assess the redox activity of Fe are often strong oxidants and reductants. The reduction and oxidation capacities that result from these reactions can be seen as a maximum redox activity and may not be entirely representative for reactions with specific radionuclides.

Examples of relevant research are:

* Clay minerals:

* Nitroaromatic compounds (e.g., Neumann et al., 2008; Neumann et al., 2011)

* Dithionite and $O_2$ or $H_2O_2$ (e.g., Stucki, 2011; Gorski et al., 2012a; Qian et al., 2023).

* Fe(oxyhydr)oxides:

* Ti-EDTA (Heron et al., 1994)

* Ascorbate (e.g., Roden, 2003)

A third method is based on biochemical processes (microbial reduction). Anaerobic microbes are often used to study the redox activity of, for example, structural Fe(III) in clay minerals. The Fe associated with a mineral is characterized before and after reduction using methods such as XAS or Mössbauer spectroscopy. These methods are limited to certain electrolyte conditions that are suitable for the microbes. Electron mediators, such as anthraquinone-2,6-disulfonate (AQDS), are often used to enhance the redox reaction occurring between a microbe and structural Fe. Many redox reactions in nature occur via microbial interaction. Therefore, these methods can give a good representation of which Fe is actually redox-active under natural conditions.

Examples of studies using microbial interaction to study redox-active Fe in minerals are:

* Clay minerals (e.g., Bishop et al., 2011; Dong et al., 2009; Dong et al., 2003; Jaisi et al., 2009; Jaisi et al., 2005; Kostka et al., 1999; Luan et al., 2014; Pentráková et al., 2013; Seabaugh et al., 2006; Stucki, 2011; Yang et al., 2012)

* Magnetite (e.g., Byrne et al., 2016; Dong et al., 2000)

* Fe(oxyhydr)oxides (e.g., Roden, 2003)

A fourth method assesses the redox activity of Fe by measuring the electron transfer (as electrical current) between an (Fe) mineral and a working electrode (possibly making use of mediators). An electric potential, reducing or oxidizing, high or low, is applied, and any resulting redox reactions with the Fe mineral are directly measured as electric current. This makes it possible to study both the extent and kinetics of a redox reaction. Measurements can be performed under well-defined conditions such as pH, ionic strength, and redox potential.

Examples of studies using this method to characterize the redox activity in Fe minerals are:
* Clay minerals (e.g., Gorski et al., 2012a; Gorski et al., 2012b; Hoving et al., 2017; Sander et al., 2015)

* Fe(oxyhydr)oxides (e.g., Aeppli et al., 2018)

* Fe sulfides (e.g., Hoving et al., 2017).

### 2.3.2.3 Characterization of redox reactions

Batch interaction experiments. In a batch interaction experiment, the $K_d$ (as a distribution solid-liquid ratio) formalism is often applied. However, it is not to be seen as being strictly “sorption.” The $K_d$ concept is based on reaction reversibility within a reactive transport conceptual model (see Chapter 3 on transport). When sorption is combined with a reduction reaction coupled with the oxidation of a mineral substrate (such as magnetite, pyrite, or mackinawite), the reduced RN species:

* has a much lower solubility as pure phase (e.g., U, Se, Sb, Tc, Mo); consequently, solid phases from the reduced (or oxidized) RN species can form;

* is however often present as single sorbed (e.g., U(IV), Tc(IV)) species (Chakraborty et al., 2010; Yalçıntaş et al., 2016);

* have a priori very limited ability to return to a solution, with the sorption reaction often being irreversible (Couture et al., 2015; Markelova, 2017).

Spectroscopic techniques. The determination and quantification of the redox state of elements and speciation in solution can be characterized using spectroscopic techniques by making use of special spectrochemical cells. Understanding the redox reactions of radionuclides on mineral surfaces is of importance to predict their mobility in the environment. For instance, process understanding on a molecular scale is the basis for geochemical modeling calculations that is applied in calculations for the safety assessment of nuclear waste repositories (Geckeis et al., 2013; Gorski et al., 2013). This understanding requires techniques such as XANES spectroscopy, which allow for the unambiguous identification of redox states because they are element-specific and sensitive to electronic and local structures. XAS covers almost all elements of the periodic table, but the penetration depth through water and mineral phases required for *in situ* redox speciation is given only for energies of a few keV; hence, the application of XAS is typically limited to certain elements, starting with the 3d transition metals. Since all atoms of a probed element contribute equally, and in a statistically representative way, for a measured XAS signal, this method is intrinsically suited to perform quantitative speciation analyses, without the need to account for redox-dependent absorption coefficients or the spectroscopically “silent” oxidation states of a given element. These advantages are offset by a relatively poor detection limit (approximately 1 ppm, depending on the element and matrix effects, often one or two orders of magnitude higher), and an intrinsically poor spectral resolution, especially for local structure analyses. In spite of these disadvantages, XAS is often the only spectroscopic method available for redox speciation, with the notable exception of Mössbauer spectroscopy, which is widely used for Fe redox speciation but limited to only a few other elements (Sn, Sb, Te, and possibly also I and Np).

Moreover, for actinides (An), the high-energy-resolution (HR-) XANES technique, at the An $M_{4,5}$-edge, provides a significant

advantage (Vitova et al., 2013) over the widely applied conventional An $L_3$-edge XANES mode, when investigating the oxidation states and electronic structures of An elements (Bahl et al., 2017). It shows high sensitivity to minor contributions from An oxidation states in mixtures, as it directly probes the An 5f unoccupied states, which are sensitive to changes in chemical bonding. Synchrotron-based X-ray methods in combination with optical spectroscopy ultimately provide information on the structure and bonding. As an example, UV-Vis spectroscopy, as a tool for molecular speciation, is sensitive to the coordination symmetry of An complexes (Meinrath, 1997).

In sorption studies with redox-sensitive An, such as Np and Pu, it is of utmost importance to achieve well-defined redox conditions. Many studies within this topic involve strict precautions to prevent the interference of minor amounts of oxygen, for example by working under an inert gas atmosphere (Banik et al., 2016; Elo et al., 2017; Marsac et al., 2015). In addition, electrochemical methods make it possible to control and, thus, systematically vary redox conditions. The application of an electrical field *in situ* can attenuate potential artefacts resulting from the oxidation and subsequent alteration of samples during sampling, storage, and transportation prior to experimental work (Pidchenko, 2016). Moreover, a mechanistic understanding of the interfacial reactions within the first minutes up to hours of contact can be obtained.

The methods of X-ray spectroscopy, optical spectroscopy, and electrochemistry can be combined by designing a spectro-electrochemical cell (Pidchenko, 2016). Specific questions that are so far unsolved could be investigated using this setup. Such questions include how an applied redox potential influences the surface speciation of redox-sensitive radionuclides, or what surface sites facilitate electron transfer processes (Banik et al., 2016; Marsac et al., 2015).

The redox state of elements in mineral phases can be determined/quantified using synchrotron methods such as high-energy-resolution fluorescence detection XANES (HERFD-XANES).

Owing to relatively long core-hole lifetimes, elements that can be probed by their K edges (up to and including period 5) allow for XANES measurements with a resolution sufficient to discriminate different oxidation states even when several oxidation states of an element are present in the same sample. This has been demonstrated, for instance, for the four common oxidation states of Se (-II, 0, IV, and VI) (Scheinost et al., 2008; Rojo et al., 2018). Actinides, however, are typically probed at L-edges with much shorter core-hole lifetimes and a correspondingly poor resolution (>5 eV). Furthermore, the evaluation of XANES is complicated in very diverse coordination environments, such as those with uranium changing from an yl-type at oxidation states VI and V to a highly symmetric cubic coordination for U(IV) (Hennig et al., 2007; Ikeda et al., 2007). Because of this, the more exotic pentavalent oxidation state of U has been intrinsically difficult to detect and can be easily overlooked in oxidation state mixtures. Owing to the advent of secondary monochromatization with, for instance, bent analyzer crystals in Rowland geometry, the spectral resolution can be significantly improved using HERFD-XANES (Kvashnina and Scheinost, 2016). However, a substantial breakthrough in the detection and quantification of various U oxidation states in
mixtures has been obtained only after applying the HERFD-XANES technique to M edges (Butorin et al., 2017; Butorin et al., 2016a; Butorin et al., 2016b; Prieur et al., 2018; Kvashnina et al., 2018; Kvashnina et al., 2017; Leinders et al., 2020). This method has also been applied to other actinides such as Pu, Am, and Np (Kvashnina et al., 2019; Epifano et al., 2019; Prieur et al., 2018; Schacherl et al., 2022a). Detection limits depend on the element investigated, the chemical composition of a sample, and the instrumental setup (including the synchrotron source) but can reach values in the order of ppm. Such state-of-the-art detection methods benefit from additional beamline developments such as cryogenic sample environments that allow reduced beam damage (Schacherl et al., 2022b).

In addition to improving the spectral resolution of XANES by employing analyzer crystals, substantial progress has been made in the statistical treatment of spectroscopic data in order to derive reliable speciation as well as to extract the pure spectra of endmember components from spectral mixtures (Rossberg et al., 2003; Rossberg et al., 2009), a method also applied to quantitatively derive the structures of the different oxidation states of Se and Tc (Yalçintaş et al., 2016; Rojo et al., 2018; Charlet et al., 2007). This method was further enhanced by coupling with Monte-Carlo simulations in order to derive the specific surface geometry of Pu(III) on magnetite (Rossberg and Scheinost, 2005; Kirsch et al., 2011). Finally, novel approaches can be used to directly derive a radial pair distribution function from EXAFS spectra, a method particularly well-suited for strongly disordered structures such as that of the Cm(III) aqua complex (Rossberg and Funke, 2010).

*Modeling of redox reactions.* Thermodynamic models are particularly relevant and required for a comprehensive prediction of An behavior (migration) in natural systems. Only few studies in the literature have attempted to quantitatively describe redox reactions in the presence of minerals, via mechanistic SCMs.

Soltermann (2014) described the sorption and oxidation of Fe(II) on Wyoming montmorillonite by using the 2SPN SC/CE model and considering the formation of ferric surface complexes (≡S$^{5}$OFe$^{2+}$) along with ferrous surface complexes (≡S$^{5}$OFe$^{+}$). The developed model is capable of describing the adsorption of Fe(II) on different montmorillonites; however, it does not provide a realistic mechanistic description, since it does not account for the intrinsic redox properties of a clay mineral itself.

Degueldre and Bolek (2009) discussed Pu sorption onto Al₂O₃, FeOOH, and SiO₂ colloids based on a simple R_d approach. Their calculations included consideration of the effect of redox potential, which was shown to have a major impact on the calculated Pu uptake data. Schwantes and Santschi (2010) proposed an SCM involving all Pu redox states (from +III to + VI) and minerals, as well as a reversible surface-mediated Pu(V)-Pu(IV) reaction, to explain the observed kinetics of Pu uptake. Experimentally determined redox potential (pe; E_h), oxygen partial pressures, or redox speciation are essential for the application of SCMs to describe sorption behavior of redox-sensitive elements. However, in most studies mentioned above, such data are often not provided owing to difficulties associated with their experimental determination (Altmaier et al., 2010; Schüring, 2000).

Marsac et al. (2015) investigated the uptake of Np(V) on illite under anaerobic conditions. Partial reduction of Np(V) to Np(IV) at the mineral surface was observed under conditions where only

Np(V) was expected to be stable in an aqueous solution. The results were explained using a TSM, based on measured redox potentials and the high stability of the inner-sphere surface complexes of tetravalent An (Np(IV)). The stability constants of such surface complexes can conveniently be estimated by analogy to other elements, for example to Th(IV), which is not redox-sensitive. This approach was also successful in the case of Pu adsorption on illite (Gorski et al., 2013; Banik et al., 2016). Based on such ideas, Banik et al. (2017) simulated the uptake of Np, specifically Np(IV) and Np(V), on illite even at high salt contents. The sorption behavior of pure oxidation states, which are required in the model, have been experimentally tested for Np(V) under oxidizing conditions and for Np(IV) under an inert gas atmosphere, by Nagasaki et al. (2017), Nagasaki et al. (2016). A model of Banik et al. (2017) was supported by these independent studies, while the mechanism of the redox reactions and the origin of the negative measured E_h-values remained unclear. To understand those findings, a SCM providing a more detailed molecular-scale description for the macroscopic quantification of the retention of redox-sensitive radionuclides under variable redox conditions should be developed. Starting from an approach by Marsac et al. (2015) and Banik et al. (2017) is possibly an appropriate option, since the modeling strategy was found to be successful for various actinides and cover a broad range of solution compositions.

*Kinetics.* Radionuclide redox kinetics in aqueous solutions depend on various parameters. These parameters include the redox potential defining the concentration of available electrons, the concentration of redox-active substances in the system under investigation, and the number of electrons that need to be transferred for a specific redox reaction. Other chemical specificities can, however, play a role. These specificities include chemical reactions as a consequence of both redox transitions (e.g., the kinetically hindered reduction of PuO₂⁺ to Pu(IV)) and the type of binding of a redox-sensitive radionuclide to an electron acceptor or donor (e.g., inner-sphere complexation or outer-sphere interaction), and, therefore, also geochemical boundary conditions. Therefore, redox kinetics require a detailed understanding of ongoing reactions. Numerous investigations have been published (e.g., Hixon and Powell, 2018; Ma et al., 2020; Romanchuk et al., 2016).

It is also worth noting that, in natural systems, Tc(VII) reduction is found to depend strongly on the Fe(II) content of rocks (Huber et al., 2017 and references therein). Rock treatment and conditioning prior to experimental studies, such as the preservation of the rock oxidation state, is consequently of utmost importance for mimicking *in situ* conditions. Reoxidation by contacting with oxygenated groundwater, simulating oxidized melt water into a granitic fracture, seems to depend on the nature of granite-bound Tc: Surface-precipitated TcO₂ seems to reoxidize faster than surface-sorbed Tc(IV). The kinetics of redox reactions have been measured for Tc(VII)/Tc(IV) pairs and correlate with the redox potential measured in solution (Kobayashi et al., 2013). Relatively rapid reduction is observed when the E_h is far below the Tc(VII)/Tc(IV) borderline. Experiments also revealed an acceleration of Tc(VII) reduction in the presence of Fe(II) solid phases, highlighting the possible relevance of specific surface interactions for electron transfer.

Cevirim-Papaioannou et al. (2018) found similar dependencies of U(VI) reduction to U(IV)O₂, in their solubility studies. Depending on the reducing agent they applied to adjust the
reducing conditions, they identified the $pe_{(exp)} - pe_{(U(VI)/U(IV))}$ difference as a driving force for reduction. This explains that, even when adding strongly reducing Sn(II) to a U(VI) solution at hyper-alkaline conditions, up to almost 2 years are required to reduce U(VI) quantitatively to U(IV)O₂, because the Eh difference is relatively low at high pH.

### 2.3.3 Solid solutions

#### 2.3.3.1 Theory and current understanding

Solid solutions are mixtures of two or more components on a molecular level existing over a range of chemical compositions. Solid solutions are ubiquitous in natural and industrial systems and play a major role in the fate and transport of elements in the earth's surface environments. Solid-solution formation can lead to a retention of harmful ions in the environment and is included in the design of engineered barriers, such as for the immobilization of radionuclides or industrially generated pollutants (Prieto et al., 2016). The disordered arrangement of different atoms over the same lattice is stabilized by the entropy, providing an energetic advantage compared to that of pure, stoichiometric endmembers. The charge difference and size of substituting atoms as well as the flexibility of a host structure, such as its Young's modulus, are important for determining the extent of solid-solution formation. The chemical composition of a solid solution can be defined in terms of fixed endmember compositions—for example, the binary solid solution $U_{1-x}Pu_xO_2$ can exist in a range of compositions specified by the mole fraction, that is, the $X_{PuO_2}$ of the endmember $PuO_2$. Solid-solution formation is typically favored between isostructural endmembers, where exchangeable atoms have the same charge and similar ionic radii. In many systems, the mixing is incomplete. In general, solid-solution formation is limited at low temperatures and is favored at higher temperatures (Bosbach et al., 2020).

A special case, relevant to nuclear waste management, is the formation of solid solutions in contact with aqueous solutions. Depending on the conditions, solid solution-aqueous solution systems (SS-AS) can form via the re-crystallization of pre-existing minerals or via co-precipitation from supersaturated solutions (Klinkenberg et al., 2014). These processes can compete with classical sorption processes such as the adsorption or formation of a pure radionuclide phase. As a consequence of the uptake of radionuclides into solid solutions, these are structurally bound into a thermodynamically stable phase, sometimes leading to significantly lowered radionuclide solubility (Bosbach et al., 2020). For example, long-term laboratory experiments conducted over several years indicated that $^{226}$Ra uptake into (Ba,Ra)SO₄ may take place in time scales relevant to deep-geological-waste disposal (Klinkenberg et al., 2014).

##### 2.3.3.1.1 Solid-solution formation and radionuclide retention.
In the safety assessment of planned nuclear disposal sites, the co-precipitation of radionuclides with sulfates and carbonates is considered, as this provides more realistic aqueous concentrations of radionuclides in the pore water of the nearby fields of nuclear waste repositories (Curti, 1999; Curti et al., 2010). For instance, barite is widely studied in the field, as it is the main scavenger of $^{226}$Ra, a long-lived radionuclide (Brandt et al., 2015; Brandt et al., 2018; Klinkenberg et al., 2018). Under repository-

relevant conditions, Ba and Ra released from radioactive waste will react with porewater containing sulfate, triggering the formation of (Ba,Ra)SO₄ (Curti et al., 2010) or even (Ba,Sr,Ra)SO₄ (Klinkenberg et al., 2018). In the case of $^{226}$Ra, an ideal solid solution with BaSO₄ irreversibly binds a contaminant in a well-defined, stable phase. Depending on the amount and the properties of the solid, the formation of solid solutions can become a very important retention mechanism and is, therefore, included in some scenarios of safety cases for deep geological repositories (NAGRA, 2002).

##### 2.3.3.1.2 Thermodynamic and kinetic aspects of solid-solution formation.
The most common tool used to rationalize the possible equilibrium states of a solid solution is the Lippmann diagram (Glynn, 2000; Lippmann, 1980). The Lippmann diagram is a phase diagram in which the solidus is a function of the solid-phase mole fraction, while the solutus is a function of aqueous activity fractions. The Lippmann model introduces the concept of total activity product, ($\sum \Pi$). This is defined as the sum of the partial activity products contributed by each component (endmember) of a solid solution. Considering the Lippmann diagram of a (Ba,Ra)SO₄ solid solution (Figure 1) as an example, the solidus is expressed ($\sum \Pi$) as the total solubility product, using Equation 12, while the solutus is expressed ($\sum \Pi$) in terms of the aqueous phase composition given in Equation 13.

$$ \begin{aligned} \sum \Pi_{eq} &= K_{RaSO_4}a_{RaSO_4} + K_{BaSO_4}a_{BaSO_4} \\ &= K_{RaSO_4}\gamma_{RaSO_4}X_{Ra} + K_{BaSO_4}\gamma_{BaSO_4}X_{Ra} \end{aligned} \quad (12) $$

$$ \sum \Pi_{eq} = \frac{1}{\frac{X_{Ra,aq}}{K_{RaSO_4}\gamma_{RaSO_4}} + \frac{X_{Ba,aq}}{K_{BaSO_4}\gamma_{BaSO_4}}} \quad (13) $$

where $a_{Ba^{2+}}$, $a_{Ra^{2+}}$, and $a_{SO_4^{2-}}$ represent the activities in the free solution; $K_{BaSO_4}$ and $K_{RaSO_4}$, the solubility product of endmembers BaSO₄ and RaSO₄, respectively (with $K_{sp}(BaSO_4) = 10^{-9.97}$ (Hummel et al., 2002) and $K_{sp}(RaSO_4) = 10^{-10.27}$ (Langmuir and Riese, 1985)) at a temperature of 298.15 K; and $X_{Ba}$ and $X_{Ra}$, the molar fraction of BaSO₄ and RaSO₄ in the solid. $\gamma_{BaSO_4}$ and $\gamma_{RaSO_4}$ are the activity coefficients of the endmembers in the solid solution. $X_{Ra,aq}$ and $X_{Ba,aq}$ are the aqueous activity fractions.

The solid-solution composition corresponding to an equilibrium aqueous composition is determined by a horizontal line connecting the solutus and the solidus. For instance, an aqueous solution with a mole fraction of Ba in solution ($XB^{a+}_{aq}$) equal to 0.7 (Figure 1) corresponds to a solid-phase composition with a mole fraction of BaSO₄ ($X_{BaSO_4}$) of 0.3.

While the Lippmann diagram provides the basis for understanding the equilibrium dynamics between a solid solution and an aqueous solution, it does not consider the effect of supersaturation. The key feature of a solid-solution formation is that trace element contaminants can co-precipitate as a solid solution, even though the solution is undersaturated with respect to a pure endmember solid containing that element. The concept of supersaturation for solid solutions is extensively described in a review by Prieto (2009), and the common tools include the stoichiometric supersaturation function, ($\Omega_{st}$) (Prieto, 2009) and “$\delta$”-functions (Astilleros et al., 2003). As an example, the supersaturation function of a hypothetical solution containing
<figure><img src="image_2.png" /><figcaption>FIGURE 1
(A) Lippmann diagram showing the equilibrium states for (Ba,Ra)SO₄ for a solid solution-aqueous solution. (B) Supersaturation of (Ba,Ra)SO₄ as a function of the mole fraction of BaSO₄ (XBaSO₄). These diagrams were constructed considering a temperature of 298.15 K and a regular mixing model with an enthalpy of mixing between Ra and Ba equal to 2,947 J/mol (Vinograd et al., 2018).</figcaption></figure>

10⁻⁴ mol/L Ba(aq), 10⁻⁴ mol/L SO₄(aq), and 10⁻⁷ mol/L Ra(aq) was computed using Equation 14 and is depicted in Figure 1B.

$$ \Omega_{st}(X_{Ba}) = \frac{(a_{Ba}^{2+})^{X_{Ba}} (a_{Ra}^{2+})^{X_{Ra}} (a_{SO_4}^{2-})}{(K_{BaSO_4} \gamma_{BaSO_4} X_{Ba})^{X_{Ba}} (K_{RaSO_4} \gamma_{RaSO_4} X_{Ra})^{X_{Ra}}} \quad (14) $$

A supersaturation Ω equal to 1 denotes equilibrium (dashed line). Solid solutions with compositions below the dashed line are expected to dissolve, whereas solid solutions with compositions above the dashed line can precipitate with the maxima representing the thermodynamically favored solid-solution composition. The example shows that, although a solution may be undersaturated with respect to pure RaSO₄, radium present in trace amounts is likely to be removed from an aqueous solution through the formation of a BaSO₄-enriched solid solution.

Because of the complexity of sulfate mineral nucleation and precipitation kinetics, highly supersaturated solutions can lead to metastable phases (Poonoosamy et al., 2016; Prieto, 2014; Weber et al., 2017), and thermodynamic equilibrium calculations might not be sufficient for a reliable prediction of ²²⁶Ra mobility (Curti et al., 2019). While the mathematical description of the thermodynamic equilibrium behavior of solid solutions is well-established, the implementations of kinetics still remains a challenging task (Noguera et al., 2016). The effect of the kinetics of solid solution precipitation is accounted for in various works by considering the nucleation rate of the solid solutions (Noguera et al., 2016; Pina et al., 2000).

### 2.3.3.1.3 Nucleation of solid solutions in porous media.
The understanding of nucleation processes in porous media is essential for the selection and parameterization of mathematical equations implemented in numerical reactive transport tools where these processes are dominant. These equations are used to assess the geochemical long-term evolution of interfaces in deep geological nuclear waste repositories in which precipitation and co-precipitation processes are recurrent. The nucleation mechanisms of minerals are described by Kashchiev and van Rosmalen (2003). In

porous media, both homogeneous and heterogeneous nucleation can occur (Poonoosamy et al., 2016; Poonoosamy et al., 2020). Homogeneous nucleation occurs in a pore solution, whereas heterogeneous nucleation occurs on the surface of existing solids that provide active centers for nucleation. It is generally acknowledged that homogeneous and heterogeneous nucleation prevail at high and low supersaturation levels, respectively (Kashchiev and van Rosmalen, 2003; Kügler et al., 2016).

The Classical Nucleation theory (CNT) is commonly applied to describe nucleation processes in porous media or confined spaces (Poonoosamy et al., 2019; Prasianakis et al., 2017). Although the applicability of CNT is highly debated (Gebauer et al., 2014), recent studies support the applicability and validity of CNT for sulfate minerals (Curti et al., 2019; Prieto, 2014) and even carbonates (Henzler et al., 2018). Although CNT is based on empirical formulations, it allows for smooth integration of experimental observations on mineral precipitation and nucleation kinetics into reactive transport simulations (Poonoosamy et al., 2019).

The nucleation rate (J) depends on supersaturation and can be calculated as follows

$$ J = r \exp \left( - \frac{\Delta G_c}{kT} \right) \quad (15) $$

where k is the Boltzmann constant, T is the absolute temperature, and r is a pre-exponential factor related to solubility. ΔG_c is the energy required for the formation of a nucleus of critical size and is given as:

$$ \Delta G_c = \frac{\beta v_0^2 \sigma^3}{(kT \ln \Omega)^2} \quad (16) $$

where v₀ is the molecular volume in the solid phase, β is a geometry factor that depends on the shape of the nucleus, and σ [J/m²] is the surface tension of the nucleating phase and the solution interface.

For a solid solution, the stoichiometric supersaturation (Ωst) is considered in Equations 15, 16. Therefore the nucleation rate for a (Ba,Ra)SO₄ solid solution is written as a function of the composition:
$$J(X_{RaSO_4}) = r(X_{RaSO_4}) \exp\left(-\frac{\Delta G_c(X_{RaSO_4})}{KT}\right) \quad (17)$$

where $r(X_{RaSO_4})$ and $\Delta G_c$ are functions of the solid-solution composition (Pina et al., 2000).

Nucleation in porous media is highly debated (Stack, 2015). Several authors (Prieto, 2014; Putnis and Mauthe, 2001; Stack et al., 2014) have demonstrated the dependence of nucleation kinetics on the pore size, with an inhibition of nucleation occurring in smaller pores. Unlike a free bulk solution, a porous medium divides a fluid into a large number of pores, which act as independent volumes with respect to nucleation. Because the nucleation probability is essentially proportional to the solution volume, precipitation in small pores is most likely hindered compared to that in a bulk solution. This point is illustrated by the basic equation for induction time (Kashchiev and van Rosmalen, 2003), given as:

$$t_i = \frac{1}{JV} \quad (18)$$

where $t_i$ is the induction time, $J$ is the nucleation rate (the rate of appearance of supercritical nuclei per unit solution volume), and $V$ is the volume of the fluid.

In the same line of thoughts, other authors (Emmanuel and Berkowitz, 2007) consider the solubility of minerals to be pore-size-dependent. This is commonly referred as the pore-size controlled solubility (PCS) effect. Under PCS consideration, the effective solubility of a mineral in a given pore volume, $S_e$, is related to its bulk solubility, $S_0$, as follows:

$$\frac{S_e}{S_0} = \exp\left(-\frac{V_s^c}{mRT} \frac{2\sigma \cos \theta}{r}\right) \quad (19)$$

where $V_s^c$ is the molar volume of the salt, $m$ is the number of ions per unit volume, $R$ is the gas constant, $T$ is the absolute temperature, $\sigma$ is the surface tension, $r$ is the radius of a cylindrical pore, and $\theta$ is the contact angle between the nucleating phase and the wall of the cylindrical pore. In most cases, $\theta > 90^\circ$, and $\exp(-V_s^c/mRT \frac{2\sigma \cos \theta}{r}) > 1$, such that the effective solubility of a mineral in a given pore volume is always higher than the solubility in bulk water, that is, $S_e > S_0$. PCS, shown in Equation 20, considers that a given solution can be supersaturated in a big pore and undersaturated in a small pore, which can give rise to a heterogeneous mineral precipitation distribution with preferential mineralization in larger pores.

Although many studies found an inhibition of nucleation in small pores, according to the latest review on mineral crystallization in confinement (Meldrum and O'Shaughnessy, 2020), it is not possible to make this a generalization. Indeed, other studies have shown that nucleation occurs uniformly across all pore sizes, without size dependence (Borgia et al., 2012), or even preferentially in fine capillaries (Hedges and Whitelam, 2012). Crystallization in confinement is still an open question and requires further investigation, as it affects the pore network alteration and transport properties of host-rock.

In recent work (Poonoosamy et al., 2021), the authors tested and validated theoretical models describing the crystallization of a sulfate solid solution against microfluidic experimental datasets. The investigation showed that i) CNT extended for solid solutions can be used to determine the composition of nucleating phases,

ii) nucleation kinetics can be a competitive growth mechanism, and iii) zonation phenomenon in solid-solution formation is governed by a complex interplay between the diffusion of reactants and the crystallization kinetics, as well as other factors such as surface tension and lattice mismatch.

### 2.3.3.2 Research methods

Typically, SS-AS systems are studied in batch precipitation, solubility, or uptake experiments. Complementary to the evolution of the chemical composition of a solution via wet chemical and radiochemical methods (e.g., ICP-MS, liquid scintillation counting, or gamma spectrometry), solids are characterized in-depth to understand the details of structural uptake (Brandt et al., 2015; Klinkenberg et al., 2018). A review on the characterization of SS-AS systems is provided by Prieto et al. (2016). Classical methods for characterization are spectroscopy, such as TRLFS, X-ray absorption (EXAFS and XANES), and Raman or IR spectroscopy. Lately, methods with high spatial resolution providing chemical information, including Secondary Ion Mass Spectrometry (SIMS), Micro-X-ray Fluorescence (μ-XRF), and electron microscopy (SEM-EDX or TEM-EDX), have been applied to unravel the distribution of minor components in host materials. X-ray and neutron diffraction have been applied to assess thermodynamic quantities such as the excess mixing volume.

Complementary, theoretical methods, including DFT (Vinograd et al., 2013) and MD (Heberling et al., 2014), have also successfully been applied to SS-AS systems.

Recently a novel microscopic method for the observation of dissolution/precipitation processes in tightly confined spaces at the nano-to-micrometer scale has been developed. It consists of a combination of microfluidic experiments (lab on a chip) and pore-scale modeling (Poonoosamy et al., 2019). Micromodels enable the in situ live monitoring of flow and reactive transport in 2D geometries and have been used, for example, to study single and multiphase flow, including the effects of pore geometry, capillary pressure, and fluid-phase saturation (Poonoosamy et al., 2019 and references therein). Combined with spectroscopic techniques and pore-scale modeling, micromodels enable the spatiotemporal visualization of mineralogical changes and ion concentrations, resulting in supersaturation determination.

### 2.3.3.3 Remaining uncertainties

Uncertainties remain with respect to the competition between processes in pore spaces, for example, the adsorption of $^{226}$Ra within a clay rock versus the formation of secondary phases. Some clay pore waters (e.g., Callovo-Oxfordian Clay and Opalinus Clay) contain significant amounts of Sr and sulfate, such that minor amounts of Ba and $^{226}$Ra may lead to the formation of a (Ba,Sr,Ra)SO$_4$ solid solution.

Despite the extensive studies of the thermodynamic properties of (Ba,Sr,Ra)SO$_4$ solid solutions, summarized by Vinograd et al. (2013), Vinograd et al. (2018), only a few studies (Pina et al., 2000) have been dedicated to assessing the nucleation of solid solutions or were limited to computational investigations (Noguera et al., 2016). The effect of pore size dependency on nucleation remains unknown.
## 2.4 Upscaling – "bottom-up approach"

As will be discussed in the chapter on transport, describing retention in transport mainly requires using a distribution coefficient ($K_d$), where the entire complexity of radionuclide/water/rock interactions is reduced to a single empirical value. Because of its simplicity, it is still the most used approach in transport modeling. However, sorption is a very complex process and heavily depends on the conditions. Therefore, the predictive capabilities of the transport calculations are only valid if the sorption ($K_d$) value is determined under similar (chemical/mineralogical) conditions.

For this reason, sorption models on pure clay minerals based on a mechanistic approach used to understand and quantify the processes controlling the uptake of radionuclides were developed. Recently, Baeyens and Bradbury (2017) and Bradbury and Baeyens (2017) developed so-called thermodynamic sorption databases for montmorillonite and illite, respectively. Such databases contain mineral-specific characteristics such as site types, site capacities, protolysis constants, selectivity coefficients, and surface complexation constants. These sorption models can be incorporated into reactive transport models that couple the influence of chemical interactions with the physical transport processes that take place in a water-saturated porous medium.

A commonly used approach to describe sorption in complex systems is the Component Additivity or Bottom-Up approach. This approach is based on the hypothesis that the uptake of radionuclides in complex mineral/groundwater systems can be quantitatively predicted from the knowledge and understanding of the mechanistic sorption processes on single minerals, and models developed to describe them.

The “bottom-up” approach, as studied in the EURAD WP5-FUTURE, considers three pathways (schemes) in order to address open questions: i) *from dispersed to compacted solids* (Are sorption data measured in dispersed systems applicable for compacted systems?), ii) *from single minerals to clay rocks* (Can we describe the sorption and transport behavior of radionuclides in mineral assemblages/clay host rocks by adding up effects described by a model developed for the main mineral components?), and iii) *from single radionuclide to multi-component chemical systems* (How does the presence of competing elements influence the mobility of the radionuclide under study?). In this project, bottom-up schemes were applied to three main groups of elements/radionuclides: weakly sorbing (anionic Se-species), moderately sorbing (Ra and Ba), and strongly sorbing (transition metals and lanthanides/actinides), and they were applied on a selection of clay phases (pure minerals: illite, smectites, kaolinite, and their mixtures (clay rocks), namely Opalinus Clay, Callovo-Oxfordian Clay, and Boda Claystone Formation (Maes et al., 2024).

Increasingly more literature have used the sum of major components to describe the behavior of a mineral assemblage. Some examples are given below. Bradbury and Baeyens (2011b) studied the sorption of Cs, Co, Ni, Eu, Th, and U on MX-80 bentonite and Opalinus Clay and obtained satisfactory predictions in applying the bottom-up approach. Marques Fernandes et al. (2015) tested the bottom-up approach for similar radionuclides on Opalinus Clay and Boda Claystone rocks. Montavon et al. (2022) successfully used this approach to describe uranium retention in Callovo-Oxfordian Clay.

The aspect of competition was extensively studied by Marques Fernandes and Baeyens (2020) for illite and montmorillonite. Orucoglu et al. (2022) investigated the competition between Pb²⁺, Co²⁺, Zn²⁺, and Mg²⁺ on montmorillonite edge surfaces. Sorption competition was investigated for different systems (see different contributions in (Maes et al., 2024)) and found to be a very important aspect for describing or understanding experimental observations from radionuclide sorption experiments. Often, deviations from a model find their origin in the presence of competing elements that are not taken into account in the model (insufficiently defined aqueous environment); hence, it is imperative to constrain as detailed as possible the aqueous system in a model.

Results obtained in this project provide convincing arguments that bottom-up schemes are applicable for the description of sorption processes, albeit not always straightforward, as relevant underlying mechanisms and sorption competition aspects must be taken into account rigorously. The transferability of data and models has been established from dispersed to compacted systems and from single to complex mineral assembly.

## 2.5 Sources of sorption data

This section presents an extended (but not complete) overview of available sources of sorption data, mainly $K_d$ based, relevant for clay and crystalline host rocks for safety assessment purposes.

### 2.5.1 Sorption on 2:1 clay minerals

Since 1995, the Paul Scherrer Institute has been conducting sorption studies on montmorillonite and developing a semi-quantitative mechanistic non-electrostatic sorption model for the uptake of (radio)contaminants. The studies started with detailed investigations of Ni and Zn sorption on montmorillonite (Baeyens and Bradbury, 1997; Bradbury and Baeyens, 1997b) and continued with a wide group of elements in Bradbury and Baeyens (2005b). The same basic non-electrostatic sorption model could be implemented for the clay mineral illite, and similar to montmorillonite, a number of (radio)contaminants were modeled (Bradbury and Baeyens, 2009a; Bradbury and Baeyens, 2009b). These studies formed an important basis for the development of sorption databases used in major Swiss safety assessments. Comprehensive compilations of the sorption data and modeling can be found in Baeyens and Bradbury (2017), Bradbury and Baeyens (2017), Baeyens and Marques Fernandes (2018), and Baeyens and Marques Fernandes (2018).

In the context of the European Project EURAD WP5-FUTURE ((Maes et al., 2024), final technical report), new sorption data on several 2:1 minerals have become available: Ra and Ba on illite and montmorillonite for the development of a sorption model (Klinkenberg et al., 2021; Marques Fernandes et al., 2023); Ra and Ba on illite and nontronite (sorption model, transfer batch to compacted state, contributions of Garcia-Guttierez et al. in Maes et al., 2024); Cd, Co, Zn, andNi on illite/Febex/(kaolinite)/mixtures (sorption model—additivity, transfer batch to compacted state, contribution of Missana et al. in Maes et al., 2024); Zn sorption reversibility on illite and montmorillonite (contribution of Dähn et al. in Maes et al., 2024); and $UO_2^{2+}$ on illite ($PO_4^{2-}$ competition for sorption, contribution of Del Nero et al. in Maes et al., 2024)).
## 2.5.2 Sorption on argillaceous rocks

### 2.5.2.1 Bentonite

To support the Swiss safety case, dedicated sorption studies have been performed on MX-80 bentonite. A study was specifically designed for the Swiss repository design, and the data are published in a technical report (Bradbury and Baeyens, 2011a). The bottom-up approach was also tested on selected datasets (Bradbury and Baeyens, 2011b).

A study by Grambow et al. (2006) presents sorption data for 11 different radionuclides (Cs, Ni, Pb, Eu(III), Am(III), Cm(III), Ac(III), Tc(IV), Th, Zr, and U(IV)) on MX-80 bentonite. The experimental sorption datasets were modeled, and the model uncertainty was assessed.

Comprehensive sorption studies have been carried out on FEBEX clay (a Spanish bentonite with a 98% smectite content), exhaustively characterized in the frame of the Spanish nuclear waste safety case, within homonymous projects (Huertas et al., 2000) and other international projects. Fundamental chemical parameters of FEBEX bentonite were analyzed to support geochemical modeling (Fernández et al., 2004). The radionuclides whose sorption was analyzed include anions such as selenite (Missana et al., 2009; Mayordomo et al., 2016) and cations of different valence: Cs (Missana et al., 2004; Missana et al., 2014), Ca, Co, Sr (Missana and Garcia-Gutierrez, 2007; Missana et al., 2008b; Missana et al., 2009; Mayordomo et al., 2019), Ga (Benedicto et al., 2014a; Missana et al., 2009; Benedicto et al., 2014b), Pu (Begg et al., 2015), Np (Benedicto et al., 2014a), U (Missana et al., 2004; Missana et al., 2009), and other elements considered as actinide analogues (Galunin et al., 2009; Galunin et al., 2010). Radionuclide sorption on FEBEX smectite was mainly described based on cation exchange and surface complexation, which were interpreted considering a two-site (weak, S<sub>w</sub>OH; strong, S<sub>s</sub>OH) NEM. The sorption irreversibility of some radionuclides onto FEBEX bentonite was also analyzed (Missana et al., 2004; Begg et al., 2015). Recently, new sorption data on FEBEX bentonite became available for selenite, Ra, Ba, Cd, Co, Zn, and Ni. Emphasis of the sorption study was on development of sorption models, transfer of sorption data from the dispersed to compacted state, sorption reversibility, and sorption competition (see contributions of Garcia Guttierez et al. and Missana et al. in the final technical report of the project FUTURE project (Maes et al., 2024)).

### 2.5.2.2 Opalinus Clay

Similar as in the case of bentonite for the Swiss safety case, sorption measurements have been carried out on Opalinus Clay and related formations. Reports on these studies were published by Lauber et al. (2000) and Baeyens et al. (2014a). The bottom-up approach was tested on a number of datasets (Bradbury and Baeyens, 2011b; Marques Fernandes et al., 2015)).

### 2.5.2.3 Boom Clay

Within the geological disposal program in Belgium (coordinated by ONDRAF/NIRAS), sorption studies have been performed by SCK CEN to complement radionuclide migration studies. Van Laer et al. (2016) gives an overview of all the sorption experiments performed on Boom Clay and the main constituting minerals illite and montmorillonite. These emphasize the influence of natural organic matter on radionuclide sorption behavior. Boom

Clay contains a considerable amount of organic matter that exhibits an important interaction (complexation and colloid formation) with a wide variety of radionuclides (transition metals, lanthanides, and actinides), changing its sorption behavior.

### 2.5.2.4 Boda Claystone

The geological disposal program in Hungary (coordinated by PURAM) considers two occurrences of the Boda Claystone Formation as potential sites (Breitner et al., 2015; Németh et al., 2016). In addition to diffusion investigations, sorption studies were performed by EK on Boda Claystone involving Cs, Sr, Co, and I (Mell et al., 2006a). Within the frame of a Swiss-Hungarian bilateral project, the bottom-up approach was tested for Cs, Co, Ni, Eu, Th, and U, in cooperation with (PSI) Paul Scherrer Institute (Marques Fernandes et al., 2015). Sorption was mostly related to illite-rich clayey matrices at the microscale (Osán et al., 2014; Kéri et al., 2016). The bottom-up approach based on illite was found to work well for Boda Claystone, despite the high hematite content due to the oxidized nature of the rock compared to the reducing properties of Opalinus Clay. Recently, new sorption data on samples of Boda Claystone became available for selenite and Ni. Sorption studies were complemented with diffusion studies, and emphasis was placed on the reversibility of sorption and competition effect (Ni vs. Co) (see contributions of Czömpöly et al. in (Maes et al., 2024)).

### 2.5.2.5 Callovo-Oxfordian argillite

Sorption studies on Callovo-Oxfordian argillite have been performed in several laboratories.

Studies for the following species have been reported: I⁻ (Montavon et al., 2014; Bazer-Bachi et al., 2006), SO₄²⁻ (Bazer-Bachi et al., 2007), Sr (Altmann et al., 2015), Cs (Chen et al., 2014b; Melkior et al., 2005), Ni (Chen et al., 2014a), Zn (Altmann et al., 2015), Pb (Orucoglu et al., 2018), Sn (Kedziorek et al., 2007; Latrille et al., 2006), Eu (Descostes et al., 2017), Pu (Latrille et al., 2006), and organic ligands and acids (EDTA, isosaccharinate (ISA), phthalate, oxalate, and polymaleic acid) (Dagnelie et al., 2018; Rasamimanana et al., 2017; Durce et al., 2014). Often, these studies were combined with migration experiments.

In Andra (2005), K<sub>d</sub> values on Callovo-Oxfordian argillite for safety relevant radionuclides are listed.

Recently, new sorption data on samples of Callovo-Oxfordian argillite became available for selenite, Ra, Ba, Eu, and U (see contributions of Garcia-Gutierrez et al., Hassan-Loni et al., and Montavon et al. in the final technical report of the FUTURE project (Maes et al., 2024)). Significant progress on the understanding of U retention in Callo-Oxfordian clay stone (Montavon et al., 2023; Montavon et al., 2022) was achieved.

## 2.5.3 Sorption on crystalline rocks and components

Compilations of sorption distribution coefficients (K<sub>d</sub> for Swedish (SKB) and Finnish (POSIVA) safety cases in crystalline bedrocks are available in Crawford (2010) and Hakanen et al. (2014). More recently, the sorption behavior of various elements has been intensively studied using rock samples provided by the Bukov Underground Research Lab (URL), Czech Republic), Onkalo (Finland), Forsmark (Sweden), and other sites with crystalline rock types potentially suitable for spent-fuel disposal.
Li et al. (2018), Li et al. (2020), and Puhakka et al. (2019) investigated the sorption behavior of Se(IV) on Grimsel granodiorite (mainly plagioclase, K-feldspar, quartz, and biotite), using artificial Grimsel groundwater and different concentrations of Se(IV). Fukushi et al. (2013) focused their study on Eu(III) sorption on granite. The distribution coefficients of Ba on Olkiluoto pegmatitic granite and veined gneiss, as well as those of Grimsel granodiorite, and their main minerals (quartz, plagioclase, K-feldspar, and biotite) were obtained using batch sorption experiments (Muuri et al., 2018; Muuri, 2019).

The sorption of Ni on fracture-filling calcite from the Bukov URL was studied and the resulting $K_d$ values at different pH levels and initial concentrations were compared. The results confirmed substantial Ni sorption, indicating its potential significance in radioactive waste repositories cf. (Parigi et al., 2022). The significant influence of pH on calcite stability and surface charge must be considered when providing $K_d$ values for performance evaluation calculations cf. (Al Mahrouqi et al., 2017).

(Fabritius et al., 2022) studied the sorption of Ra on biotite minerals and on whole rock samples using reference groundwater from the Finnish and Swedish spent nuclear fuel repositories at Olkiluoto and Forsmark, respectively. The rock types included veined gneiss, mica gneiss, diatexitic gneiss, granitic pegmatoid, granodiorite, monzogranite, and tonalite, and different water types were used, ranging from fresh to saline. The experimental and PHREEQC-interpreted distribution coefficient ($K_d$) results will be used in the safety case calculations of the Finnish and Swedish spent fuel repositories to further improve the understanding and reduce the uncertainties regarding the behavior of Ra in crystalline rocks. A multi-site SCM was developed to describe the Se sorption on biotite (Li et al., 2020).

Batch sorption experiments focused on rock-forming minerals, partially also present as fracture-filling minerals, including muscovite, orthoclase, and quartz. Various aspects of investigations comprise surface complexation modeling, TSM modeling, and data compilation for databases (e.g., Noseck et al., 2012; Noseck et al., 2018; Videnská et al., 2013; Stamberg et al., 2014; Videnska et al., 2015; Stockmann et al., 2017; Britz, 2018; Richter et al., 2016).

An application of the smart $K_d$ concept (see Section 3.1.3.3) in reactive transport calculations for porous sedimentary systems is described in detail by Noseck et al. (2012).

Newly formed sheet silicates and other precipitates in fractures, such as sulfates (e.g., baryte, celestite) or carbonates (e.g., calcite), show remarkable retention capacities for some radionuclides, such as Ra (Brandt et al., 2015; Chagneau et al., 2015; Klinkenberg et al., 2018; Curti, 1999), rare earth elements (REEs), and U(IV) (Drake et al., 2018).

A mechanistic explanation for heterogeneous sorption on nanotopographic muscovite surfaces due to the availability of energetically favorable sorption sites has been provided by Schabernack et al. (2023), who used kinetic Monte-Carlo (kMC) simulation techniques. The energetic differences of surface sorption sites available at nanotopographic structures such as steps, pits, and terraces were investigated, and the energies of adsorption and desorption reactions were obtained from DFT calculations.

### 2.5.4 Sorption databases

Publications of sorption databases ($K_d$ values) used in performance assessments are very sparse in the open literature, and they are mainly in the form of technical reports issued by waste management organizations (WMOs; e.g., ANDRA, https://www.andra.fr/sites/default/files/2018-04/dossier-options-surete-apres-fermeture_0.pdf). Examples are sorption databases on illite, bentonite, and argillaceous rocks that have been used in the Swiss safety cases for marl formations (Bradbury and Baeyens, 1997a), as well as those on bentonite (Bradbury and Baeyens, 2003b), illite (Bradbury and Baeyens, 2017), and Opalinus Clay and related formations (Bradbury and Baeyens, 2003a; Baeyens et al., 2014b).

RES³T (https://www.hzdr.de/db/res3t.login) is an acronym for the “Rossendorf Expert System for Surface and Sorption Thermodynamics” and is an open-source digitized thermodynamic sorption database comprising mineral-specific datasets for sorption modeling (e.g., SCPs, thermodynamic data, and surface complexation reactions). The database covers 145 minerals, approximately 140 sorbing ligands, 2,150 SSA measurements, 1,900 surface site datasets, and 6,700 surface complexation reactions, and it is based on literature references to ensure re-traceability and transparency.

## 3 Transport processes in host rocks

This chapter describes the current predictive models for the migration of radionuclide (RN) as dissolved species (colloid-related transport behavior and the transport of gases are not discussed here) in a fully saturated environment, the role and influence of structural properties, as well as the experimental and numerical approaches utilized for clay and crystalline host rocks. In this context, the main processes, such as diffusion, advection, and retention in clay and crystalline rocks, are characterized, including the main components and heterogeneities along the migration path and the understanding of coupled transport and retention processes.

### 3.1 Introduction to diffusion, advection, and retention processes

#### 3.1.1 Basics

For a more in-depth reading, we refer to several text books (Crank, 1975; Dullien, 1991; Mitchell, 1993; Grathwohl, 1998; Cussler, 2009; Appelo and Postma, 2004).

In porous (fully saturated) media, the net flux (J) of a solute is equal to the sum of all the fluxes induced by the different gradients that may be present in the considered system: chemical potential, hydraulic head, electrical potential, and temperature; this is described by Onsager reciprocal relations. However, the two most dominant transport processes for radionuclide migration are:

* Diffusion (created by a gradient of chemical potential and, hence, concentration), and;

* Advection (created by a gradient of hydraulic pressure).
Diffusive flux is the net movement of dissolved species present in an aqueous solution because of a difference of chemical potential (with the concentration gradient as driving force). Fick's first law describes the relationship between a diffusive flux and a concentration gradient

$$J = -D_e \frac{\partial C}{\partial x} \qquad (20)$$

where J is the solute flux [mol/m²s], $D_e$ is the effective diffusion coefficient [m²/s] of the solute in the porous medium, and C is the aqueous-phase solute concentration [mol/m³].

The effective diffusion coefficient ($D_e$) of each dissolved species directly depends on the diffusion coefficient in pure water, $D_0$, its diffusion-accessible porosity in medium $\eta$ [-], and geometric factor G [-] depending on the microstructure of the porous medium:

$$D_e = \eta G D_0 \qquad (21)$$

The geometric factor is a lumping parameter that averages the complex morphology of a pore structure and typically consists of two terms: tortuosity $\tau$ [-] and constrictivity $\delta$ [-]. (Please note that the G factor is sometimes defined as $\tau^2/\delta$ in papers cited in this report (e.g., in Appelo et al. (2010); Glaus et al. (2015a)). Tortuosity accounts for the effects of an increased path length of a molecule diffusing through the water-filled pores of a porous medium and is, therefore, larger than one. Constrictivity accounts for pore narrowing and widening but is often arbitrarily set equal to one. Both constrictivity and tortuosity cannot be measured independently. Only the lumping geometric factor G can be estimated from a diffusion experiment if the porosity is known:

$$G = \frac{\delta}{\tau^2} \qquad (22)$$

The molecular diffusion coefficient, $D_0$, of a dissolved species in pure water depends essentially on its chemical nature (e.g., atomic or molecular weight, hydrated radius, steric effects, and electrical charge), temperature, and the viscosity of the fluid. Typical values for the molecular diffusion coefficients of ions and molecules in water are in the range of $10^{-9}$ and $10^{-10}$ m²/s (Li and Gregory, 1974).

Fick's second law describes the transient evolution of flux:

$$\frac{\partial C}{\partial t} = \frac{D_e}{\eta} \frac{\partial^2 C}{\partial x^2} \qquad (23)$$

The radionuclide transport in a porous medium is influenced by retention processes at the water-rock interface (see previous chapter on retention processes). When these processes are rapid, linear, and reversible, their effect on the transport regime can be accounted for by a retardation factor, R [-], which can be related to the solid-liquid distribution coefficient $K_d$ [m³/kg] (often derived from sorption experiments, see Chapter 2 on retention processes):

$$R = 1 + \frac{\rho_d K_d}{\eta} \qquad (24)$$

where $\rho_d$ is the dry bulk density [kg/m³].

By rewriting Fick's second law to include reversible retardation, one obtains:

$$\frac{\partial C}{\partial t} = \frac{D_e}{(\eta + \rho_d K_d)} \frac{\partial^2 C}{\partial x^2} = \frac{D_e}{\eta R} \frac{\partial^2 C}{\partial x^2} = D_a \frac{\partial^2 C}{\partial x^2} \qquad (25)$$

Other frequently used diffusion parameters (note that the correct name of each diffusion coefficient must always be clearly and explicitly stated along with its definition to avoid confusion; simply referring to a “diffusion coefficient” or to “diffusivity” is insufficient and can lead to serious errors when transferring/using data) are the pore diffusion coefficient, $D_p$, and the apparent diffusion coefficient, $D_a$, which are related to each other in following way:

$$D_a = \frac{G}{R} D_0 \text{ with } GD_0 = D_p \text{ and } D_e = \eta GD_0 = \eta D_p = \eta R D_a \qquad (26)$$

The product $\eta R$ is also often referred to as the “rock capacity factor.” The rock capacity factor $\alpha$, relates the concentration in a porous medium, $C_{bulk}$, to the concentration in a solution, C, i.e., the sum of the amount of the species in a solution and in a solid normalized to the solution and solid volumes (Tournassat and Steefel, 2015). Hence the term refers to a “storage capacity” of a porous medium for a certain species.

$$C_{bulk} = \alpha C = \eta RC = (\eta + \rho_d K_d)C \qquad (27)$$

The transport of dissolved species by advection is the result of the movement of pore water under the effect of a hydraulic pressure gradient. The hydraulic pressure gradient causes the water to move out of a porous medium at an average velocity called the Darcy velocity, $V_{darcy}$ [m/s].

The Darcy law describes the flowrate (Q) of water forced through a porous medium of cross section S [m²] and length L [m] by a hydraulic gradient ($\Delta P/L$):

$$Q = K S \frac{\Delta P}{L} \qquad (28)$$

The Darcy velocity corresponds to the water flowrate Q divided by the cross-section S of the considered porous volume:

$$V_{darcy} = \frac{Q}{S} = K \frac{\Delta P}{L} \qquad (29)$$

where Q is the flowrate [m³/s], $\Delta P$ is the hydraulic pressure difference that is exerted over the length of the porous medium (expressed in m of water column height), and K is a proportionality constant called hydraulic conductivity [m/s]. While the Darcy velocity represents the velocity of water flowing out of a porous medium, the pore water velocity ($V_p$) takes into account the fact that, inside the porous medium, water flows only through the pore space. Both parameters are related via the water porosity or total porosity, $\eta$:

$$V_p = \frac{V_{darcy}}{\eta} \qquad (30)$$

Of note, as soon as an advective flow is present, diffusive flow also becomes affected by mechanical and hydrodynamic dispersion and becomes a dispersion coefficient ($D_{disp}$), which is a function of the advective velocity ($a$ is the longitudinal dispersivity [m]):

$$D_{disp,p} = D_p + aV_{darcy} \qquad (31)$$
Mechanical dispersion (hydrodynamic dispersion) results from ground/pore water moving at rates both greater and less than the average linear velocity (Fetter, 1992). This is caused by: 1) fluids moving faster through the center of the pores than along the edges, 2) fluids traveling shorter pathways and/or splitting or branching to the sides, and 3) fluids traveling faster through larger pores than through smaller pores.

Because solute-containing water does not travel at the same velocity everywhere, “mixing” occurs along flow paths. This mixing is called mechanical dispersion and results in a distribution of solutes at the advancing edge of flow. The dispersion coefficient in the transport equation accounts for the combined effects of mechanical dispersion and molecular diffusion, both of which cause spreading from highly concentrated areas towards less concentrated areas.

The ratio of the advective transport rate over the diffusive transport rate is known as the Péclet number, $P_e$:

$$P_e = \frac{LV}{D} \qquad (32)$$

where L is the characteristic length [m] (travelled length of the fluid), and V and D represent the water velocity and the diffusion coefficient, respectively. The Péclet number is used as an indicator of the dominant transport mechanism. Equation 32 is the general form of the equation for the Péclet number, and depending on the underlying assumptions in the solute transport equation, it can have different forms (Huysmans and Dassargues, 2004).

The advection-diffusion transport equation reads.

$$\frac{\delta C}{\delta t} = \frac{D_{disp,p}}{R} \frac{\partial^2 C}{\partial x^2} - \frac{V_{darcy}}{\eta R} \frac{\delta C}{\delta x} \qquad (33)$$

Taking into account that $V_{darcy}/\eta = V_p$, and introducing the notation of “apparent” coefficients to incorporate the retardation effect, we obtain:

$$\frac{\delta C}{\delta t} = D_{disp,a} \frac{\partial^2 C}{\partial x^2} - V_a \frac{\delta C}{\delta x} \qquad (34)$$

The classical diffusion-retardation equation considers the concentration of solutes to be the same over the entire pore space and accounts for linear reversible sorption processes. The latter are mainly linked to the main sorption mechanisms in clays, i.e., ion exchange and surface complexation. The equation disregards i) the presence of EDLs that has an effect on the concentration profiles of ions within pores, and ii) non-reversible retention processes (e.g., irreversible chemisorption, surface precipitation, incorporation, solid-solution formations, and precipitation).

## 3.1.2 Numerical modeling of coupled physical and chemical processes in transport

Numerical modeling of coupled physical and chemical processes of fluid flow, solute transport, and chemical reactions is vital for predicting the migration of RN in host rocks. In general, such numerical models include reactive-transport equations (Steefel and Lasaga, 1994) for both solute transport and chemical reactions, as well as the conservation of momentum for fluid flow (Steefel et al., 2005). A multitude of reactive transport codes are available, see

Steefel et al. (2014) for a non-exhaustive list of reactive transport codes such as PHREEQC, PHREEQC-PHAST, HPx, HYTEC, CRUNCHFLOW, ORCHESTRA, PFLOTRAN, OpenGeoSys, TOUGHREACT, and GEM2MT.

### 3.1.2.1 Reactive multi-species diffusion models (mainly applicable to clays)

For each species, diffusion is coupled to speciation reactions in solution, sorption reactions on mineral surfaces, as well as mineral reactions.

One-dimensional reactive multi-species diffusion in solution may be described for i species in solution with concentration $C_i$ as

$$\frac{\partial C_i}{\partial t} = D_p \frac{\partial^2 C_i}{\partial x^2} - R_i - r_i, \quad i = 1, 2, \dots, n_{aqueous} \qquad (35)$$

where $D_p$ is the pore diffusion coefficient [m²/s]; $R_i$ [mol/m³s] denotes a source/sink term due to equilibrium reactions between the solution and non-moving species sorbed on mineral surfaces according to

$$R_i = \sum_{j=1}^{n_{sorbed}} v_{ij} \frac{1 - \eta}{\eta} \rho_s \frac{\partial S_j}{\partial t}, \quad i = 1, 2, \dots, n_{aqueous} \qquad (36)$$

where $n_{aqueous}$ and $n_{sorbed}$ [-] denote the number of species in solution and sorbed, respectively, and $v_{ij}$ [-] are stoichiometric coefficients. $r_i$ [mol/m³s] is a source/sink terms due to mineral reactions.

The mass action law or the assumed equilibrium chemistry in solution and between solutions and both minerals and mineral surfaces yields constraints that define the source/sink term for each species and represents the mechanistic understanding of surface and/or ion-exchange reactions of the following form

$$\equiv \mathbf{s}\mathbf{c}_1 + \mathbf{c}_2 \Leftrightarrow \equiv \mathbf{s}\mathbf{c}_2 + \mathbf{c}_1 \qquad (37)$$

Mineral reactions can be kinetically or by chemical equilibrium controlled and have the following form

$$r_i = \sum_{m=1}^{N_m} v_{im} r_m, \quad i = 1, 2, \dots, n_{aqueous} \qquad (38)$$

where $r_m$ [mol/m³s] is the net rate of precipitation ($r_m > 0$) or dissolution ($r_m < 0$) of a mineral per unit volume of rock, and $v_{im}$ [-] is the number of moles of species i per mole of mineral m present in the rock.

This can be further extended by accounting for effects of electrical potentials that develop during the movement of charged species (Nernst-Planck formalisms).

### 3.1.2.2 Reactive transport including advective fluid flow (important for crystalline rocks)

At the pore scale, crystalline rock can be considered as a binary domain consisting of a solid material and the pore space. The transport of (solute) species occurs in pore spaces and reactive processes such as the adsorption/desorption of RN, and mineral dissolution/precipitation take place at fluid-solid interfaces. The Navier-Stokes equation (NS) or the Stokes equation is employed as a momentum equation to describe fluid flow in a pore space (Li et al., 2008; Molins, 2015). Interfacial reactions are usually treated as a boundary condition at fluid-solid boundaries (Molins, 2015).
Pore-scale simulation equations.

(1) Fluid flow using NS equation:

$$ \frac{\partial \mathbf{u}}{\partial t} + (\mathbf{u} \cdot \nabla)\mathbf{u} = -\frac{1}{\rho}\nabla p + \nu\Delta\mathbf{u}, \quad (39) \\ \nabla \cdot \mathbf{u} = 0, $$

where $\mathbf{u}$ is a velocity vector, $p$ is pressure, $\rho$ is the fluid density, and $\nu$ is the kinematic viscosity of the fluid.

(2) Solute transport and chemical reactions; reactive transport equations:

$$ \frac{\partial c_i}{\partial t} + \nabla \cdot (\mathbf{u}c_i - D_i \nabla c_i) = R_i, \quad (40) $$

where $c_i$ is the molar concentration of species $i$ in solution, $D_i$ is the dispersion/diffusion coefficient, and $R_i$ is the reaction term.

(3) Surface reactions:

$$ -D_i \nabla c_i = \xi_{im} r_m, \quad (41) $$

where $r_m$ denotes the surface reaction rate and $\xi_{im}$ is the stoichiometric coefficient of component $i$ in each surface reaction $m$.

At the continuum scale, a porous medium is typically conceptualized as a continuum with bulk parameters that characterize its physical and chemical properties. At this scale, mobile and immobile phases are assumed to coexist at each point in space and time. The Darcy's equation is applied as a momentum equation, while reactive transport is described by the mass balance equation of each species (Steefel and Lasaga, 1994). Reactive transport in fractured porous media can be simulated using the coupled Stokes-Brinkman equation and reactive-transport equations (Yuan et al., 2016; Yuan and Zhang, 2019), which allows for efficiently modeling the alterations of rock properties caused by chemical reactions:

$$ \mu K_{perm}^{-1} \mathbf{u} + \nabla p - \mu^* \Delta \mathbf{u} = 0 \quad (42) \\ \nabla \cdot \mathbf{u} = 0, $$

where $K_{perm}^{-1}$ is a permeability tensor, $\mathbf{u}$ represents the physical velocity of the fluid in the fractures and the Darcy velocity in porous media, $\mu$ is the dynamic viscosity of the fluid, and $\mu^*$ is the effective dynamic viscosity.

### 3.1.3 Accounting for retention processes in transport: Empirical and thermodynamic approaches

There are different approaches for describing sorption and desorption processes at water-rock interfaces:

a) Empirical relationships, such as e.g., Langmuir, Freundlich, or linear sorption isotherms (Kd), and

b) Conceptual TSMs that describe ion-surface associations based on mechanistic approaches to characterize coordinative properties of surface complexes using, for example, the DDLM or more sophisticated model theories (Ochs et al., 2012).

Approaches from Langmuir and Freundlich, as well as TSMs, assume a thermodynamic equilibrium between a solid and surrounding electrolytes. However, compared to empirical approaches, mechanistic TSMs are based on the law of mass action.

#### 3.1.3.1 Empirical – constant $K_d$ approach

In performance and safety assessment (PA/SA) calculations, sorption is usually described using linear distribution coefficients (also see Section 2.1) constant in space and time, i.e., simple $K_d$ values (NAGRA, 2002; POSIVA, 2012b; SKB, 2011), and applied in the classical advection-diffusion equation as a retardation factor, $R$.

Since the sorption of many elements and radionuclides strongly depends on the geochemical properties of the surrounding solution and sediments, the approach of constant $K_d$ values is justified only for constant hydrogeological and geochemical conditions expected in various parts of a GDS.

For transport simulations for PA/SA purposes, typically, databases with experimental $K_d$ values (and uncertainty ranges) are compiled, and the databases are often extended for different conditions, geological strata, and other factors.

Within large time frames or in dynamic systems close to an engineered barrier system, a variability of geochemical conditions can be expected due to spatial heterogeneities and environmental variations. Therefore, it may be desirable to describe the transport and/or sorption of potential contaminants as a function of important environmental parameters such as pH or ionic strength.

#### 3.1.3.2 Process understanding-Thermodynamic approaches

Reactive transport simulations typically rely on transport codes coupled with a geochemical code, capable of using TSMs based on mass action law equations:

$$ \equiv sc_1 + c_2 \Leftrightarrow \equiv sc_2 + c_1 \quad (43) $$

These coupled models are useful for acquiring process and system understanding, as they are able to account for different chemical conditions, in contrast to the constant $K_d$ approach.

Different TSMs have been elaborated on to describe the mineral-fluid interface during sorption processes (Ochs et al., 2012) (also see Section 2.2.2). TSMs that are used to calculate defined surface complexation processes of, for example, metal cations (Me) are called SCMs. All TSMs are based on the law of mass action and mole balance equations (Davis and Kent, 1990a), but they involve different descriptions of EDLs, i.e., electrostatic interaction term (Davis and Kent, 1990a; Westall and Hohl, 1980). In non-EDL models (noEDLMs) electrostatic interaction is not considered; thus, a DDL is neglected. Each TSM (except noEDLMs) considers different electrostatic planes that can contribute to sorption processes.

For use in PA/SA, when it comes to simulating the impact of changing geochemical conditions over very long time frames and/or under variability in composition over larger distances, a fully coupled reactive transport approach may result in unfeasible computation times owing to the complexity of mechanistic model approaches (Noseck et al., 2014). However, it enables the complementing of sorption databases by underpinning a selected
value range, as well as the binding of sorption values for current/
changing conditions or alternative scenarios.

For application in safety assessment (Noseck et al., 2012; Noseck
et al., 2018), the “smart K<sub>d</sub> concept” (http://www.smartkd-concept.
de/) has been proposed and developed, and it has been implemented
in transport modeling. It combines state-of-the-art mechanistic
TSM theories with the traditional empirical, constant, linear, and
reversible K<sub>d</sub> approach and presents a computation-time-efficient
approach to simulate long-term contaminant transport over large
areas and long distances.

It allows for the interpolation of distribution coefficients
(referred to as smart K<sub>d</sub> values) for any observed (or hypothetic)
combination of environmental boundary conditions, treating
complex systems in a more rigorous way. This approach is based
on SCMs and considers further ion exchange, solid-phase equilibria,
and complexation in aqueous solutions (Stockmann et al., 2017).

An application of the smart K<sub>d</sub> concept in reactive transport
calculations for porous sedimentary systems is described in detail by
Noseck et al. (2012).

3.2 Transport in clays

3.2.1 Structural aspects of clays

Argillaceous sediments are deposited in various marine
(Opalinus, Callovo-Oxfordian, and Boom Clay) and lacustrine
(Boda Claystone) environments. Because of water-rock
interactions and microbial activity, physico-chemical changes
occur in clay sediments after their deposition, and the sediments
undergo diagenesis by being progressively transformed in clay rocks.
The diagenesis of clay sediments depends on many factors, such as
the depositional environment, types of clay minerals, amount of
organic matter, pore water chemistry, geothermal gradient, and
burial history. Increasing burial depth leads to sediment
compaction and chemical diagenesis, which results in a
systematic decrease of porosity, water expulsion, and progressive
transformation reaction of smectite to illite (Velde and Meunier,
2008; Weaver, 1989).

The principal minerals in argillaceous sediments are clay
minerals, with additional quartz, feldspars, and calcite, small
amounts of pyrite, Fe-carbonates, and apatite. In addition, solid
or immobile organic matter in the form of more or less evolved
kerogen may be present. When argillaceous sediments are
sufficiently compacted, such that they are no longer soft or
plastic but still not fissile (i.e., not cleavable along weak planes as
shales or metamorphic slates), they are called clay rock or mudstone.
Depending on its mineralogy (smectite content), a clay sediment
may exhibit swelling properties in the presence of water; this is
important in applications for nuclear waste disposal, as it provides
self-sealing properties. Argillaceous sediments contain variable
amounts of water present both in the clay mineral structures and
in the nanometer pore spaces. High water contents (>20 vol%) result
in soft and plastic clay (e.g., Boom Clay), whereas clay rocks/
mudstones (Opalinus Clay, Callovo-Oxfordian Clay, and Boda
Claystone) are drier (<20 vol% pore water) and stiffer. Further
compaction and water expulsion from the pores result in further
alignment of clay mineral grains such that the rocks become fissile,
forming shales and, eventually, slates.

Clay minerals are phyllosilicates (from the Ancient Greek word
φύλλον "phyllon", meaning leaf), and the principal clay minerals in
argillaceous host rocks are smectite (or "swelling clay"), illite,
kaolinite, and chlorite. These phyllosilicates are composed of
jointed octahedral (O) layers (such as that containing Al³⁺, Mg²⁺,
or Fe²⁺) and tetrahedral (T) layers (those containing Si⁴⁺ or Al³⁺)
(Figure 2). Illite and montmorillonite (main smectite) are composed
of TOT layer sequences, while kaolinite is composed of TO
sequences. Chlorites are TOT layer sequences with brucite-like
interlayers with compositions near to (Mg,Fe,Al)₃(OH)₆.

Due to isomorphic substitution where a higher-valency element
is replaced by a lower-valency element (Si⁴⁺ by Al³⁺ in tetrahedral
layers, or Al³⁺ by Mg²⁺ and Fe²⁺ in octahedral layers), permanent
electrical charges are formed. These negative charges become
charge-balanced by (hydrated) cations in the interlayers between
two TOT layers. Due to the higher negative charge of illite TOT
layers, cations (mainly K⁺) are fixed strongly in interlayers without a
hydration shell. In contrast, montmorillonite, which has a lower
charge, has interlayer cations that maintain their hydration shells.
The TO layers of kaolinite are uncharged. An excellent overview of
the surface properties of clay minerals is given by Tournassat et al.
(2015). For an extensive overview on the properties of clays, the
reader is referred to Mitchell (1993), Velde and Meunier (2008), and
Meunier and Fradin (2005).

3.2.2 Main processes controlling transport in clays

3.2.2.1 Diffusion-dominated transport

In compact clay systems (natural host rocks, and clay-based
engineered barriers) used for the geological disposal of radioactive
waste, diffusion is the dominant transport mechanism owing to the
very low hydraulic conductivities of these rocks (<10<sup>-11</sup> m/s).

Clay host rocks may experience tectonic disturbances, but due to
the presence of swelling clays and water, they tend to self-seal
relatively quickly (within weeks to years). The higher the amount
of water and smectites, the prompter is the self-sealing response
(plastic Boom Clay vs. Callovo-Oxfordian Clay and Opalinus Clay
host rocks). After tectonic disturbances, the intrinsic permeability of
fractured clay rock is expected to return to that of an intact non-
fractured rock (OECD/NEA, 1996).

Advective transport in clays due to the presence or creation of
faults and fractures (via natural processes or excavation processes) is
considered negligible in terms of geological time scales because of
self-sealing properties (Bernier et al., 2007; Bock et al., 2010;
Li, 2013).

3.2.2.2 Retention predominantly via ion exchange and
surface complexation

Clay minerals such as montmorillonite and illite (including
illite-smectite mixed layers) are important in the disposal of
radioactive waste, since the former is the major clay mineral
component of near-field bentonite (buffer and backfill) and the
latter is present at significant levels in many argillaceous host rock
formations under consideration.

The main sorption mechanisms of aqueous species on clay
minerals are cation exchange and surface complexation (see
Section 2.1 (Classical retention processes) for more details). The
sorption behavior of a broad series of elements has been well
investigated.
<figure><img src="image_3.png" /><figcaption>FIGURE 2
Representation of the structures of different types of clay minerals (Source, https://www.soils4teachers.org/mineralogy, with permission from SSSA).</figcaption></figure>

In addition to these classical retention processes, redox-induced precipitation is also an important retention mechanism for redox-sensitive radionuclides (see Section 2.2.2 for more details).

For transition metals, another retention mechanism has been observed in clays. At higher loadings of transition metals, the surface precipitation of neo-formed phyllosilicates or the formation of LDHs may occur (also see Section 2.3.1).

The presence of carbonate and sulfate minerals in natural clays may also lead to solid-solution formation, causing strong retention for some radionuclides (e.g., Ba, Sr, and Ra) (see Section 2.3.3, Solid solutions).

### 3.2.2.3 Special effect of negatively charged clay surfaces on radionuclide transport behavior – double layer effects

Because of negatively charged clay mineral surfaces, anionic, cationic, and neutral species exhibit distinct diffusion behavior.

* Anionic species become repelled by the negatively charged surface of clay minerals, which limits their accessibility to certain pores (smallest size; access restriction by constrictivity effects) and reduces their movement to a more limited zone of bigger pores. This phenomenon is known as “anion exclusion” or “Donnan exclusion.” The higher the negative charge, the stronger is the decrease in anion-accessible porosity (see the contribution of Van Laer et al. in (Maes et al., 2024)).

* Cationic species, on the contrary, become attracted by negative clay surfaces, leading to interaction processes that enhance the sorption/uptake in host rocks. These interaction processes are generally expected to lead to a slowed down diffusive transport and are called retardation processes.

* Neutral species are hardly affected by the electrostatic forces.

Recent studies (EC CATCLAY project, Altmann et al. (2015); and EURAD WP5 FUTURE project (Maes et al., 2024) have, however, revealed that the uptake of cations in clays not only leads to an increase in the rock capacity factor (cf. Equation 26) but may also have an impact on diffusive properties by altering $D_e$ values (potential explanations for this peculiarity will be outlined in the following section). For this reason, the effect of nuclide uptake by the solid phase does not inevitably lead to retardation in all cases but is rather the result of the combined impact on R and $D_e$ (cf. Equation 25).

The classical description of the diffusive processes of solutes interacting with surfaces according to Equation 25 is based on the assumption that species bound to a surface are immobile, viz. they can only change location after a desorption step (stationary and mobile phase concept borrowed from chromatography). Sorption, thus, contributes only as a “capacity factor” increasing the retardation of solutes. Once the sorption capacity is reached, or sorption sites are saturated, the amount of solute transported per unit cross section area and time, viz. the diffusive flux, is exactly the same as that for the case of an inert solute that does not interact with a clay surface. The only driving force for diffusion is the concentration gradient of dissolved species present in the bulk aqueous phase (cf. Equation 20, Fick’s first law). In such cases, $D_e$ values of different solutes would only depend on their respective $D_0$ values and on $\eta G$ factors (Equation 21). The latter factor is often assumed to be identical for all aqueous species in a given porous medium (e.g. (Aldaba et al., 2014; Jacops et al., 2017b),). However, this assumption is invalid for charged solutes in the nano-scale pore networks of negatively charged clay minerals and clay rocks (cf. Figure 3). For an overview of the literature, the reader is referred to reviews such as (Bourg et al., 2003; Miller and Wang, 2012;
<figure><img src="image_4.png" /><figcaption>**FIGURE 3**
Texture of a clay rock, illustrating the various water properties and the heterogeneity of the different types of clay minerals at the nanometer scale.</figcaption></figure>

Dark blue colors indicate water domains with bulk water properties, light blue colors indicate water domains, which are influenced by the negatively charged surfaces of clay minerals and their electrical double layers (EDLs). Depending on the amount of fixed charges of the latter, the overall charge equivalents of positively and negatively charged solutes may vary within the different water domains. The green color represents water domains in the interlayers of swelling clays. Of note, the fractional distribution of clay and water domains does not represent the reality.

Shackelford and Moore, 2013; Altmann et al., 2015; Tinnacher et al., 2016). The $D_e$ values derived for anionic species (Equation 17) are generally lower than those of uncharged reference species such as HTO or $D_2O$, while those of cationic species are generally higher, even after normalization for individual different $D_0$ values (see Table 4). In addition, for non-charged species (dissolved gasses), unexpected phenomena, such as a functional dependence of the “geometric factor” on the molecular size, has been demonstrated (Jacops et al., 2017a; Jacops et al., 2017b).

While potential differences in geometry factors between anionic and neutral species could be a valuable explanation for the lower $D_e$ values of anionic species in the sense of a restricted accessibility of anionic species to the pore space of negatively charged clay surfaces (called “anion exclusion” (Wigger and Van Loon, 2017)), such an explanation would fail for the enhanced $D_e$ values of cationic species. It has been demonstrated (Glaus et al., 2017) for different clay materials that the diffusion-accessible porosity of uncharged species is equal to the total water porosity. If the enhanced $D_e$ values were related to geometry effects, the only reasonable explanation would be found in the reduced tortuosity factors, $\tau^2$ (cf. Equation 22). Consequently one would have to assume shorter diffusion paths (Shackelford and Moore, 2013) for cationic species compared to neutral species, which cannot be justified from a physical point of view.

In view of the conceptual difficulties in explaining the different diffusive behaviors of differently charged species by geometry effects solely, it appears to be more sensible to question the assumption of the concentration gradients of bulk aqueous species being the only driving force for diffusion. In fact, it has been demonstrated for the

<table><caption>**TABLE 4** Normalization of the effective diffusion coefficient ($D_{e,i}$) of an aqueous species i with the effective diffusion coefficient ($D_{e,HTO}$) of tritiated (HTO) or deuterated (HDO) water and the respective diffusion coefficients in bulk water ($D_0$). Dependence on the electrical charge of the species.</caption><thead><tr><td>Dissolved species</td><td>D<sub>e,i</sub>/D<sub>e,HTO</sub>·D<sub>0,HTO</sub>/D<sub>0,i</sub></td></tr></thead><tbody><tr><td>Anions</td><td>&lt;1</td></tr><tr><td>HTO or HDO</td><td>= 1</td></tr><tr><td>Cations</td><td>&gt;1</td></tr></tbody></table>

diffusion of cationic species in compacted smectite clay that the direction and magnitude of observed diffusive fluxes depend on the concentration gradients of cationic species bound to the cation-exchange sites present in interlayer spaces, rather than those in the aqueous phase (Glaus et al., 2013).

The process-based physical models for such phenomena are still under development and should take into account structural characteristics of clay material, aside from those of diffusing solutes. Although no direct experimental evidence exists for a Brownian type of motion of surface-associated species, it is reasonable to assume that a similar relationship exists, according to Equation 20 (Fick's first law), between the concentration gradient of these species and their respective diffusive fluxes, as for the species in the aqueous phase. Conceptually, the overall observed diffusive flux can be described as being the result of the superposition of the individual fluxes in two phases, one representing a bulk aqueous species and the second representing a mobile surface species (cf. Equation 24).
Within the scope of such a concept, Fick's first law, shown in Equation 20, can be expanded to expressions such as

$$J = -D_e \frac{\partial C_f}{\partial x} - D_s \frac{\partial C_s}{\partial x} \qquad (44)$$

where the subscript f denotes species in the bulk such as the ("free") aqueous phase, and s denotes the mobile surface species, with $D_s$ being the respective effective diffusion coefficient of these surface species; hence, the term "surface diffusion" was developed (a summary of quantitative expressions related to $D_s$ from key references can be found in (Maes et al., 2017a; Maes et al., 2017b)). Depending on the particular characteristics of a given case, it is possible to simplify Equation 24, viz. by assuming (i) that either one of the terms on the right-hand side of Equation 24 can be neglected and (ii) by replacing the concentration of the sorbed species by using well-known thermodynamic equilibrium expressions between $C_s$ and $C_f$. A few examples from the literature are mentioned in the following sections. Experimental evidence for the validity of expressions such as Equation 24 had already been found many decades ago (e.g., (van Schaik et al., 1966)).

Glaus et al. (2007) showed that the diffusion across the bulk aqueous phase of cationic species bound to a clay mineral surface via cation exchange (e.g., Na⁺ and Sr²⁺) in compacted montmorillonite was negligible, and diffusion was dominated by cationic species present in the cation-exchange sites. Cation-exchange selectivity according to the Gaines-Thomas convention was used to relate $C_s$ and $C_f$. The diffusion of cationic species across compacted montmorillonite was, thus, described using the assumption of clay porosity being represented as a single phase. The same approach was applied by Birgersson and Karnland (2009), who used an electrostatic Donnan-type distribution between the clay phase and the bulk aqueous phase to describe the chemical equilibration of differently charged species. These authors could successfully describe the diffusive behavior of negatively and positively charged species in compacted smectite minerals and bentonites.

Gimmi and Kosakowski (2011) compiled a broad collection of literature data for the diffusion of cationic species in clay minerals, clay rocks, and argillaceous soils. They used a relationship, derived from Equation 24, between observed $D_e$ values and pertinent $R_d$ (or $K_d$) values to evaluate $D_s$ values. Therefore, they could demonstrate that the resulting $D_s$ values comprised rather limited parameter ranges for the cationic species tested in their study, confirming that Equation 24 is generically applicable.

All these studies were, however, mainly related to cationic species interacting with clay surfaces via cation exchange or via electrostatic interactions. Many representatives from the transition-metal and lanthanide or actinide series were shown to form strong surface complexes via chemical bonding with surface groups present at the so-called edge surfaces of clay minerals. Owing to the presumed covalent character of this type of surface bonding, it is hardly expected that such complexes exhibit an appreciable surface mobility. Recent studies by Glaus et al. (2015a), Glaus et al. (2020), Montoya et al. (2018), and Altmann et al. (2015), who investigated the diffusion of Co²⁺, Zn²⁺, Eu³⁺, and Cs⁺ in compacted illite as a function of a broad variety of solution parameters, such as pH, ionic strength, or the concentration of the respective stable isotope elements, also gave evidence for surface diffusion phenomena

being active in the case of these cations. Significantly larger $D_e$ values were indeed observed for these elements than could be expected from direct proportionalities with the corresponding $D_0$ values of an uncharged reference species. It was demonstrated that the characteristic parameter dependencies of the $D_e$ values were not related to the formation of surface complexes at the clay edge surfaces but rather to the species bound to planar surfaces. As an example (Glaus et al., 2015a), the diffusion depths in compacted illite of bi-valent transition metals significantly increased with decreasing ionic strength, while the total concentrations of these metals in clay at the interface between the clay and solutions were almost identical for all ionic strengths. This was a clear indication that the $D_e$ values were affected by the variation of ionic strength, rather than by sorption. It is well known that the uptake processes of cations at planar surfaces are susceptible to variations in the ionic strength, while the formation of surface complexes at edge surfaces is much less dependent on ionic strength (e.g., Bradbury and Baeyens, 2005b). It was further demonstrated that the $R_d$ values derived from diffusion experiments exhibited the same parameter dependencies as those measured in static batch sorption experiments, in which the tracer distribution between a liquid and the clay phase was measured using dilute clay suspensions.

In contrast to the very small pore sizes of smectites (order of 1 nm), the pore sizes of compacted illite may be larger by a factor of approximately 5–10 (Sanchez et al., 2008), and the planar surfaces may, thus, be in contact with bulk-type water. The distribution of charged species near such charged surfaces is frequently described using the concept of an EDL (Figure 4).

The neutralization of negative surface charges occurs via two different types of cationic species. The first type is located directly at the clay surface and is assumed to be bound as an inner-sphere complex; thus, it is considered as immobile. These surface species are frequently referred to as Stern layer species. The second type is a diffuse swarm of cations, near the surface, and assumed to be mobile. Anions, accompanied by an excess of diffuse swarm cations, may also enter the EDL domain to a reduced extent. Depending on the charge density of the solid phase and various solution parameters, cationic species are, thus, enriched in EDLs compared to their concentrations in bulk aqueous phase. If a concentration gradient parallel to planar surfaces exists, the diffusive flux of cationic species in a charged clay will, thus, be larger compared to that in uncharged clay, owing to the higher population of cations in the EDL and the resulting larger concentration gradients parallel to the clay surface. Similar to cation exchange in the interlayers of smectites, the unexpectedly large fluxes of transition metal cations can be understood as an enrichment phenomenon at the planar illite surfaces.

For quantitative purposes and for simplifying the numerical procedures for diffusion calculations, the Donnan approach has been proposed to represent the local distribution of the electrical potential near planar surfaces (Appelo et al., 2010; Appelo and Wersin, 2007; Altmann et al., 2015; Glaus et al., 2015a; Glaus et al., 2020; Gimmi and Altmann, 2018). Figure 4 shows a schematic representation of the equilibrium distribution of cationic and anionic species between free pore water and EDL water for a Donnan approach using a parallel clay pore geometry.

Similarly, this concept was used for a circular pore geometry (Glaus et al., 2015a; Glaus et al., 2020). Notably, the choice of the
<figure><img src="image_5.png" /><figcaption>FIGURE 4
Schematic representation of the Donnan approach for an electrical double layer (EDL) at planar clay surfaces. Species in blue (ions in free water) and green (ion swarm in the Donnan volume) are assumed to be mobile, while red cationic species are present as immobile surface complexes in the Stern layer. If a concentration difference of a tracer cation species between the left and right side of a clay is present, the individual flux contributions ($J_i$) would result from the local concentration gradients parallel to the clay surface in the different domains.</figcaption></figure>

model geometry is a rather arbitrary and an uncritical decision. In the Donnan approach, the extension of an EDL is finitely defined by the Donnan thickness, $d_{DL}$, exhibiting a constant potential (the Donnan potential, $\psi_D$) across the Donnan volume. The equilibrium distribution of cationic and anionic species between free pore water and water in the Donnan volume is governed by $\psi_D$, which is basically derived from the electrical potential at the surface, taking into account the shielding effect of cations in the Stern layer. $\psi_D$ and $d_{DL}$ are adjustable parameters in this model.

The total flux for each species $i$ ($J_{tot,i}$) is the sum of the fluxes in the free pore water ($J_{free,i}$) and the Donnan layer ($J_{DL,i}$) (Appelo et al., 2010; Altmann et al., 2015; Glaus et al., 2015a) (in the case of smectite-type clays with accessible interlayers, an extra flux, $J_{IL}$, is present and needs to be added to the total flux:

$J_{IL,i} = -f_{IL}\eta G_{IL} \frac{C_{IL,CEC}}{z_i} \frac{\partial \beta_i}{\partial x} D_0$, where $\beta_i$ is the molar or equivalent fraction, and $C_{IL,CEC}$ is the CEC expressed in interlayer water (mol/L), Appelo et al., 2010):

$$J_{free,i} = -f_{free}\eta G \frac{\partial C_i}{\partial x} D_0 \quad (45)$$

$$J_{DL,i} = -f_{DL}\eta G \frac{C_{DL,i}}{C_i} q_n \frac{\partial C_i}{\partial x} D_0 \quad (46)$$

$$J_{tot,i} = -\eta G \left( f_{free} + f_{DL} \frac{C_{DL,i}}{C_i} q_n \right) \frac{\partial C_i}{\partial x} D_0 = -D_{e,i} \frac{\partial C_i}{\partial x} \quad (47)$$

According to the assumed mobility of the cationic species (bulk aqueous and diffuse layer species are mobile; Stern layer species are immobile), the $D_{e,i}$ values for the situation depicted in Figure 4 were calculated as follows:

$$D_{e,i} = \eta G \left( f_{free} + f_{DL} \frac{C_{DL,i}}{C_i} q_n \right) \quad (48)$$

where $\eta$ [-] is the total water-accessible porosity, G [-] is the geometry factor, $f_{free}$ [-] is the porosity fraction containing free pore water, $f_{DL}$ [-] is the fraction of the Donnan layer porosity, $D_0$ [m²/s] is the diffusion coefficient in bulk water, and $q_n$ [-] is the viscosity factor. The ratio $c_{DL}/c$ represents the concentration ratio of cationic species between the diffuse layer and the bulk aqueous phase. $q_n$ reflects changes in the surface mobility compared to the mobility in the bulk pore aqueous phase. The application of expressions such as Equation 28 in transport models is limited to situations in which the various parameters are independent at

concentration c. This can reasonably be assumed to be the case in many examples of tracer diffusion in electrolyte media.

Within the concepts of such electrostatic approaches for the description of surface uptake reactions and diffusion processes, not only the behavior of counter-ions, but also those of the cations (anions, in the case of negatively charged clay surfaces), are implicitly included. A valuable overview of different models used to describe the behavior of anions in charged clays can be found in Tournassat and Appelo (2011). The picture shown in Figure 4 should not be confused with the anion exclusion concept, in which an assumption is made that part of the porosity is completely inaccessible to anions (e.g., Van Loon et al., 2007). The volume fraction of free pore water in the Donnan approach is always larger than the equivalent exclusion volume in the latter approach. The results of both approaches can, however, be mutually converted for the purpose of comparison.

The description of the diffusion of charged species in charged porous argillaceous media may, thus, be based on the use of lump-defined effective diffusion coefficients, comprising the characteristic geometric properties of the porous medium and the presence of multiple chemical concentration gradients, despite the comprehension that it may actually be the result of the superposition of the diffusion of individual chemical species in multiple phases. Because the local concentration gradients of these species cannot always be measured using experimental methods, the use of thermodynamic equilibrium models for the description of lump-defined effective diffusion coefficients is indispensable for this purpose. Currently, reactive transport codes as PHREEQC (Appelo and Wersin, 2007), CRUNCHFLOW (Steefel and Tournassat, 2021), ORCHESTRA (Jenni et al., 2021), and FLOTRAN (Krejci et al., 2021) have implemented double layer effects on diffusion in their codes, using different approaches for the spatial dependence of surface potential. Evidently, the diffusion model represented in Equation 28 uses an idealized picture of the chemical species and equilibria involved. This also applies for assumptions regarding geometry. The model does not describe the effects of the orientation of clay particles with respect to the direction of diffusion. In the idealized picture, the concentration gradients responsible for tracer diffusion are assumed to be parallel to the clay surfaces, while the respective concentration gradients perpendicular to the surface determine the chemical species equilibration at those surfaces and, thus, the chemical enrichment or depletion factors.
A broad variety of influencing factors have an effect on the parenthetical expression in Equation 28. These can be roughly subdivided into factors affecting geometric aspects and those affecting the chemical effects of a model. Among the latter, the composition of an aqueous solution in equilibrium with the clay phase may have a direct influence on the ratio of concentration in the diffuse layer and in the bulk aqueous phase ($c_{DL,i}/c_i$). For example, it has been demonstrated that the electrolyte concentration, the presence of other metal cations competing with the cation of interest, or pH may directly affect this distribution ratio (Glaus et al., 2015a; Glaus et al., 2020; Glaus et al., 2021; Glaus et al., 2015b; Glaus et al., 2020; Glaus et al., 2021). A recent study within the frame of the EURAD WP 5 FUTURE work package has further shown that the cationic composition of a clay surface produces similar effects. The effective diffusivity of transition-metal cations, such as $Co^{2+}$ and $Zn^{2+}$, depends distinctly on the cation loading, in compacted illite, with the highest diffusivities being measured for Li- and Na-illite, while the K- and Cs-illite forms show lower values (contribution of Zerva et al. in (Maes et al., 2024)). This behavior can be interpreted in terms of a different shielding of the fixed negative lattice charges, which leads to different electrostatic effects related to changes in the extension of an electrical double layer.

Alternative assumptions for the equilibrium distribution of cations between the diffuse layer and the bulk aqueous phase can be made, such as a Poisson-type distribution instead of a discrete Donnan-type distribution. Such assumptions result in similar diffusion models (cf. e.g. (Ochs et al., 2001; Tachi and Yotsuji, 2014)). Furthermore, it has been suggested that Stern layer species may exhibit surface diffusion effects (Krejci et al., 2021; Krejci et al., 2023). The various ambiguities in the model assumptions can only be resolved to a limited extent by macroscopic diffusion experiments. Ambiguities remain, for instance, with respect to the individual species concentrations and their molecular mobilities. More elaborated methods, such as spectroscopic, diffractive, and microscopic techniques (e.g., Strawn and Sparks, 1999; Bostick et al., 2002; Park et al., 2008; Malikova et al., 2010; Ferrage, 2016) are required to bound critical model assumptions. This is also the case for the use of atomistic simulation methods, which can give further insights into the detailed structural characteristics and processes related to the surface reactivity of clays (Kosakowski et al., 2008; Cygan et al., 2009; Bourg et al., 2017; Bourg and Sposito, 2010; Churakov and Gimmi, 2011; Holmboe and Bourg, 2014; Tinnacher et al., 2016; Tournassat et al., 2016a; Churakov, 2013; Gimmi and Churakov, 2019).

The broad variety of model concepts for describing the dynamics of diffusive transport processes in clay minerals and clay rocks render the definition of an effective diffusion coefficient increasingly more complex. In fact, a variety of different terminologies can be found in the literature, for the exact definition of effective diffusion coefficients, which sometimes makes it difficult to compare the results of individual works. Valuable attempts to provide an overview of this situation can, however, also be found in the literature (e.g., Shackelford and Daniel, 1991; Shackelford and Moore, 2013). Despite the "unphysical" nature of such coefficients, they are frequently applied for engineering purposes, because concentration gradients in the bulk aqueous phase are the only sound basis for transport

simulations. In order to attain an improved process-based understanding of the diffusion of charged species in charged porous argillaceous media, it is, however, imperative to improve the physical basis of the understanding of surface-induced diffusion phenomena.

### 3.2.2.4 Effect of particle size and orientation on diffusion

Clay mineral particles exhibit an anisometric shape associated with a wide size range (from tens of namometers to tens of micrometers).

Regarding the influence of particle size on diffusion properties, decreasing the particle size leads to an increase in the ratio between external (basal and lateral surfaces) and internal (interlayer) surfaces (Reinholdt et al., 2013). For non-swelling clay minerals, a decrease in particle size, in turn, leads to an increase in double layer effects (cf. 3.2.2.3), which affects the diffusion properties of cationic and anionic species. For swelling low-charged clay minerals such as smectite, this effect can be enhanced by the osmotic swelling of particles and the partial or total delamination of individual smectite layers, depending on the ionic strength and interlayer cation type. In this case, low ionic strength values and monovalent interlayer cations can promote the osmotic swelling process, leading to an increase in the Donnan volumes associated with anion exclusion and an increase in the retardation process for cation diffusion (cf. 3.2.2.3). In the specific case of water diffusion, the role played by external surfaces of clay particles on the overall diffusion process is limited, in line with the limited and very local influence of the electrostatic force between clay interfaces and water molecules (Bourg and Sposito, 2011). Indeed, Tertre et al. (2018) (among others) showed similar effective diffusion coefficients for water for porous media with the same porosity and made of uncharged kaolinite or charged illite clay particles. Furthermore, these authors demonstrated the crucial role of the osmotic swelling process in the water diffusion of swelling clay materials. In the case of porous media made of highly charged swelling vermiculite particles, without osmotic swelling, the effective diffusion coefficients, $D_e$, for water were similar to those obtained for porous media made of kaolinite or illite particles at a given interparticle porosity, irrespective of the particle size fractions used (i.e., 0.1–0.2, 1–2, or 10–20 µm). In the specific case of vermiculite samples, reductions in the water molecule mobility in interlayer spaces have, thus, been considered as dead volume for overall water diffusion and confirmed using various simulation approaches (Tyagi et al., 2013; Asaad et al., 2021; Ferrage et al., 2023). Furthermore, the situation is very different in porous media made of low-charge smectite, which exhibit osmotic swelling, for which the definition of internal (interlayer) and external (surface-related) surfaces becomes rather vague. In this case, the swelling of interlayer volumes, occupying a large part or all of the pore space, coupled with a reduction in interlayer water mobilities, generally leads to an overall decrease in the effective water diffusion coefficient (Bourg and Tournassat, 2015).

Owing to the anisometric particle shapes of clay minerals (Asaad et al., 2022; Bardot et al., 1998; Reinholdt et al., 2013; Tournassat et al., 2003; Sayed Hassan et al., 2006)), their mutual arrangement in porous media is most often associated with the development of anisotropy in pore network organizations. This effect most often leads to contrasting diffusion properties between the normal and in-