Migration of Contaminants in Fractured-Porous Media in the Presence of Colloids: Effects of Kinetic Interactions

The presence of colloids in groundwater can enhance contaminant migration by reducing retardation effects. The field experiments performed in a water-conducting feature of a shear zone at the Grimsel Test Site in Switzerland indicate that the sorption processes for contaminants on mobile and immobile colloids are kinetically controlled and that colloid filtration occurs. To enable the modelling of those experiments, an appropriate model of colloid-facilitated contaminant transport in fractured-porous media is developed. The physical system is modelled as a single planar fracture with adjacent fully saturated porous rock matrix. Contaminants can diffuse into the rock matrix but colloids cannot. In the mathematical model, the 1D advective contaminant transport along the fracture is coupled with contaminant diffusion into the rock matrix perpendicular to the fracture. Radioactive decay and sorption processes for contaminants in the rock matrix (linear equilibrium), on the fracture surface (linear equilibrium and linear kinetic reactions) and on mobile and filtered colloids (linear kinetic approach), are considered. The model for colloid transport includes filtration of colloids in the fracture and their remobilization. A useful approach is developed that can be applied to adequately describe a natural system (crystalline rock) with our double-porosity model (single fracture integrated into porous rock). Numerical solutions are obtained using an implicit finite difference scheme and realized in the transport code COFRAME. The transport code is validated by the comparison of calculated results with field experiment data. Some sets of simulations are performed to study the effect of kinetics for interaction processes.


Introduction
Among the host rocks which are discussed as options for the disposal of radioactive waste in deep geological repositories, crystalline rock formations such as granite or granodiorite have been considered. Crystalline formations provide a high mechanical stability, advantageous reducing hydrochemical conditions and a high sorption capacity for cations due to the occurrence of minerals as, e.g. feldspars or mica. However, in the envisaged repository depth (500-1000 m below surface), granite is fractured. In order to ensure the safe containment of the waste, especially in the case of high-level radioactive waste, the following barriers are considered in designing final repositories (SKB 2011): -In the KBS-3 concept, the waste is packed into copper containers to provide a high corrosion resistance. -As buffer and backfill material, dry compacted bentonite is favoured. In contact with water, bentonite swells and thereby creates a very low permeable diffusion barrier delaying radionuclide release into the fracture. In addition, bentonite provides a high sorption capacity for many radionuclides.
However, colloids (clay particles) might be released from the bentonite block when bentonite comes into contact with water (Schäfer et al. 2012;Baik et al. 2007;Missana et al. 2003;Pusch 1999;Gardiner et al. 2001). These colloids can be transported further through the fracture system by the groundwater. Such a scenario is particularly considered at the end of an ice age during glacial retreat, when glacial melt waters with very low ionic strength and high pH values (>9) might intrude into the formation down to repository depth (SKB 2011). Such water properties encourage the formation of stable colloids. But natural groundwater also contains colloids which are generated by weathering processes, conversion of organic substances, dissolution and precipitation processes as well as hydrochemical and hydraulic changes in the environment.
Colloids are defined as particles with diameters between approximately 1 and 1000 nm. The surface properties of colloids dominate compared to the bulk properties, i.e. colloids have an enormous specific sorption capacity. Consequently, colloids can sorb significant amounts of contaminants, which are usually rather immobile due to strong sorption on solid material. Because of their charge and size, colloids and colloid-bound contaminants can be transported without sorption, tend to linger in the middle of the fracture, where the fluid velocity at laminar flow is maximum, and might not diffuse into the rock matrix. Hence, the average colloid velocity can be higher than the average fluid velocity. As a consequence, the transport of contaminants adsorbed on such colloids might be significantly increased. For these reasons, the impact of colloids has to be taken into account in the assessment of the long-term safety of a repository for radioactive waste.
The fact that colloids can act as carrier for contaminants has been shown in numerous laboratory and field tests (Bates et al. 1992;Kersting et al. 1999;Puls and Powell 1992;Huber et al. 2012;Smith 1993;Kim et al. 1996;Baik and Hahn 1997). Field experiments performed within the Colloid Formation and Migration (CFM) project under near-natural flow conditions at the Grimsel Test Site (GTS) in Switzerland indicate that desorption processes for contaminants from mobile and immobile colloids are kinetically controlled (Blechschmidt and Lanyon 2008). This is true for the fluid transport times of several hours to several days that have been applied in these dipole experiments. These experiments further show that the colloid filtration proceeds according to a linear kinetic approach. At least in one of the dipole experiments, colloid remobilization seems to be important. In order to adequately describe these experiments, an appropriate model of colloid-facilitated contaminant transport is required.
The colloid-facilitated contaminant transport within fractured-porous media with kinetically controlled interaction processes has been previously modelled using different assumptions and approaches. In Ibaraki and Sudicky (1995), a very comprehensive model of coupled transport of radionuclides and colloids in discretely fractured-porous media is presented. According to that model, contaminants can be transported by the water both in the fracture and in the porous matrix, whereas the colloid transport occurs only in the fracture. The equilibrium as well as kinetic approaches for sorption of contaminants is considered. The interaction of colloids with the fracture surface is described by a linear kinetic reaction neglecting the colloid remobilization. The radioactive decay is taken into account, however without consideration of nuclide chains. In Li et al. (2004), three cases of colloid transport in a single planar fracture surrounded by saturated porous rock have been modelled based on different assumptions and solved analytically: (1) colloid transport only within the fracture considering an irreversible filtration onto the fracture surface (without remobilization), (2) colloid transport only within the fracture considering the mechanism of filtration as well as remobilization and (3) colloid transport including diffusion of colloids into the rock matrix and irreversible filtration onto the fracture surface. The transport of radionuclides is modelled under the assumption that the sorption reactions of radionuclides on colloids and solid material are in equilibrium. The radioactive decay is included in the model, however disregarding the contribution coming from parent nuclides.
In the model developed in this work, the key processes which affect the contaminant transport in a fractured-porous medium in the presence of colloids are treated explicitly. Linear kinetic reactions were used for describing the sorption processes of contaminants on colloids. The interaction of colloids with the fracture surface is modelled by a linear kinetic reaction including both filtration and remobilization terms. The justification and limitation of the approaches implemented in the model to describe the interaction processes are discussed. Radioactive decay is considered comprehensively, taking into account the complete nuclide chains. Based on the classical double-porosity model, a fractured-porous media are modelled as a singular planar fracture, which is situated in water-saturated porous rock. However, fractures in a natural crystalline formation largely vary in size, properties and structure. They are connected with each other and can form large-scale networks, especially in the case of a large-scale model area. In this work, a useful approach is developed that can be applied to adequately describe a natural system (crystalline rock) with our double-porosity model (single fracture integrated into porous rock). The kinetic approaches used in the model for description of sorption and filtration/remobilization processes are consistent with the interaction processes observed during the short-term CFM experiments at GTS mentioned above and can as well be applied for the simulation of the equilibrium sorption conditions by using appropriately high reaction rates (possibly relevant for long-term safety assessments of deep geological repositories in crystalline formations).
The developed model is realized in the transport code COFRAME (COlloid facilitated transport in FRActured-porous MEdia). COFRAME is a computation module of the RepoTREND program package (Reiche 2012) that is used for the assessment of the longterm safety of a repository for radioactive waste. An important issue of this work was to investigate the effect of kinetics on contaminant sorption on colloids and on the fracture surface as well as on interactions of colloids with the fracture surface. For this reason, a number of simulations varying the appropriate model parameters were performed with COFRAME. This work provides a model and a tool, which could be used to analyse, better understand and explain the results obtained in the field experiments performed at GTS. In this paper, the results of the modelling and simulation of one of these experiments, namely CFM RUN 10-01, are presented.

Physical System
Most rock types, especially crystalline rock formations, consist of a system of fractures, sheared zones and other mechanical defects which are generally interconnected over large distances. In this paper, all these defects are summarized by the term fracture. Fractures can be open as well as partially or totally filled with loose or solid substances. But all fractures have one property in common: Along a fracture, a groundwater flow is possible (fractures which are completely isolated from the conducting fracture network are not considered in this work). Hence, fractures act as main path for the advective transport of species dissolved in water.
Compared to this, the groundwater flow through the undisturbed rock matrix (solid rock between fractures) of a crystalline rock formation can be considered as negligible. But the water volume in the fully saturated rock matrix is often much larger than the water volume in fractures. Over a long timescale, due to molecular diffusion and sorption processes, significant amounts of solutes can be retarded and accumulated by the rock matrix.
There are several model concepts which adequately describe the physical system as well as the relevant physical and chemical processes. The model of the fractured-porous medium developed in this work is based on the so-called double-porosity model (Barenblatt et al. 1960). According to this approach, a fractured-porous medium is interpreted as the combination of two separate but interacting systems with different porosities: the porous rock matrix and the fracture network.
According to our approach, the geosphere is modelled as a single planar fracture with the aperture 2b (m) and the length L (m), thereby 2b L (see Fig. 1). The fracture and the pores of the rock matrix are completely saturated with water. Because of the symmetry, only one half of the fracture and the adjoining half of the rock matrix must be considered. It is assumed that the groundwater flow only occurs within the fracture. Radionuclides or other contaminants may diffuse through the fracture surface into stagnant matrix pore water, but colloids cannot do so due to their larger size. Further simplifications and additional model assumptions are: 1. Groundwater flow is assumed to be laminar. Because of this assumption, Darcy's law can be applied. 2. The dispersion of radionuclides and colloids in the fracture in flow direction can be described by Fick's law. 3. Transversal dispersion is neglected, and the complete mixing of contaminants and colloids across the fracture occurs immediately. This assumption allows a 1D modelling of the contaminant and colloid transport in the fracture. 4. The contaminant transport along the fracture is significantly faster than the transport in the rock matrix. Thus, the matrix diffusion can be restricted to a 1D problem in space, directed normal to the groundwater flow. 5. The matrix diffusion is limited by the maximum penetration depth x m (m). The introduction of this parameter helps to optimize the numerical model and to improve the program performance. This parameter should be set to an appropriate value and thus is not a model restriction.

Fig. 1
Schematic representation of colloid-facilitated transport in a fractured-porous medium (modified after Baek and Pitt 1996) The modelled system is illustrated in Fig. 1. The z-axis is the direction of the groundwater flow along the fracture, and the x-axis is directed towards the rock matrix. The line x = 0 is the symmetric axis along the fracture. The inflow boundary is assumed at z = 0. The formulation of the transport equations is based on mass balance over a representative elementary volume (Bear 1979;Helmig 1997). The following processes are taken into account (see Figs. 1, 2, exchanging processes between components are labelled Q): -advective transport, longitudinal mechanical dispersion and molecular diffusion of contaminants and colloids in the fracture, -contaminant exchange between the fracture and the adjoining rock matrix due to molecular diffusion Q fp , -limited molecular diffusion of contaminants within the rock matrix in the direction perpendicular to the fracture axis, -sorption of contaminants within the rock matrix Q pr , -sorption of contaminants onto the fracture surface Q fr , -sorption of contaminants on mobile and immobile colloids Q frm and Q fri , -filtration and remobilization of colloids Q c due to their interaction with the fracture surface as well as thereby induced filtration and remobilization of contaminants adsorbed on colloids Q cr , -radioactive decay under consideration of nuclide chains.
In the modelling system, contaminants and colloids can be present as the following components ( Fig. 2): -dissolved mobile contaminants in the fracture with concentration C fr (mol/m 3 ), -dissolved mobile contaminants in the pore water of the rock matrix C pr (mol/m 3 ), -mobile contaminants adsorbed on mobile colloids C cr (mol/m 3 ), Fig. 2 Exchange processes between the components that represent radionuclides and colloids in the modelling system -immobile contaminants adsorbed on the fracture surface S fr (mol/m 2 ), -immobile contaminants adsorbed on rock matrix pore surfaces S pr (mol/kg), -immobile contaminants sorbed on filtered colloids S cr (mol/m 2 ), -mobile colloids C c (mol/m 3 ), -immobile colloids filtered on the fracture surface S c (mol/m 2 ).
Additional assumptions about colloids are: -All colloids have the same size (monodisperse) and properties. -All colloids are injected into the system by the groundwater through the inflow boundary (z = 0). -Contaminants sorbed on the colloids do not impact on the colloid properties.

Fluid Flow Model
The approach developed here for an adequate mapping of the natural system to obtain the resulting physical model as a single planar fracture with adjacent rock matrix is based on a model proposed in NAGRA (1994). When creating a model of the real area, it must be taken into account which part of the fluid-accessible cross section actually allows advective transport. A real domain has the cross-sectional area A (m 2 ) (e.g. the cross-sectional area of the repository for radioactive waste) and the length L (m).
The cross-sectional area A is intersected with a number of small-scale water-conducting channels (traces). Each channel has a width g i (m) and an aperture 2b (average value). The total trace length of channel per unit area g (m/m 2 ) is defined as: This domain is modelled as a single planar fracture situated in porous rock (Fig. 3). The model fracture parameters are as follows: length L, aperture 2b and width g i . This ensures that the groundwater volume flow rate Q (m 3 /year) passing the fracture corresponds to the flow rate through the real transport domain. The cross-sectional area in the model should be equal A = 2x max g i with 2x max the height of the cross section. In this way, it is guaranteed that the real ratio of the part of the cross section accessible for advective flow to the part not accessible is kept in the model.

The Darcy velocity (filter velocity) q (m/year) is defined as:
The effective flow velocity u (m/year), i.e. the average transport velocity of the groundwater and accordingly the velocity of the dissolved contaminants, is coupled with the Darcy velocity q by the effective flow-through fraction of the volume of the relevant media. The effective flow-through volume of the modelled transport domain is 2b · g A · L · n f (groundwater flows only along the fracture), where n f (−) is the porosity of the fracture filling (n f = 1 in case of an unfilled fracture). Accordingly, the fraction of the flow-through volume in the total volume A · L is 2bgn f . This factor may be regarded as a kind of the porosity of water-conducting channels in the considered transport domain. Hence, the effective flow velocity u in the fracture can be calculated as: One more comment on the approach presented here: In general, the hydraulic fracture aperture 2b used in the model differs from the average mechanical aperture of the channels. Many different aspects affect the transportation ability of the fracture system. All those aspects must be identified and taken into account in determining the fracture aperture 2b as well as in determining other parameters. More details on these issue are contained in Reiche et al. (2014).

Modelling of Sorption Processes
The following approaches for describing the sorption processes are used in the model of the transport system presented in this work.

Sorption of Contaminants Dissolved in the Pore Water on the Pore Surfaces
The process of adsorption of contaminants dissolved in the pore water on the pore surfaces of the rock matrix can be assumed to be in equilibrium. This assumption is reasonable because the contaminant flow in the rock matrix is very slow so that the adsorption can be considered as a fast process compared to the transport processes. For the description of the reversible equilibrium adsorption reaction, the linear isotherm is used: where K pr (m 3 /kg) is the distribution coefficient for contaminants in the rock matrix.

Sorption of Dissolved Contaminants on the Fracture Surface
Two alternative concepts for the adsorption of contaminants dissolved in the fracture water on the fracture surface are realized: -linear equilibrium reaction: -linear kinetic reaction: where K fr (m) is the distribution coefficient for contaminants on the fracture surface and k fr (year −1 ) is the reaction rate for the kinetic reaction of contaminants adsorbed on the fracture surface.

Sorption of Dissolved Contaminants on Colloids
The sorption processes for contaminants dissolved in fracture water with colloids (both mobile and filtered) are treated as reversible linear kinetic reactions (Ibaraki and Sudicky 1995): where K frm (m 3 /mol) and K fri (m 3 /mol) are the distribution coefficients for contaminants on the mobile and immobile colloids, k frm (year −1 ) and k fri (year −1 ) are the reaction rates for the kinetic reaction of contaminants adsorbed on the mobile and immobile colloids, respectively. The kinetic approach is justified for processes that are slow compared to the transport processes, where an equilibrium between the dissolved and adsorbed phase has not been established. The use of the linear isotherm in the kinetic term means that colloids can adsorb an unlimited amount of contaminants. The application of this isotherm is reasonable for low contaminant concentrations.

Interaction of Colloids with Fracture Surface
The filtration of colloids onto the fracture surface is modelled as reversible, i.e. the colloids filtered on the fracture surface may be released back into the solution. This filtration/remobilization process is described by the following linear kinetic reaction (Oswald and Ibaraki 2001;Li et al. 2004): where λ f (m −1 ) and R mb (year −1 ) are the filtration and the remobilization coefficients, respectively. u c (m/year) is the average colloid velocity. The filtration and remobilization coefficients indicate the ability of the solid material to attach colloids and in general depend on colloid size and properties. Colloids considered in this work are assumed having the same size and properties so that the filtration and remobilization coefficients can be represented by constant values. The exchange term Q cr expressing the filtration and remobilization of colloid-bound contaminants due to the interaction of colloids with the fracture surface is consistent with the Eq. (9) and given as: 3 Mathematical Model

General Transport Equation
In the modelled transport system, mobile and immobile contaminants and colloids are represented by eight different components P (Fig. 2). It is assumed that the dispersion mechanism Table 1 Dispersion coefficients D P and average velocities u P for each component P It is also assumed that such parameters as porosity, velocity, dispersion and diffusion coefficients are constant within the modelling area during the total simulation time. Based on these assumptions, the general 1D-governing equations for component P can be formulated as follows: where t (year) is the time, y (m) is the coordinate in space, n (−) is the porosity, u P (m/year) and D P (m 2 /year) are the average velocity and the dispersion coefficient of the component P, respectively. Q P (mol/m 3 /year) or (mol/m 2 /year −1 ) is, respectively, the volume-and surface-related source-sink term. P (mol/m 3 ) and (mol/m 2 ) is the volume-and surfacerelated concentration of the appropriate component, respectively.
The dispersion coefficient D P is typically expressed by the combination of the molecular diffusion and mechanical dispersion: where D m P (m 2 /year) is the molecular diffusion coefficient and α P (m) is the longitudinal dispersion length of the component P. Table 1 provides an overview of the velocities and dispersion coefficients for each component P. D f (m 2 /year) is the dispersion coefficient for contaminants dissolved in the fracture water, and D c (m 2 /year) is the dispersion coefficient for mobile colloids and for contaminants adsorbed on mobile colloids. The molecular diffusion is the only transport mechanism in the rock matrix. Therefore, the dispersion coefficient for contaminants in the rock matrix is reduced to the molecular diffusion coefficient D p (m 2 /year).
In the following, the transport equation for each component P is established, based on the general transport Eq. (11) and considering the above assumptions.

Diffusive Contaminant Flux Across the Fracture-Matrix Interface
The diffusive contaminant flux Q fp (mol/m 2 /year) normal to the fracture-matrix interface can be expressed by Fick's law as: where n p (−) is the matrix porosity.

Colloid Transport Along the Fracture
The colloid transport along the fracture can be described as: The hydrodynamic dispersion coefficient for colloids is defined as: where D mc (m 2 /year) and α c (m) are the molecular diffusion coefficient and the longitudinal dispersion length of mobile colloids, respectively.

Contaminant Transport in the Rock Matrix
Due to the equilibrium sorption concept (4) in the rock matrix, the concentration of the adsorbed contaminants S pr can be expressed in terms of the concentration of the dissolved contaminants C pr . This allows the description of the contaminant transport in the rock matrix with only one differential equation as follows, taking into account the decay of radionuclides under consideration of nuclide chains: where n p (−) is the matrix porosity, λ (year −1 ) is the decay constant for radionuclides, k is the index for parent nuclides of the relevant radionuclide. R p (−) is the retardation factor defined in the conventional way as: where ρ (kg/m 3 ) is the density of the rock matrix.

Colloid-Facilitated Radionuclide Transport Along the Fracture
The equation which describes the transport along the fracture of contaminants dissolved in the fracture water and exchange between the relevant components can be formulated as: The hydrodynamic dispersion coefficient D f for dissolved contaminants is defined as: where D mf (m 2 /year) and α f (m) are the molecular diffusion coefficient and the longitudinal dispersion length of contaminants dissolved in the fracture water, respectively. The appropriate equation for contaminants adsorbed on mobile colloids can be written as: The differential equation describing the amount of contaminants adsorbed on filtered colloids is: In the case of the linear kinetic approach for the sorption of contaminants on the fracture surface (6), the component S fr is described by a separate equation: If the equilibrium sorption reaction (5) for contaminants with the fracture surface can be assumed, then by introducing the retardation factor R fr (−) the transport Eq. (19) can be reformulated to:

System of Equations
Equations (14) and (15) for mobile C c and immobile S c colloids are independent of contaminant concentrations. Therefore, these equations can be solved first. The results can then be used as input parameters for solving the system of -five coupled differential Eqs. (17), (19), (21), (22) and (23) for unknown concentrations C pr , C fr , C cr , S cr and S fr in the case of the linear kinetic approach for the sorption of contaminants on the fracture surface or -four coupled differential Eqs. (17), (25), (21), and (22) for unknown concentrations C pr , C fr , C cr and S cr in the case of the linear equilibrium sorption reaction for contaminants with the fracture surface.

Initial Conditions
The groundwater can contain natural background concentrations of mobile contaminants C r0 (mol/m 3 ) and mobile colloids C c0 (mol/m 3 ). The modelling system is initially in equilibrium. This assumption requires in particular for colloids that remobilization of colloids occurs R mb > 0 in case colloids are filtered on the fracture surface λ f > 0. Furthermore, in equilibrium state, the concentration of contaminants in the pore water of the rock matrix is assumed to be equal to the concentration of contaminants dissolved in the fracture water. With those considerations in mind, the initial values of all components can be derived from two known (easy measurable) background concentrations mentioned above as follows: In the case of a kinetic approach for sorption of contaminants on the fracture surface, the initial condition for a concentration of contaminants adsorbed on the surface fracture is:

Boundary Conditions
For the inflow boundary, the following assumptions for mobile components have been made: -Contaminants and colloids (both natural and from an artificial source) are transported into the system by the groundwater through the inflow boundary (no sources of contaminants and colloids within the system). -The input functions are known over the total simulation time.
In this case, the Dirichlet (or first-type) boundary condition can be applied for the inflow boundary. For colloids, it specifies the concentration values on the boundary and can be given as follows: where C c (t) (mol/m 3 ) is the concentration of colloids from a source at the inflow boundary z = 0 at the time t, where C r (t) (mol/m 3 ) is the concentration of mobile contaminants from a source at the inflow boundary z = 0 at the time t, and p (−) is the element-specific multiplication factor that defines which fraction of contaminants coming from a source at the inflow boundary at the time t is dissolved in the groundwater. The parameter p(t) can be calculated automatically by the program based on the equilibrium assumption: or specified by the user as element-specific constants p = const. At the outflow boundary, the transmission boundary condition is applied: Furthermore, considering that contaminants can pass the rock matrix only to the depth x m , the following boundary conditions for the component C pr can be applied:

Numerical Solution
Generally, it is impossible to solve the equations for colloids and contaminants described in the previous section analytically. In this study, the transport equations coupled with the initial and boundary conditions are numerically solved using the Crank-Nicolson implicit finite difference method (Wikipedia; Crank and Nicolson 1947). This method, which is based on the trapezoid rule in time, is unconditionally stable and second-order accurate in time. Furthermore, the numerical method used employs an upwind difference approximation for the first derivative in the spatial coordinate and a central difference approximation for the second derivative in the spatial coordinate. The discretization of the colloid transport Eq. (14) leads to a system of algebraic equations with a tridiagonal coefficient matrix, which can be directly solved by the Thomas algorithm (Wikipedia; Thomas 1949). Coupled Eqs. (17), (19), (21) or (17), (25), (21) describing the contaminant transport lead to an unsymmetric sparse linear system, which in this study is solved by the unsymmetric multifrontal (Davie 2004) method.
The absolute stability of this numerical method allows a high flexibility regarding the choice of the time step size. The time step size is determined only by the requirements on the accuracy of results. The concept of using variable time stepping was developed to achieve the best performance by maintaining the required results accuracy. The realized algorithm uses an initial estimated time step and reduces it until the required resolution in time is obtained whenever needed or increases it so that the admissible accuracy is not exceeded. The decision of when to reduce or to increase is made adaptively using the numerical solution that was calculated in the previous time step. The maximum relative difference is determined in the following way: where P m i and P m−1 i are concentrations at node i calculated at current m and previous m − 1 time step, respectively.
The user defines the tolerance limits for changes in the concentration d min and d max . The time step size will be reduced according to: and appropriately increased: However, the time step size can be reduced or increased only until the permissible lower t min or appropriate upper t max limit. These limits change depending on the current time t in the following way: where p max and p min are multiplication constants. The similar approaches for the time step size control are described in Ibaraki and Sudicky (1995) and Forsyth and Sammon (1986).

Model Verification
Verification of the numerical model described above and realized in the transport code COFRAME is performed using the following existing analytical solutions: -Analytical solution of the transport equation for mobile colloids C c along the fracture that is initially colloid-free C c (z, 0) = 0 (z > 0), by a permanent constant colloid inflow C c (0, t) = C c,0 = const (t ≥ 0) and the condition at the outflow boundary ∂C c (∞,t) ∂z = 0 (t ≥ 0) for the special case without remobilization R mb = 0 and without colloid diffusion into the rock matrix given by Abdel-Salam and Chrysikopoulos (1994), -Analytical solution of the transport equations for radionuclides along the fracture C fr and in the porous rock matrix C pr for a system without colloids C c (t) = 0 obtained by Tang et al. (1981) in the case of the permanent constant radionuclide inflow C fr (0, t) = C fr,0 = const (t ≥ 0) into an initial contaminant-free model domain C fr (z, 0) = 0, C pr (x, z, 0) = 0 (x > 0, z > 0) with the following boundary conditions on the outflow and in the rock matrix: , -Analytical solution given in Tien and Li (2002), Baek and Pitt (1996), Ahn (1988) for radionuclides dissolved in the fracture water C fr in the presence of the constant colloid concentration C c (t) = C c,0 = const (t ≥ 0) whereby the colloid filtration and remobilization as well as all sorption processes for radionuclides are described by a linear equilibrium reaction. The inflow boundary condition for C fr is C fr (z = 0, t) = C fr,0 exp (−λt), and the further initial and boundary conditions are defined as in the case above.
In all cases, the numerical solutions performed by COFRAME showed an excellent match with the analytical solutions (Reiche et al. 2014).

Analysing the Effect of Kinetic Interactions
An important aspect of this work is to improve the understanding of kinetically controlled interaction processes between contaminants, colloids and solids during transport. In order to investigate the effect of kinetics as well as the role of the relevant parameters, a number of simulations varying the appropriate model parameters were performed with COFRAME.

Effect of Colloid Remobilization
Changing of environment parameters (such as flow rate, pH value or ionic strength) can lead to the remobilization of colloids filtrated on the fracture surface. Some numerical calculations have been performed with COFRAME to study how the effect of colloid remobilization impacts the system behaviour. First, a number of test cases were calculated in which only the remobilization coefficient R mb was verified in the range of 0 to 1 year −1 . The value of the filtration coefficient  λ f = 0.1 m −1 has been chosen in such a way that, if no remobilization occurs, then the majority of mobile colloids will be filtered on the fracture surface before they reach the outflow boundary of the model. Colloids are introduced into the system at the time t = 0 through the inflow boundary with a constant concentration C c,0 . The molecular diffusion of colloids in the fracture is neglected (D mc = 0), in order to focus on the interaction of colloids with the fracture surface and to avoid the influence of other effects where possible. All parameters used are listed in Table 2. Figure 4 shows the resulting breakthrough curves of colloids at the distance z = 50 m. The normalized concentration of mobile colloids C c /C c,0 in the case without remobilization (R mb = 0) is about 0.01. The higher the remobilization coefficient, the more filtered colloids are released, and consequently, the concentration of mobile colloids in the fracture increases.
Next, a set of simulations with different reaction rates was performed to investigate the effect of kinetics by the filtration/remobilization processes. The reaction rates have been varied in such a way that the ratio of filtration-to-remobilization reaction rates remains constant. The varied parameters are listed in Table 3. Further parameters used in the simulations are listed in Table 2. Figure 5 presents the spatial distributions of normalized concentration of mobile colloids C c /C c,0 at different times for the cases a to d. A fraction of colloids is filtered on the fracture surface. The higher the reaction rate, the faster mobile colloids change into the immobile phase (until equilibrium is reached), which further causes stronger retardation of the front.

Effect of Kinetic Adsorption for Radionuclides on Colloids and the Fracture Surface
If sorption processes of contaminants on colloids and on the fracture surface are slow compared to the timescale of other transport processes, then the amounts of dissolved and adsorbed contaminants are not in equilibrium. An equilibrium state is reached when the rates of sorption and desorption reactions are fast enough and exactly equal. In the transport model presented in this study, the kinetic sorption processes are described by linear kinetic sorption reactions (7), (8), (9). To analyse the effect of kinetics, a set of simulations with varying reaction rates was performed.
The values of most system parameters specified for the test system are within the range of the values which are characteristic for the dipole tests performed at GTS within the CFM project (Sect. 5.3). The molecular diffusion in the fracture as well as into the porous rock matrix and the radioactive decay are neglected (D mc = D mf = 0, D p ≈ 0, λ = 0) to focus on the effect of kinetics. The sorption parameters and filtration/remobilization coefficients are chosen in a way that considerable amounts of mobile as well as filtered colloids are available and none of the sorption processes in the fracture dominates.
The simulated system contains only natural background colloids with a concentration C c,0 = 0.1 kg/m 3 (this value is within the range that is typical of natural colloid concentrations in groundwater systems according to field experiments, e.g. Nuclear Science and Technology 1990;Lieser et al. 1990). For the parameters applied here, a significant concentration of filtered colloids S c,0 = 0.1 kg/m 2 is accumulated in the system. No additional colloids are inserted into the system at t ≥ 0, i.e. C c (t) = const and S c (t) = const. Colloids have the same velocity as contaminants u c = u = 1 m/year. The inflow of contaminants into the model domain begins at the time t = 0 and remains constant over the total simulation time (C r (t, z = 0) = const). The distribution of contaminants between colloids C cr and solution C fr at the inlet to the fracture z = 0 satisfies the equilibrium conditions and is calculated according to (28) and (30), respectively. The dispersion length of colloids and contaminants is α c = α f = 0.1 m. The fracture aperture is 2b = 0.001 m. The reaction rates have been varied in such a way that the ratios of sorption-to-desorption and filtration-to-remobilization reaction rates remain constant.
Colloids are present in a system in mobile and immobile phase. Mobile colloids enhance and filtered colloids retard the transport of contaminants. The impact of these competing processes depends on system parameters which have to be taken into account in the analysis of the results.
First, the effect of kinetics for sorption of contaminants on colloids is investigated on the basis of cases a and b (Table 4). Here, sorption of contaminants on the fracture surface is neglected (K fr = 0). Further parameters are listed in Table 4. Figure 6 shows the spatial distributions of the normalized concentration for mobile radionuclides (C fr + C cr ) /C r(z=0) at different times for cases a and b. For comparison, a case without colloids is plotted. In cases a and b, a fraction of radionuclides is adsorbed on filtered colloids, and the concentration of mobile radionuclides in the fracture increases accordingly slower compared to the case without colloids. For the system parameters used in the simulations, it can be concluded that the higher the reaction rate, the faster mobile radionuclides change into the immobile phase, which further causes stronger retardation of the front.
In general, it can be deduced that colloids cause a retardation effect in systems in which contaminants can be adsorbed on colloids but that no sorption of contaminants occurs on the fracture surface.
Further, the same system is considered, but in cases c, d and e, the sorption of contaminants on the fracture surface is modelled much faster than the sorption of contaminants on colloids so that an equilibrium for sorption on the fracture surface described by Eq. (5) can be assumed. The distribution coefficients K frm , K frm and K fr are chosen in such a way that none Table 4 Parameters for cases a-e  Table 4. In Fig. 7, the spatial distribution of the normalized concentrations for mobile radionuclides (C fr + C cr ) /C r(z=0) at different times for cases c, d and e is presented. For comparison, the case without colloids is also plotted. Here, mobile colloids act as a vehicle which accelerates the migration of radionuclides because contaminants dissolved in the fracture water adsorb on mobile colloids. Therefore, contaminants are to a lesser extent immobilized by sorption on the fracture surface or on filtered colloids (case c). With increasing rates of kinetic reactions (sorption/desorption of contaminants as well as filtration/remobilization of colloids), the increasing impact of filtered colloids reduces the facilitating effect of the mobile colloids (case d and e; see also cases a and b in Sect. 5.2.1). Because of the presence of colloids, the total concentration of mobile contaminants increases. The slower the rates of kinetic sorption reactions on colloids, the slower are contaminants sorbed on filtered colloids, leading to an increasing concentration of mobile contaminants. Figure 8 illustrates the total concentration of mobile contaminants N total and separately the concentrations of its components-concentration of contaminants dissolved in the fracture water N fr and adsorbed on mobile colloids N cr -at the time t = 50 year for cases c and e (all concentrations in Fig. 8 are normalized to C r (t, z = 0)). Concentrations of contaminants Integrated contaminant amount in the system for cases c and e as function of elution time for the following components: contaminants dissolved in the fracture water C fr , contaminants adsorbed on mobile C cr and filtered S cr colloids dissolved in the fracture water N fr in both cases are similar at that time. In case e, the amount of filtered colloids is significantly higher than in case c. Due to the higher reaction rates of colloids with the fracture surface (filtration/remobilization) as well as by sorption reactions, in case e, a considerably higher fraction of contaminants than in case c is sorbed on filtered colloids. This is illustrated even better in Fig. 9 which shows the amounts of radionuclide components C fr , C cr and S cr in cases c and e accumulated in the system over time. Next, the same system is considered, but kinetically controlled sorption processes on the fracture surface are assumed according to (6). The distribution coefficient for contaminants   Table 5. In cases f to h as well as in cases i to k, the rates are varied in the same way as in cases c to e, with the reaction rate for the kinetic reaction of contaminants adsorbed on the fracture surface k fr remaining constant in cases f to h and being varied in cases i to k. Results presented in Figs. 10 and 11 show that the sorption kinetics of the contaminants on the fracture surface mainly impact the concentration of the dissolved contaminants C fr and therefore the total mobile concentration (in the given system especially in the first time, if the equilibrium has not yet been reached). It has nearly no influence on the velocity of the contaminant front.

Application to a Dipole Tracer Test of the CFM Project at the GTS
The intention for developing the presented model was to provide a tool describing the colloid-facilitated transport of contaminants for the long-term safety assessment of deep geological repositories in crystalline formations. This model should cover the relevant interaction processes between contaminants, colloids and solids. At the Grimsel Test Site, dipole tests have been performed addressing colloid-facilitated transport in a fracture zone under near-natural flow conditions as expected in a granite host rock for a deep geological repository.
Primarily, the developed model should be able to depict the key processes observed during these experiments. A first step of the model validation was accomplished by application to one of these dipole tests. The results are presented below. Justifications for all assumptions used to model the experiment (e.g. reduction in the 3D-problem to a 1D-problem or the parameter derivation from independent experiments) are not treated in this paper. The emphasis was placed on the model itself and not on the detailed description of the experimental set-up. More information about the dipole tests performed at GTS can be found in Blechschmidt and Lanyon (2008), Huber et al. (2012), and new results are planned to be published soon.

Experiment Overview
In order to assess the impact of bentonite colloids on the radionuclide transport in granite, a series of laboratory and field experiments at the underground rock laboratory at GTS, Switzerland (www.grimsel.com), is performed within the international partner project CFM under the auspices of NAGRA, Switzerland (Blechschmidt and Lanyon 2008). One set of experiments are dipole tests in the shear zone with FEBEX bentonite (Ca-bentonite) colloids and radionuclides or homologues under defined flow conditions. In these experiments, a defined chemical solution with colloid-bound radionuclides or homologues is injected into a dipole field in the shear zone and breakthrough curves are obtained at the outflow. The fluid velocity and therewith the transport time is controlled by the outflow rate . From these experiments, information and data for the interaction of colloids with the fracture and for the kinetics of radionuclide colloid interaction can be derived.
In this study, for the model validation, one dipole experiment, namely CFM RUN 10-01, is selected . The dipole experiment is simulated with the transport code COFRAME. For the simulation, it is assumed that all colloids are of similar size having similar properties. This assumption is supported to some extent by the colloid size distribution before injection into the dipole and at different points in time during extraction from the dipole, which was determined by LIBD s-curve analysis . The size distributions are the same, indicating no size chromatography effect. Injections of colloids and tracers at the inflow boundary of the model area are applied. CFM RUN 10-01 was performed as point dilution test by circulating the tracer solution through an injection interval in a closed flow loop and allowing the natural flow in the shear zone to displace the tracer solution from the interval (Geckeis et al. 2011). This resulted in an injection function, where the tracer concentration exponentially decreased. The injection function was detected for the uranine tracer and applied in the model to the colloids and homologues regarding their total injected mass . The concentrations in the injection solution were as follows: Colloids = 99 mg/l, Eu = 7.87 × 10 −8 mol/l, Hf = 6.39 × 10 −8 mol/l, Tb = 7.16 × 10 −8 mol/l, Th = 6.41 × 10 −8 mol/l. The parameters describing the characteristics of the model are summarized in Table 6. For the experiment, the conservative tracer uranine, the stable homologues Tb-159, Hf-178 and natural Eu (mixture of 47.8 % Eu-151 and 52.2 % Eu-153) and the radionuclide Th-232 (T 1/2 = 1.405×10 10 year) are considered. Sorption coefficients of radionuclides on colloids and on fracture material are derived from batch experiments, see e.g. in Huber et al. (2012). Desorption rates from colloids are derived by fitting the experimental breakthrough curves. For contaminant sorption on the fracture surface, a high kinetic rate is assumed. Sorption parameters are summarized in Table 7.

Breakthrough Curves
The breakthrough curves for the experiment plotted in all simulation curves are scaled by a factor of 0.8 to take into account the recovery of the conservative tracer of 80 %. Breakthrough curves for uranine are plotted in Fig. 12. Therewith, the contaminant dispersion length and the cross section of the modelling area were adjusted.
In Fig. 13, the breakthrough curves for the colloids are shown. In order to fit the shape of the experimental curve, a reversible, kinetically controlled colloid filtration and remobilization    Fig. 14 Breakthrough curves for homologues from experiment CFM RUN 10-01 and modelling. The curves are scaled by an additional factor of 0.8 in order to account for an irreversible colloid filtration with a filtration coefficient of 0.04 m −1 , since they represent the colloid-bound fraction of the homologues were assumed. Further, the reduction in the recovery is caused by an additional irreversible colloid filtration. The respective data are summarized in Table 8 In a final step, the breakthrough curves for the homologues are simulated. They are shown in Fig. 14. Interaction coefficients of homologues with colloids and fracture are taken from batch experiments. The reaction rates k frm and k fri of the homologues desorbing from the colloids are adjusted to fit the experimental curves. The resulting simulation curves coincide well with the experimental curves.
In conclusion, the field experiment results can be explained by the following processes. As described above, the homologues are equilibrated with the colloid-bearing solution in advance and are quantitatively bound (∼99 % adsorbed) to the colloids, when injected into the dipole. During the transport through the fracture, the homologues desorb from the colloids and subsequently adsorb to the fracture-filling material. The velocity of these processes is determined by the detachment rates k frm and k fri of the homologues from the colloids, since the adsorption reaction of the homologues to the fracture-filling material is assumed to be fast.
The results show that the experiment can be successfully simulated with the new COFRAME code taking into account the experimental boundary conditions and using plausible assumptions for the interaction data of contaminants, colloids and fractured medium. Further applications to other CFM field tests as well as for colloid-facilitated radionuclide transport under repository conditions are planned and will be published in the near future. and interactions of colloids with the fracture surface are treated as linear kinetically controlled processes. For the sorption of contaminants on the fracture surface, two approaches are used: linear equilibrium and linear kinetic reactions. The numerical solution based on the Crank-Nicolson implicit finite difference method was realized in the COFRAME transport code and verified using existing analytical solutions for colloid and contaminant migration in a fracture-porous medium.
A useful approach was suggested in this work for adequate mapping of a natural system (fractured crystalline rock) onto its model consisting of a single planar fracture with adjacent rock matrix.
To verify the code and to investigate the effect of kinetics, a number of simulations were performed with COFRAME. Colloids can be present in a system as mobile and immobile components. In systems where contaminants do not sorb on the fracture surface, immobile colloids act as an additional retarding factor for contaminants. In systems where contaminants do sorb on the fracture surface, on the contrary, the presence of colloids can significantly accelerate the transport of contaminants and also lead to an increase in the total concentration of mobile contaminants compared to the case without colloids.
The reaction rates of the interacting processes for contaminants as well as for colloids affect the transport characteristics, which is best indicated by the shape of the spatial distribution curve (especial its tailing). The reaction rate for the sorption of contaminants on the fracture surface mainly impacts the concentration of mobile contaminants. The reaction rates for sorption of contaminants on colloids and for interacting of colloids with the fracture surface affect both accelerating and retarding effects on the contaminant migration. However, the impact and the extension of each of these competing effects do not only depend on reaction rates but also on the system parameters (e.g. on the accumulated amount of filtered colloids or distribution coefficients). For the system parameters used in the simulations (with significant amounts of mobile as well as immobile colloids), it can be concluded that the higher the rates of kinetic reactions, the faster mobile radionuclides change into the immobile phase due to their sorption on filtered colloids, which further causes stronger retardation of the front.
The developed model and the numerical solution were further validated by application to the data from a dipole test performed in the frame of the CFM project at GTS in Switzerland. The experiment was simulated with COFRAME. The modelled breakthrough curves coincide well with the experimental curves. The computed results also confirm that the experimental results can be best explained using a kinetic approach for desorption of contaminants from colloids with subsequent fast sorption on the fracture-filling material and for interaction of colloids with the fracture.