Plan
Comptes Rendus

Hydrology, Environment
Multiscale modelling of transport in clays from the molecular to the sample scale
Comptes Rendus. Géoscience, Volume 346 (2014) no. 11-12, pp. 298-306.

Résumé

We report some recent applications of multiscale modelling to the transport of ions, water and CO2 in clays. On the one hand, simulations on different scales allow us to investigate the physicochemical processes underlying the geochemical and transport behaviour of these fluids in the interparticle pores and at the surface of clay minerals. We discuss more specifically the insights gained from molecular simulations into the acidity of surface edge sites, ion exchange and the behaviour of clay interlayers in contact with a CO2 reservoir. On the other hand, upscaling the descriptions from the molecular level to the macroscopic scale without forgetting the fundamental role of interfaces on the mesoscopic scale provides a means to capture complex phenomena such as electrokinetic couplings. We illustrate the complementarity of molecular dynamics, lattice-based mesoscopic simulations and Pore Network Models to address this issue.

Métadonnées
Reçu le :
Accepté le :
Publié le :
DOI : 10.1016/j.crte.2014.07.002
Mots clés : Clay minerals, Multiscale simulation, Acidity, Ion exchange, Sorption, Electrokinetic effects

Benjamin Rotenberg 1, 2 ; Virginie Marry 1, 2 ; Mathieu Salanne 1, 2 ; Marie Jardat 1, 2 ; Pierre Turq 1, 2

1 Sorbonne Universités, UPMC Univ. Paris 06, UMR 8234 PHENIX, 75005 Paris, France
2 CNRS, UMR 8234 PHENIX, 75005 Paris, France
@article{CRGEOS_2014__346_11-12_298_0,
     author = {Benjamin Rotenberg and Virginie Marry and Mathieu Salanne and Marie Jardat and Pierre Turq},
     title = {Multiscale modelling of transport in clays from the molecular to the sample scale},
     journal = {Comptes Rendus. G\'eoscience},
     pages = {298--306},
     publisher = {Elsevier},
     volume = {346},
     number = {11-12},
     year = {2014},
     doi = {10.1016/j.crte.2014.07.002},
     language = {en},
}
TY  - JOUR
AU  - Benjamin Rotenberg
AU  - Virginie Marry
AU  - Mathieu Salanne
AU  - Marie Jardat
AU  - Pierre Turq
TI  - Multiscale modelling of transport in clays from the molecular to the sample scale
JO  - Comptes Rendus. Géoscience
PY  - 2014
SP  - 298
EP  - 306
VL  - 346
IS  - 11-12
PB  - Elsevier
DO  - 10.1016/j.crte.2014.07.002
LA  - en
ID  - CRGEOS_2014__346_11-12_298_0
ER  - 
%0 Journal Article
%A Benjamin Rotenberg
%A Virginie Marry
%A Mathieu Salanne
%A Marie Jardat
%A Pierre Turq
%T Multiscale modelling of transport in clays from the molecular to the sample scale
%J Comptes Rendus. Géoscience
%D 2014
%P 298-306
%V 346
%N 11-12
%I Elsevier
%R 10.1016/j.crte.2014.07.002
%G en
%F CRGEOS_2014__346_11-12_298_0
Benjamin Rotenberg; Virginie Marry; Mathieu Salanne; Marie Jardat; Pierre Turq. Multiscale modelling of transport in clays from the molecular to the sample scale. Comptes Rendus. Géoscience, Volume 346 (2014) no. 11-12, pp. 298-306. doi : 10.1016/j.crte.2014.07.002. https://comptes-rendus.academie-sciences.fr/geoscience/articles/10.1016/j.crte.2014.07.002/

Version originale du texte intégral

1 Introduction

Understanding the transport of water, ions, gas or oil through clay barriers is essential in the contexts of hydrology, petroleum and gas reservoir engineering, as well as the geological disposal of radioactive waste (ANDRA, 2005; Bradbury and Baeyens, 2003) or CO2 sequestration (DePaolo and Cole, 2013; Gaus, 2010). Clay rocks display a complex multiscale structure, from the regular stacking of aluminosilicate lamellae on the microscopic scale to their disordered assembly into particles on mesoscopic scales (10−9 to 10−6 nm) with corresponding interlayer and interparticle porosities, up to the macroscopic assembly of clays with other minerals (e.g., carbonates or quartz). To this hierarchy of length scales correspond a variety of physicochemical phenomena occurring over times spanning several orders of magnitude, many of them being associated with the structural charge of the mineral surface: chemical reactivity of surface sites, electrostatic interactions between the surface and ions in the interstitial solution, specific ion effects, wetting, swelling, electrokinetic couplings, etc. The electromagnetic signature of materials, in response to applied electromagnetic or acoustic waves, is at the basis of many logging tools and of seismo-electric exploration (see, e.g., Pride, 1994).

From the modelling point of view, this challenge set by multiple length and time scales can only be addressed within the framework of a multiscale strategy. On the one hand, this means adopting various complementary descriptions of the same system to capture different types of complexities. The reorganization of electrons around nuclei during chemical reactions requires using quantum descriptions, but their computational cost limits their use to a few hundred atoms (Boek and Sprik, 2003; Liu et al., 2013; Suter et al., 2008). Molecular simulations allow us to accurately describe clay layers and interlayers as well as the interface between clay particles and solutions over a few nanometres (Boek et al., 1995; Ferrage et al., 2011; Greathouse and Cygan, 2006; Hensen and Smit, 2002; Marry et al., 2002; Michot et al., 2012; Rotenberg et al., 2010a; Sposito et al., 1999). Mesoscopic simulations are then needed to model phenomena in pores with sizes in the range 10–100 nm (Boek and Venturoli, 2010; Rotenberg et al., 2010b; Stukan et al., 2010). On the other hand, one should make the link explicit between the various levels of description, exploiting more accurate models to calibrate coarser ones retaining the relevant information on larger scales. This “bottom-up” strategy includes, for example, using ab initio calculations to parameterize force fields (Cygan et al., 2004), or molecular simulation to derive simple kinetic models (Carof et al., 2014; Dufrêche et al., 2010; Rotenberg et al., 2007a,b). Even though all the references cited above are applications to clay minerals, this strategy is of course very general and has been successful in many other contexts.

In the present work, we review some recent theoretical and numerical developments in our group, using molecular and mesoscale simulations to address a number of issues related to the transport of ions, water and other fluids in clays. In Section 2, we first illustrate the interest of various simulation techniques to gain insight into geochemistry problems: ab initio molecular dynamics (MD) to predict the acidity of clay edge sites, classical MD to understand the thermodynamics of ion exchange, and Grand-Canonical Monte Carlo (GCMC) to investigate the interaction of clay interlayers with CO2 reservoirs. In Section 3, we then show the interest of a multiscale simulation strategy to describe the transport of water and ions in clays, including electrokinetic couplings, from the molecular to the macroscopic sample scale, using a combination of molecular simulation, lattice-based mesoscopic simulation (Lattice Boltzmann Electrokinetics) and Pore Network Models.

2 Geochemistry: insights from simulation on the molecular scale

In order to predict the transport and retention of mobile species in complex, interfacial materials such as clays, one first needs to understand the underlying microscopic mechanisms. Here we examine three examples of such processes for which molecular simulation proved very useful in complementing experimental approaches.

2.1 Acidity of clay edges

Lateral surfaces of clay particles may provide sorption sites contributing significantly to the retention of heavy metal cations. This sorption depends largely on the protonation state of the various edge sites, which is controlled by their pKa and the pH of the solution (Tournassat et al., 2013). Fig. 1 illustrates silanol SiOH and aluminol AlOH and AlOH2 sites on the (0 1 0) surface of a dioctahedral TOT clay. Titration experiments only provide a global measure of the surface state and do not allow one to assess the state of individual sites (which however matters for cation sorption) (Bourg et al., 2007). In practice, an efficient strategy consists in introducing a priori or experimental structural information in semi-empirical models such as MUSIC (Hiemstra et al., 1996; Machesky et al., 2008) in order to reproduce the titration curve (Tournassat et al., 2004).

Fig. 1

Lateral surfaces of clay particles display edge sites arising from the reaction of broken SiO and AlO bonds with water molecules. The resulting groups, here silanol and aluminol on the (0 1 0) surface of a dioctahedral TOT clay such as montmorillonite, may protonate or deprotonate depending on the pH of the interstitial solution. The pKa of each site can be determined separately from ab initio simulations, whereas experimental titration data can only provide global information. Colour code: yellow = Si, green = Al, red = O, white = H. In the left part water molecules are in grey and the spheres represent Na+ (blue), Cs+ (orange) and Cl (pink) ions. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of the article.)

In this context, molecular simulation is an attractive alternative for evaluating the acidity of individual sites. The pKa is related to the reaction free energy for the transfer of the surface proton to a water molecule, e.g., for a silanol group:

SiOH+H2OaqSiO+H3Oaq+.(1)
This involves the breaking and formation of chemical bonds so that DFT-based simulations are necessary (even though reactive force fields are now available, they remain less accurate, in particular for proton transfer reactions). Since the relaxation of bond lengths on edge surfaces plays an important role on the acid–base properties (White and Zelazny, 1988), it has first been proposed to introduce refined structures from ab initio simulations into bond-valence-based models, both in terms of bond lengths (Bickmore et al., 2003) and of solvation structure (coordination numbers and bond lengths) (Bickmore et al., 2006; Machesky et al., 2008). Using ab initio molecular dynamics (AIMD), Churakov studied water confined between edges of pyrophyllite and evidenced transient proton exchange between surface groups mediated by surface water molecules (Churakov, 2007). He also proposed a first attempt to provide an acidity scale of the surface groups from their deprotonation energy in vacuum (Churakov, 2006), which neglects entropic contributions and the solvent reorganization after deprotonation.

By combining AIMD with thermodynamic integration (Sulpizi and Sprik, 2010), we computed the free energy for the deprotonation reaction (1) of SiOH, AlOH and AlOH2 sites on the (0 1 0) surface of pyrophyllite (Tazi et al., 2012). SiOH is slightly more acidic than AlOH2, with pKa values close to 7 for both sites, whereas AlOH does not deprotonate in water (pKa ∼ 22). In addition, the reorganization free energy, arising from bond length relaxation and solvent reorganization after deprotonation, is an important contribution to the overall free energy and cannot be neglected. For example, the strong stabilization of SiO after deprotonation is explained by the arrival of two water molecules donating H-bonds. Following the same strategy, Liu et al. recently extended this study to assess the role of AlIII→MgII substitutions near edges in montmorillonite (Liu et al., 2013). It was found that Mg substitutions increase the pKa of the neighbouring SiOH group by 2–3 pKa units. Overall, these studies suggest that the sites mainly responsible for heavy cation sorption should be the deprotonated SiOH and AlOH2 sites.

2.2 Ion exchange

One of the most important processes governing cation retention by clays is cation exchange, between ions in solution and interlayer counterions, as illustrated in Fig. 2. The exchange of interlayer Na+ counterions by Cs+ in solution is thermodynamically favourable (ΔrG < 0) and exothermic (ΔrH < 0). It was generally admitted that this was the result of favourable interactions (of various origins) between the Cs+ cation and the clay surface. This point of view was challenged by Teppen and Miller (2006), who computed by molecular simulation the free energy associated with the replacement of Na+ by larger alkaline cations, and found that this contribution was in fact unfavourable for the exchange.

Fig. 2

(Color online.) Cation exchange reaction between the clay interlayer and an aqueous solution. The contributions of both phases to the overall enthalpy and Gibbs free energy can be analyzed separately in molecular simulations.

The following thermodynamic cycle separates the contribution of ion replacement in both the clay interlayer (with a fixed number of water molecules per cation) and the aqueous phase (“aq” subscript):

(2)

This cycle involves the exchange reaction with ions in the gas phase (“g” subscript), where there are no interactions, used as a reference state. The associated thermodynamic quantities are the hydration enthalpy of both cations ΔhHNa and ΔhHCs, the standard ionic exchange reaction with the gas phase ΔrHg0 and that with the aqueous phase ΔrH0. Using molecular simulations, we showed that replacing Na+ by Cs+ in the interlayer is in fact an endothermic process (ΔrHg0>0), and that the overall exchange is exothermic (ΔrH0 < 0), only because of the hydration enthalpy difference ΔΔhH = ΔhHCs − ΔhHNa ∼ 133 kJ mol−1 associated with the replacement of Cs+ by Na+ in the aqueous phase. This conclusion, as well as the evolution of ΔrH0 with water content, was confirmed by microcalorimetry experiments (Rotenberg et al., 2009). Similarly, for the exchange of Cs+ traces in Na-montmorillonite we found that the contribution of the clay phase to the exchange is positive and large (ΔrG ∼ +120 kJ mol−1) and that the overall exchange is thermodynamically favourable only due to the replacement of Cs+ by Na+ in the aqueous phase.

While the thermodynamics of ion exchange can be studied by considering only the initial and final bulk states, understanding the exchange process requires taking into account explicitly the interface between the interlayer and the aqueous solution via the lateral surface. Using the setup illustrated in the left part of Fig. 1, we found that for Na-montmorillonite in the bilayer state water and cations can exchange via the edges with no or small activation barriers, whereas anions are excluded from the interlayer (we estimated the associated free energy cost to ∼23 kJ mol−1 in that case) (Rotenberg et al., 2007c). For larger distances between the surfaces and larger salt concentrations, this anion exclusion becomes less dramatic and some salt may enter the nanopores. In order to assess systematically the effect of salt concentration and interlayer distance on the ratio between the “internal” and “external” salt concentrations (Donnan effect), we introduced a continuous solvent model, parameterized on molecular simulation (Jardat et al., 2009). The correlations between ions modulate the salt exclusion, and two regimes were found: at low reservoir concentrations, the Coulomb attraction between ions increases the amount of salt in the interlayer space compared to the Poisson–Boltzmann predictions, while at high concentrations excluded volume (due to the finite size of the ions) dominates.

2.3 Contact with a CO2 reservoir

Several options are currently considered for the geological sequestration of CO2. One of them is to inject supercritical CO2 in saline aquifers: This fluid would partly be dissolved and trapped physically (microbubbles) or chemically (carbonate precipitation) in such porous formations, but some would eventually rise by buoyancy as a separate phase until it encounters an impermeable caprock. It is then important to predict the evolution of the often clay-rich caprock due to the contact with the reservoir (Gaus, 2010). In this context too, molecular simulation gives access to the microscopic mechanisms at play (Hamm et al., 2013). Among the many possible scenarios, we have considered the possibility of clay interlayer swelling or dehydration. In the latter case, this would lead to the shrinkage of clay particles, resulting in fractures increasing the permeability of the caprock.

To that end, we performed Grand-Canonical Monte Carlo simulations of sodium montmorillonite clays in contact with H2O and CO2 reservoirs. The number of molecules in the interlayer fluctuate as a result of matter exchange with the reservoirs which set the chemical potentials of the fluid (see Fig. 3). For given thermodynamic conditions in the reservoirs we monitor the interlayer composition and pressure for each interlayer distance, from which we determine the swelling free energy. The thermodynamically stable states then correspond to the minima of this free energy. For the considered thermodynamic conditions (T = 348 K and P = 125 bar) we found that from the initial state with H2O only, the contact with the CO2 reservoir does not change the location of these minima (in other words, we should not expect swelling or shrinkage) (Botan et al., 2010). We further investigated the structure and dynamics of CO2 intercalated in montmorillonite interlayers, demonstrating in particular the slowing down of all interlayer species in the presence of CO2. Experiments have confirmed the presence of intercalated CO2 in hydrated montmorillonite interlayers (Giesting et al., 2012a,b) and recent molecular dynamics studies provided further insights into the interlayer properties (Cygan et al., 2012; Myshakin et al., 2013).

Fig. 3

(Color online.) Grand-Canonical Monte Carlo simulations. H2O and CO2 molecules are exchanged between the interlayer of a sodium montmorillonite clay and the reservoirs, which set the chemical potentials of both species. While in the reservoirs, under the considered pressure and temperature, a H2O-rich phase and a CO2-rich phase coexist, one finds only the intercalation of CO2 molecules inside the hydrated interlayers, which corresponds to the minima of the swelling free energy as a function of the basal spacing (right). See Botan et al. (2010) for more information.

3 Transport, from the molecular to the sample scale

Modelling the transport of fluids in charged multiscale materials such as clay is a challenge because it is limited by the smallest pores, with sizes in the 1–100-nm range, where interfacial effects such as wetting or electrokinetic couplings play a dominant role. Under confinement down to the nanometre scale, one should also expect the departure from continuous hydrodynamics (Bocquet and Charlaix, 2010). Once continuous transport equations on the pore scale are known (e.g., Stokes), they can be upscaled to the corresponding ones (e.g., Darcy's law) on the macroscopic sample scale, using numerical or analytical homogenization approaches. Here we illustrate how simulations on three different scales allow one to assess the relevance of continuous models, to solve the electrokinetic equations on the pore scale, and to upscale the electrokinetic couplings in complex heterogeneous materials.

3.1 Hydrodynamics in clay nanopores

Using MD simulations, we investigated the applicability of continuous hydrodynamics (Navier–Stokes equation) in clay nanopores, with pore widths ranging from 2 to 9 nm (Botan et al., 2011). Fig. 4 illustrates the simulated setup for sodium montmorillonite clay surfaces separated by 3.5 nm. The solvent density displays strong oscillations near the surface, reflecting the layering of water over 2–3 molecular diameters. Depending on their size and valency, cations may form either inner-sphere or outer-sphere complexes at the surface (Greathouse and Cygan, 2006; Marry et al., 2008). In the absence of substitutions in the mineral layers (hence of counterions), the basal surfaces are hydrophobic, even though the presence of hydroxyl groups may render the competition between water–water and water–surface interactions more delicate, as in the case of talc (Rotenberg et al., 2011).

Fig. 4

Non-equilibrium molecular dynamics. Under the effect of an applied force, mimicking a pressure gradient, the steady-state flow profile (black symbols) may deviate from the prediction of continuous hydrodynamics, in particular due to the layering of the solvent density (blue dashed line) at the surfaces. When the distance between sodium montmorillonite clay surfaces exceeds 3 nm, the velocity profile far from the surface is reasonably well described by a Poiseuille flow; however only if a slip boundary condition is introduced, with a slip length derived from molecular simulations. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of the article.)

Reprinted with permission from Botan et al. (2010). Copyright 2011 American Chemical Society.

Once the pore composition is known for given thermodynamic conditions (e.g., from prior Grand-Canonical Monte Carlo simulations, see Section 2.3), it is possible to simulate the effect of a pressure gradient by applying a force along the surface on the fluid. The steady-state velocity profile is then measured and compared to the prediction of continuous hydrodynamics (parabolic Poiseuille flow profile in this case). For sufficiently large pores, this profile is parabolic near the centre of the pore, and the effective viscosity corresponding to the measured curvature is comparable to the bulk viscosity. However, deviations are observed near the surface, as an expected result of layering. More importantly, the profile far from the surface is reproduced by the Navier–Stokes equation only if a slip boundary condition is used. This in turn allows the determination of the corresponding slip length from molecular simulation. The latter is very small (∼2.1Å), but its effects on the overall flow through the pore can be large for pore widths of a few nanometres.

More recently, we followed the same strategy to investigate the electro-osmotic flow induced by an electric field along the clay surfaces (Botan et al., 2013; Rotenberg and Pagonabarraga, 2013). Here again, in order to reproduce the molecular simulation results with continuous equations, it is essential to introduce appropriate slip boundary conditions at the clay/solution interface, as had been previously observed in equilibrium MD (Dufrêche et al., 2005). This study also confirmed the limitations of the Poisson–Boltzmann equation for the description of the ionic profiles and of the resulting electric force profile in the presence of the applied electric field.

3.2 Electrokinetic couplings on the pore scale

While molecular simulation provides the most accurate description of the clay/solution interface, its computational cost prevents its use for large systems that are necessary to capture the effect of the complex pore structure over 10–100s of nanometres. To that end, one relies on simpler descriptions. The standard model for electrokinetic couplings on the pore scale is thus the so-called Poisson–Nernst–Planck (PNP) model, which couples (i) the Nernst–Planck equation for the transport of ions under the combined effects of advection, diffusion and migration due to the local electric field, (ii) Poisson's equation to determine the electrostatic potential from the electric charge distribution and (iii) Navier–Stokes (NS) equations for the evolution of the fluid under local pressure and electric potential gradients, taking into account viscous dissipation.

These equations can be solved numerically using finite element or volume methods. For example, Adler and co-workers used this direct numerical resolution in various complex systems (random packings, reconstructed and fractured porous media) (Coelho et al., 1996; Marino et al., 2001), demonstrating in particular a universal electrokinetic behaviour if appropriate rescaled quantities are introduced (Gupta et al., 2006, 2008). Recently, alternative methods have been proposed to simulate electrokinetic effects starting from a more fundamental description of the fluid than the PNP and NS equation (Pagonabarraga et al., 2010). For example, Capuani et al. proposed a hybrid lattice based approach (Lattice Boltzmann Electrokinetics, LBE) to capture the coupling of hydrodynamic flow with ion transport and the simulation of electrokinetic effects in colloidal suspensions (Capuani et al., 2004; Pagonabarraga et al., 2005). We applied this method to porous media and extended it to charged liquid–liquid interfaces (Rotenberg et al., 2008, 2010b) as well as to capture the effect of salt concentration gradients (osmosis) (Obliger et al., 2013).

As an example of the application of LBE in the perspective of a numerical homogenization strategy, Fig. 5a reports the coefficients of the transfer matrix describing cation and anion fluxes under an applied potential gradient, in a cylindrical channel of fixed radius and surface charge density, as a function of salt concentration. These fluxes include both the direct effect of the electric field on the ions (Nernst–Einstein contribution, proportional to the ion concentrations and their mobility) and the electro-osmotic contribution (advection by the fluid flow). The former drives cations and anions in opposite directions, while the latter, in the present case of a negative surface charge, increases the cation flux and mitigates that of anions. Such numerical simulations also allow one to assess the validity of approximate analytical theories obtained by linearizing the Poisson–Boltzmann equation, which can be used in models on larger scales (see Section 3.3). Fig. 5b reports the relative error of this approximate solution compared to the non-linear one for a fixed channel radius, as a function of the surface charge density and of the salt concentration.

Fig. 5

Lattice Boltzmann Electrokinetics (LBE). (a) Coefficients of the transfer matrix describing the cation (•) and anion (▪) fluxes under an applied potential (ψ) gradient, in a cylindrical channel of radius 5 nm and surface charge density σ = −0.08 e nm−2, as a function of salt concentration ρsalt. The total electrical permeability K±ψ=K±eo+K±NE (red and blue solid lines for cations and anions, respectively) as well as the Nernst–Einstein K±NE (orange and green dashed lines for cations and anions, respectively) are indicated. LBE simulations (symbols) are compared to the analytical results obtained with the linearized Poisson–Boltzmann (PB) equation (lines). (b) Relative error on the electro-osmotic permeability for cations K+eo due to the linearization of the PB equation, compared to the non-linear one (as determined by LBE simulations), as a function of the surface charge density σ and of the salt concentration ρsalt. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of the article.)

Reprinted with permission from Obliger et al. (2013). Copyright 2013 by the American Physical Society.

3.3 Electrokinetics: from the pore scale to the sample scale

Once transport properties are known on the pore scale, it is sometimes possible to upscale rigorously from a mathematical point of view the PNP equations using the homogenization approach (Allaire et al., 2010), even introducing mechanical couplings (Moyne and Murad, 2006a,b) or extensions to capture non-ideality of the electrolyte (Allaire et al., 2014). While very general in principle, most studies of electrokinetic couplings consider simple geometries (slits or cylinders) with dimensions and surface charge densities estimated from the macroscopic properties of the real system. Direct numerical simulation on complex porous networks is typically limited, due to the high computational cost, to system sizes in the range of tens to hundreds of nanometres. In order to investigate electrokinetic couplings on larger scales, including the effect of the heterogeneity of the material, we have recently proposed a simplified description based on the Pore Network Model (PNM). Such a model, originally developed by Fatt (1956) to predict multiphase flow properties in porous media, describes the porosity as a network of pores connected by channels, as illustrated in Fig. 6a. It has been extensively used and extended by petrophysicists in various situations, such as capillarity and multiphase flow through porous media (Békri et al., 2005; Blunt, 2001; van Dijke and Sorbie, 2002), or mineral dissolution and precipitation in the context of CO2 sequestration (Algive et al., 2010).

Fig. 6

Pore Network Model. (a) Transport properties on the network scale are determined by the transport matrix of each channel (red) connecting pores (yellow). (b) Channel diameter distribution, here of the Weibull form, for several average (d¯) and standard deviation (δ) in nm. (c) Evolution of the electro-osmotic coupling coefficient K0V describing the solvent flow under an applied electric field, as a function of salt concentration in the reservoir c, for various channel size distributions (d¯,δ). The surface charge density is fixed to −0.05 e nm−2 and lines indicate the result for a network composed of identical channels (δ = 0). (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of the article.)

Panels (b) and (c) reprinted with permission from Obliger et al. (2014). Copyright 2014 by the American Physical Society.

In a nutshell, the PNM approach amounts to solving a set of conservation equations on the nodes of the network (in analogy with Kirchhoff's law for a network of resistors), on the basis of local fluxes through the channels connecting the nodes, under the effect of external macroscopic gradients. In the case of electrokinetic couplings, this requires introducing a transport matrix for each channel relating the solvent, cation and anion fluxes to the pressure, concentration and potential gradients across the channel. In addition, the matrix for each channel depends not only on its diameter and surface charge density, but also on the salt concentration in the pores at both ends, via the Donnan equilibrium (see Section 2.2), so that the problem to be solved numerically is a non-linear one (Obliger et al., 2014).

From the knowledge of the transfer matrix on the channel scale, either analytically using approximations or numerically, e.g., from LBE simulations as described in Section 3.2, and of the channel diameter distribution, one can finally determine the macroscopic transfer matrix on the network scale. Examples of such distributions are shown on Fig. 6b, together with the effect of the sample heterogeneity on the macroscopic electro-osmotic coupling coefficient. In general, the coefficients of this matrix qualitatively behave as their microscopic counterpart for a channel with the average diameter. However, the combined effects of electrokinetic couplings on the local scale and of heterogeneity result in a decrease of the overall transport coefficients. This decrease is more pronounced for large surface charge densities and low salt concentrations. In addition, for a given average diameter, more heterogeneous samples result in stronger effects, due to the presence of smaller pores. With this new numerical tool at hand, we are currently examining the case of the Callovo-Oxfordian argilite in the context of the geological disposal of radioactive waste in France.

4 Conclusion

We have illustrated here some recent applications of multiscale modelling to the transport of ions, water and CO2 in clays. On the one hand, simulations on different scales allow us to address the physicochemical processes underlying the geochemical and transport behaviour of these fluids in the interparticle pores and at the surface of clay minerals. On the other hand, upscaling the descriptions from the molecular level to the macroscopic scale without forgetting the fundamental role of interfaces on the mesoscopic scale provides a means to capture complex phenomena such as electrokinetic couplings. This multiscale approach, followed in our group mainly in the case of transport, is also very fruitful for the description of other properties of clays, such as mechanical ones (Brisard and Dormieux, 2010; Carrier et al., 2014). Extensions to more complex fluids would also provide an efficient route for the modelling of hydrocarbons in oil and gas shales.

In the future, one should also benefit from recent numerical (Tyagi et al., 2013) and experimental (Brisard et al., 2012; Levitz, 2007) developments for the generation of realistic numerical samples for the description of real materials. As a recent example, Robinet et al. recently simulated the diffusion of solutes in 3D-images of a Callovo-Oxfordian clay-rich rock obtained by SEM and micro-CT experiments to investigate the effect of mineral distribution (Robinet et al., 2012). Multiscale experiments using NMR also provide an ideal tool to investigate the multiscale dynamics of mobile species in such complex materials (Porion et al., 2013).

Acknowledgements

The authors are grateful to many colleagues and students who have contributed to this work over the last ten years, in particular Alexandru Boţan, Sami Tazi, Amaël Obliger, Magali Duvail, Antoine Carof, Natalie Malikova, Jean-François Dufrêche, Rodolphe Vuilleumier, Jean-Pierre Hansen, Daniel Coelho, Benoît Noetinger, and Samir Bekri.


Bibliographie

[Algive et al., 2010] L. Algive; S. Békri; O. Vizika Pore-network modeling dedicated to the determination of the petrophysical-property changes in the presence of reactive fluid, SPE J. (2010), p. 15

[Allaire et al., 2014] G. Allaire; R. Brizzi; J.F. Dufrêche; A. Mikelic; A. Piatnitski Role of non-ideality for the ion transport in porous media: derivation of the macroscopic equations using upscaling, Physica D (2014) | arXiv

[Allaire et al., 2010] G. Allaire; A. Mikelic; A. Piatnitski Homogenization of the linearized ionic transport equations in rigid periodic porous media, J. Math. Phys., Volume 51 (2010), p. 123103

[ANDRA, 2005] ANDRA Évaluation de la faisabilité du stockage géologique en formation argileuse, France (2005)

[Békri et al., 2005] S. Békri; C. Laroche; O. Vizika Pore network models to calculate transport and electrical properties of single or dual-porosity rocks, SCA, 2005, p. 2005

[Bickmore et al., 2003] B.R. Bickmore; K.M. Rosso; K.L. Nagy; R.T. Cygan; C.J. Tadanier Ab initio determination of edge surface structures for dioctahedral 2:1 phyllosilicates: implications for acid–base reactivity, Clays Clay Miner., Volume 51 (2003), pp. 359-371

[Bickmore et al., 2006] B.R. Bickmore; K.M. Rosso; C.J. Tadanier; E.J. Bylaska; D. Doud Bond-valence methods for pKa prediction. II. Bond-valence, electrostatic, molecular geometry, and solvation effects, Geochim. Cosmochim. Acta, Volume 70 (2006), pp. 4057-4071

[Blunt, 2001] M.J. Blunt Flow in porous media pore-network models and multiphase flow, Curr. Opin. Colloid Interface Sci., Volume 6 (2001), pp. 197-207

[Bocquet and Charlaix, 2010] L. Bocquet; E. Charlaix Nanofluidics, from bulk to interfaces, Chem. Soc. Rev., Volume 39 (2010), pp. 1073-1095

[Boek and Sprik, 2003] E.S. Boek; M. Sprik Ab initio molecular dynamics study of the hydration of a sodium smectite clay, J. Phys. Chem. B, Volume 107 (2003), pp. 3251-3256

[Boek and Venturoli, 2010] E.S. Boek; M. Venturoli Lattice-Boltzmann studies of fluid flow in porous media with realistic rock geometries, Comput. Math. Appl., Volume 59 (2010), pp. 2305-2314

[Boek et al., 1995] E.S. Boek; P.V. Coveney; N.T. Skipper Monte Carlo molecular modeling studies of hydrated Li-, Na-, and K-smectites: understanding the role of potassium as a clay swelling inhibitor, J. Am. Chem. Soc., Volume 117 (1995), pp. 12608-12617

[Botan et al., 2013] A. Botan; V. Marry; B. Rotenberg; P. Turq; B. Noetinger How electrostatics influences hydrodynamic boundary conditions: Poiseuille and electro-osmostic flows in clay nanopores, J. Phys. Chem. C, Volume 117 (2013), pp. 978-985

[Botan et al., 2010] A. Botan; B. Rotenberg; V. Marry; P. Turq; B. Noetinger Carbon dioxide in montmorillonite clay hydrates: thermodynamics, structure, and transport from molecular simulation, J. Phys. Chem. C, Volume 114 (2010), pp. 14962-14969

[Botan et al., 2011] A. Botan; B. Rotenberg; V. Marry; P. Turq; B. Noetinger Hydrodynamics in clay nanopores, J. Phys. Chem. C, Volume 115 (2011), pp. 16109-16115

[Bourg et al., 2007] I.C. Bourg; G. Sposito; A.C.M. Bourg Modeling the acid–base surface chemistry of montmorillonite, J. Colloid Interface Sci., Volume 312 (2007), pp. 297-310

[Bradbury and Baeyens, 2003] M. Bradbury; B. Baeyens Near Field Sorption Data Bases for Compacted MX-80 Bentonite for Performance Assessment of a High-Level Radioactive Waste Repository in Opalinus Clay Host Rock, Paul Scherrer Institut, Switzerland, 2003

[Brisard et al., 2012] S. Brisard; R.S. Chae; I. Bihannic; L. Michot; P. Guttmann; J. Thieme; G. Schneider; P.J.M. Monteiro; P. Levitz Morphological quantification of hierarchical geomaterials by X-ray nano-CT bridges the gap from nano to micro length scales, Am. Mineral., Volume 97 (2012), pp. 480-483

[Brisard and Dormieux, 2010] S. Brisard; L. Dormieux FFT-based methods for the mechanics of composites: a general variational framework, Comput. Mater. Sci., Volume 49 (2010), pp. 663-671

[Capuani et al., 2004] F. Capuani; I. Pagonabarraga; D. Frenkel Discrete solution of the electrokinetic equations, J. Chem. Phys., Volume 121 (2004), pp. 973-986

[Carof et al., 2014] A. Carof; V. Marry; M. Salanne; J.P. Hansen; P. Turq; B. Rotenberg Coarse graining the dynamics of nano-confined solutes: the case of ions in clays, Mol. Simul., Volume 40 (2014), pp. 237-244

[Carrier et al., 2014] B. Carrier; M. Vandamme; R.J.M. Pellenq; H. Van Damme Elastic properties of swelling clay particles at finite temperature upon hydration, J. Phys. Chem. C (2014) | DOI

[Churakov, 2006] S.V. Churakov Ab initio study of sorption on pyrophyllite: structure and acidity of the edge sites, J. Phys. Chem. B, Volume 110 (2006), pp. 4135-4146

[Churakov, 2007] S.V. Churakov Structure and dynamics of the water films confined between edges of pyrophyllite: a first principle study, Geochim. Cosmochim. Acta, Volume 71 (2007), pp. 1130-1144

[Coelho et al., 1996] D. Coelho; M. Shapiro; J.F. Thovert; P.M. Adler Electroosmotic phenomena in porous media, J. Colloid Interface Sci., Volume 181 (1996), pp. 169-190

[Cygan et al., 2004] R.T. Cygan; J.J. Liang; A.G. Kalinichev Molecular models of hydroxide, oxyhydroxide, and clay phases and the development of a general force field, J. Phys. Chem. B, Volume 108 (2004), pp. 1255-1266

[Cygan et al., 2012] R.T. Cygan; V.N. Romanov; E.M. Myshakin Molecular simulation of carbon dioxide capture by montmorillonite using an accurate and flexible force field, J. Phys. Chem. C, Volume 116 (2012), pp. 13079-13091

[DePaolo and Cole, 2013] D.J. DePaolo; D.R. Cole Geochemistry of geologic carbon sequestration: an overview, Rev. Mineral. Geochem., Volume 77 (2013), pp. 1-14

[van Dijke and Sorbie, 2002] M.I.J. van Dijke; K.S. Sorbie Pore-scale network model for three-phase flow in mixed-wet porous media, Phys. Rev. E, Volume 66 (2002), p. 046302

[Dufrêche et al., 2005] J.F. Dufrêche; V. Marry; N. Malikova; P. Turq Molecular hydrodynamics for electro-osmosis in clays: from Kubo to Smoluchowski, J. Mol. Liq., Volume 118 (2005), pp. 145-153

[Dufrêche et al., 2010] J.F. Dufrêche; B. Rotenberg; V. Marry; P. Turq Bridging molecular and continuous descriptions: the case of dynamics in clays, An. Acad. Bras. Cienc., Volume 82 (2010), pp. 61-68

[Fatt, 1956] I. Fatt The network model of porous media. Part I. Capillary characteristics, Petrol. Trans. AIME (1956), pp. 144-159

[Ferrage et al., 2011] E. Ferrage; B.A. Sakharov; L.J. Michot; A. Delville; A. Bauer; B. Lanson; S. Grangeon; G. Frapper; M. Jimenez-Ruiz; G.J. Cuello Hydration properties and interlayer organization of water and ions in synthetic Na-smectite with tetrahedral layer charge. Part 2. Toward a precise coupling between molecular simulations and diffraction data, J. Phys. Chem. C, Volume 115 (2011), pp. 1867-1881

[Gaus, 2010] I. Gaus Role and impact of CO2-rock interactions during CO2 storage in sedimentary rocks, Int. J. Greenh. Gas Control, Volume 4 (2010), pp. 73-89

[Giesting et al., 2012a] P. Giesting; S. Guggenheim; A.F. Koster van Groos; A. Busch Interaction of carbon dioxide with Na-exchanged montmorillonite at pressures to 640 bars: implications for CO2 sequestration, Int. J. Greenh. Gas Control, Volume 8 (2012), pp. 73-81

[Giesting et al., 2012b] P. Giesting; S. Guggenheim; A.F. Koster van Groos; A. Busch X-ray diffraction study of K- and Ca-exchanged montmorillonites in CO2 atmospheres, Environ. Sci. Technol., Volume 46 (2012), pp. 5623-5630

[Greathouse and Cygan, 2006] J.A. Greathouse; R.T. Cygan Water structure and aqueous uranyl(VI) adsorption equilibria onto external surfaces of beidellite, montmorillonite, and pyrophyllite: results from molecular simulations, Environ. Sci. Technol., Volume 40 (2006), pp. 3865-3871

[Gupta et al., 2006] A.K. Gupta; D. Coelho; P.M. Adler Electroosmosis in porous solids for high zeta potentials, J. Colloid Interface Sci., Volume 303 (2006), pp. 593-603

[Gupta et al., 2008] A. Gupta; D. Coelho; P. Adler Universal electro-osmosis formulae for porous media, J. Colloid Interface Sci., Volume 319 (2008), pp. 549-554

[Hamm et al., 2013] L.M. Hamm; I.C. Bourg; A.F. Wallace; B. Rotenberg Molecular simulation of CO2- and CO3-brine-mineral systems, Rev. Mineral. Geochem., Volume 77 (2013), pp. 189-228

[Hensen and Smit, 2002] E.J.M. Hensen; B. Smit Why clays swell, J. Phys. Chem. B, Volume 106 (2002), pp. 12664-12667

[Hiemstra et al., 1996] T. Hiemstra; P. Venema; W.H.V. Riemsdijk Intrinsic proton affinity of reactive surface groups of metal (hydr)oxides: the bond valence principle, J. Colloid Interface Sci., Volume 184 (1996), pp. 680-692

[Jardat et al., 2009] M. Jardat; J.F. Dufrêche; V. Marry; B. Rotenberg; P. Turq Salt exclusion in charged porous media: a coarse-graining strategy in the case of montmorillonite clays, Phys. Chem. Chem. Phys., Volume 11 (2009), pp. 2023-2033

[Levitz, 2007] P. Levitz Toolbox for 3D imaging and modeling of porous media: relationship with transport properties, Cement Concr. Res., Volume 37 (2007), pp. 351-359

[Liu et al., 2013] X. Liu; X. Lu; M. Sprik; J. Cheng; E.J. Meijer; R. Wang Acidity of edge surface sites of montmorillonite and kaolinite, Geochim. Cosmochim. Acta, Volume 117 (2013), pp. 180-190

[Machesky et al., 2008] M.L. Machesky; M. Predota; D.J. Wesolowski; L. Vlcek; P.T. Cummings; J. Rosenqvist; M.K. Ridley; J.D. Kubicki; A.V. Bandura; N. Kumar; J.O. Sofo Surface protonation at the rutile (1 1 0) interface: explicit incorporation of solvation structure within the refined MUSIC model framework, Langmuir, Volume 24 (2008), pp. 12331-12339

[Marino et al., 2001] S. Marino; M. Shapiro; P. Adler Coupled transports in heterogeneous media, J. Colloid Interface Sci., Volume 243 (2001), pp. 391-419

[Marry et al., 2002] V. Marry; P. Turq; T. Cartailler; D. Levesque Microscopic simulation of structure and dynamics of water and counterions in a monohydrated montmorillonite, J. Chem. Phys., Volume 117 (2002), p. 3454

[Marry et al., 2008] V. Marry; B. Rotenberg; P. Turq Structure and dynamics of water at a clay surface from molecular dynamics simulation, Phys. Chem. Chem. Phys., Volume 10 (2008), pp. 4802-4813

[Michot et al., 2012] L.J. Michot; E. Ferrage; M. Jimnez-Ruiz; M. Boehm; A. Delville Anisotropic features of water and ion dynamics in synthetic Na- and Ca-smectites with tetrahedral layer charge. A combined quasi-elastic neutron-scattering and molecular dynamics simulations study, J. Phys. Chem. C, Volume 116 (2012), pp. 16619-16633

[Moyne and Murad, 2006a] C. Moyne; M.A. Murad A two-scale model for coupled electro-chemo-mechanical phenomena and Onsager's reciprocity relations in expansive clays. I. Homogenization analysis, Transp. Porous Media, Volume 62 (2006), pp. 333-380

[Moyne and Murad, 2006b] C. Moyne; M.A. Murad A two-scale model for coupled electro-chemo-mechanical phenomena and Onsager's reciprocity relations in expansive clays. II. Computational validation, Transp. Porous Media, Volume 63 (2006), pp. 13-56

[Myshakin et al., 2013] E.M. Myshakin; W.A. Saidi; V.N. Romanov; R.T. Cygan; K.D. Jordan Molecular dynamics simulations of carbon dioxide intercalation in hydrated Na-montmorillonite, J. Phys. Chem. C, Volume 117 (2013), pp. 11028-11039

[Obliger et al., 2013] A. Obliger; M. Duvail; M. Jardat; D. Coelho; S. Békri; B. Rotenberg Numerical homogenization of electrokinetic equations in porous media using lattice-Boltzmann simulations, Phys. Rev. E, Volume 88 (2013), p. 013019

[Obliger et al., 2014] A. Obliger; M. Jardat; D. Coelho; S. Békri; B. Rotenberg Pore network model of electrokinetic transport through charged porous media, Phys. Rev. E, Volume 89 (2014), p. 043013

[Pagonabarraga et al., 2005] I. Pagonabarraga; F. Capuani; D. Frenkel Mesoscopic lattice modeling of electrokinetic phenomena, Comput. Phys. Commun., Volume 169 (2005), pp. 192-196

[Pagonabarraga et al., 2010] I. Pagonabarraga; B. Rotenberg; D. Frenkel Recent advances in the modelling and simulation of electrokinetic effects: bridging the gap between atomistic and macroscopic descriptions, Phys. Chem. Chem. Phys., Volume 12 (2010), pp. 9566-9580

[Porion et al., 2013] P. Porion; A.M. Faugère; A. Delville Multiscale water dynamics within dense clay sediments probed by 2H multiquantum NMR relaxometry and two-time stimulated echo NMR spectroscopy, J. Phys. Chem. C, Volume 117 (2013), pp. 26119-26134

[Pride, 1994] S. Pride Governing equations for the coupled electromagnetics and acoustics of porous media, Phys. Rev. B, Volume 50 (1994), pp. 15678-15696

[Robinet et al., 2012] J.C. Robinet; P. Sardini; D. Coelho; J.C. Parneix; D. Prêt; S. Sammartino; E. Boller; S. Altmann Effects of mineral distribution at mesoscopic scale on solute diffusion in a clay-rich rock: example of the Callovo-Oxfordian mudstone (Bure, France), Water Resour. Res., Volume 48 (2012), p. W05554

[Rotenberg and Pagonabarraga, 2013] B. Rotenberg; I. Pagonabarraga Electrokinetics: insights from simulation on the microscopic scale, Mol. Phys., Volume 111 (2013), pp. 827-842

[Rotenberg et al., 2007a] B. Rotenberg; V. Marry; J.F. Dufrêche; E. Giffaut; P. Turq A multiscale approach to ion diffusion in clays: building a two-state diffusion–reaction scheme from microscopic dynamics, J. Colloid Interface Sci., Volume 309 (2007), pp. 289-295

[Rotenberg et al., 2007b] B. Rotenberg; V. Marry; J.F. Dufrêche; N. Malikova; E. Giffaut; P. Turq Modelling water and ion diffusion in clays: a multiscale approach, C. R. Chim., Volume 10 (2007), pp. 1108-1116

[Rotenberg et al., 2010a] B. Rotenberg; V. Marry; N. Malikova; P. Turq Molecular simulation of aqueous solutions at clay surfaces, J. Phys.: Condens. Matter, Volume 22 (2010), p. 284114

[Rotenberg et al., 2007c] B. Rotenberg; V. Marry; R. Vuilleumier; N. Malikova; C. Simon; P. Turq Water and ions in clays: unraveling the interlayer/micropore exchange using molecular dynamics, Geochim. Cosmochim. Acta, Volume 71 (2007), pp. 5089-5101

[Rotenberg et al., 2009] B. Rotenberg; J.P. Morel; V. Marry; P. Turq; N. Morel-Desrosiers On the driving force of cation exchange in clays: insights from combined microcalorimetry experiments and molecular simulation, Geochim. Cosmochim. Acta, Volume 73 (2009), pp. 4034-4044

[Rotenberg et al., 2008] B. Rotenberg; I. Pagonabarraga; D. Frenkel Dispersion of charged tracers in charged porous media, EPL, Volume 83 (2008), p. 34004

[Rotenberg et al., 2010b] B. Rotenberg; I. Pagonabarraga; D. Frenkel Coarse-grained simulations of charge, current and flow in heterogeneous media, Faraday Discuss., Volume 144 (2010), pp. 223-243

[Rotenberg et al., 2011] B. Rotenberg; A.J. Patel; D. Chandler Molecular explanation for why talc surfaces can be both hydrophilic and hydrophobic, J. Am. Chem. Soc., Volume 133 (2011), pp. 20521-20527

[Sposito et al., 1999] G. Sposito; N.T. Skipper; R. Sutton; S.h. Park; A.K. Soper; J.A. Greathouse Surface geochemistry of the clay minerals, Proc. Natl. Acad. Sci. U.S.A., Volume 96 (1999), pp. 3358-3364

[Stukan et al., 2010] M.R. Stukan; P. Ligneul; J.P. Crawshaw; E.S. Boek Spontaneous imbibition in nanopores of different roughness and wettability, Langmuir, Volume 26 (2010), pp. 13342-13352

[Sulpizi and Sprik, 2010] M. Sulpizi; M. Sprik Acidity constants from DFT-based molecular dynamics simulations, J. Phys.: Condens. Matter, Volume 22 (2010), p. 284116

[Suter et al., 2008] J.L. Suter; E.S. Boek; M. Sprik Adsorption of a sodium ion on a smectite clay from constrained ab initio molecular dynamics simulations, J. Phys. Chem. C, Volume 112 (2008), pp. 18832-18839

[Tazi et al., 2012] S. Tazi; B. Rotenberg; M. Salanne; M. Sprik; M. Sulpizi Absolute acidity of clay edge sites from ab-initio simulations, Geochim. Cosmochim. Acta, Volume 94 (2012), pp. 1-11

[Teppen and Miller, 2006] B.J. Teppen; D.M. Miller Hydration energy determines isovalent cation exchange selectivity by clay minerals, Soil Sci. Soc. Am. J., Volume 70 (2006), p. 31

[Tournassat et al., 2004] C. Tournassat; E. Ferrage; C. Poinsignon; L. Charlet The titration of clay minerals. II. Structure-based model and implications for clay reactivity, J. Colloid Interface Sci., Volume 273 (2004), pp. 234-246

[Tournassat et al., 2013] C. Tournassat; S. Grangeon; P. Leroy; E. Giffaut Modeling specific pH dependent sorption of divalent metals on montmorillonite surfaces. A review of pitfalls, recent achievements and current challenges, Am. J. Sci., Volume 313 (2013), pp. 395-451

[Tyagi et al., 2013] M. Tyagi; T. Gimmi; S.V. Churakov Multi-scale micro-structure generation strategy for up-scaling transport in clays, Adv. Water Resour., Volume 59 (2013), pp. 181-195

[White and Zelazny, 1988] G.N. White; L.W. Zelazny Analysis and implications of the edge structure of dioctahedral phyllosilicates, Clays Clay Miner., Volume 36 (1988), pp. 141-146


Commentaires - Politique