Benoît Claire Benjamin

Guillaume

Manu Lucie

Keywords:

Une thèse n'est jamais le travail d'une seule personne et celle-ci ne fait pas exception. J'ai eu la chance durant ces travaux d'être encadré par une équipe efficace et toujours présente pour répondre à mes questions et qui m'a énormément apporté tant sur le plan scientifique qu'humain. Je voudrais ainsi remercier Teddy Fen-Chong mon maître de thèse pour sa patience, ses bons conseils et son efficacité, Patrick Dangla pour toutes ces discussions acharnées, Jean-Michel Pereira pour son calme et pour m'avoir tenu compagnie durant toutes ces conférences et réunions, Antonin Fabbri pour sa bonne humeur et Arnault Lassin qui a repris le projet SALTCO au pied levé tout en répondant toujours présent pour les questions de géochimie. Je tiens aussi à remercier Laurent André qui a lui aussi toujours pris le temps de répondre à mes nombreuses interrogations sur le transport réactif.

Ces travaux été possibles grâce aux financement de l'ADEME à travers le projet SALTCO (no. 1094C0013), et de la chaire CTSC (Captage Transport et Stockage du CO 2 ) que je remercie sincèrement.

Je souhaite également remercier les deux rapporteurs de cette thèse Jan Carmeliet et Lionel Mercury ainsi que les autres membres du jury, Gilles-Pijaudier-Cabot et François Renard pour m'avoir fait l'honneur de lire ce manuscrit et d'avoir participé à ma soutenance.

Cette thèse a aussi été l'histoire d'une longue lutte avec ce prototype de percolation réactive. Je tiens à remercier David Hautemayou pour son soutien technique et Cédric Mezière pour m'avoir laissé lui emprunter tous ses outils. J'ai eu aussi la chance de travailler dans un environnement particulièrement chaleureux, et pour cela je tiens à remercier toutes les personnes de Képler, en particulier les thésards : mes cobureaux Adrien, Mathilde et Antoine,

Résumé

Dans un contexte de réduction internationale des émissions de gaz à effet de serre, les techniques de Captage Transport et Stockage de CO 2 (CTSC) apparaissent comme une solution à moyen terme particulièrement efficace. En effet, les capacités de stockage géologique pourraient s'élever jusqu'à plusieurs millions de tonnes de CO 2 injectées par an, soit une réduction substantielle des émissions atmosphériques de ce gaz. Une des cibles privilégiées pour la mise en place de cette solution sont les aquifères salins profonds. Ces aquifères sont des formations géologiques contenant une saumure dont la salinité est souvent supérieure à celle de la mer la rendant impropre à la consommation. Cependant, cette technique fait face à de nombreux défis technologiques ; en particulier la précipitation des sels, dissous dans l'eau présente initialement dans l'aquifère cible, suite à son évaporation par le CO 2 injecté. Les conséquences de cette précipitation sont multiples, mais la plus importante est une modification de l'injectivité, c'est-à-dire des capacités d'injection. La connaissance de l'influence de la précipitation sur l'injectivité est particulièrement importante tant au niveau de l'efficacité du stockage et de l'injection qu'au niveau de la sécurité et de la durabilité du stockage. Le but de ces travaux de thèse est de comparer l'importance relative des phénomènes négatif (colmatage) et positif (fracturation) consécutifs à l'injection de CO 2 et à la précipitation des sels. Au vu des nombreux résultats de simulations et de modélisation dans la littérature décrivant le colmatage de la porosité, il a été décidé de porter l'accent sur les effets mécaniques de la cristallisation des sels et la possible déformation de la roche mère. Une modélisation macroscopique et microscopique, tenant compte de deux modes possibles d'évaporation induits par la distribution spatiale de l'eau résiduelle a donc été développée afin de prédire le comportement mécanique d'un matériau poreux soumis à un assèchement par injection de CO 2 . Les résultats montrent que la pression de cristallisation consécutive à la croissance d'un cristal en milieu confiné peut atteindre des valeurs susceptibles localement de dépasser la résistance mécanique du matériau, soulignant ainsi l'importance de ces phénomènes dans le comportement mécanique global de l'aquifère. Sur le plan expérimental, les travaux ont porté sur l'utilisation d'un nouveau prototype de percolation réactive afin de reproduire le comportement d'une carotte de roche soumise à l'injection et ainsi obtenir l'évolution des perméabilités dans des conditions similaires à celle d'un aquifère. annual mean, in red the 5-year running mean, in green are the uncertainty estimates) . .

8.5

Crystal growing in a spherical pore of radius r p and invading the small entry channels of radius r s . The curvature of the crystal in the channel depends on the current supersaturation (adapted from [START_REF] Scherer | Crystallization in pores[END_REF]). . . . . . . . . 8.6 Particular case of a crystal growing in a cylindrical pore, with tips radius of curvature corresponding to the pore radius (adapted from [START_REF] Scherer | Crystallization in pores[END_REF]). . . . .

List of symbols

General introduction

In the context of increasing environmental concern, CCS (Carbon Capture, transport and geological Storage) appears to be an attractive solution to cope with global warming and increasing carbon dioxide content in the atmosphere. The aim of CCS is to capture carbon dioxide at the exit of the emitters and to store it in the subsurface (offshore or onshore) where it will no longer contribute to greenhouse effect [START_REF] Pacala | Stabilisation wadges: solving the climate problem for the next 50 years with current technologies[END_REF]. There are several targets for CCS such as depleted oil and gas reservoirs, unmined coal veins or deep saline aquifers. Deep saline aquifers are huge geological formations containing saline water, unfit for human consumption. The idea is then to inject the carbon dioxide within the porosity of the aquifer and to store it for an indefinite amount of time. This solution seems particularly interesting as storage capacities of deep saline aquifers can raise to several millions of tons of carbon dioxide stored per year [START_REF] Bachu | Sequestration of CO2 in geological media in response to climate change: Capacity of deep saline aquifers to sequester CO2 in solution[END_REF]. According to the estimations of the International Panel on Climate Change (IPCC [START_REF] Pachauri | Climate change 2007: Synthesis report[END_REF]), in order to keep the atmospheric concentration of carbon dioxide around the values of 2005, the reduction of emission has to be as important as 15 to 25 Gt CO 2 /year, making CCS one of the key solution for reducing carbon dioxide emissions. However, numerous studies on CCS have shown that injection of supercritical carbon dioxide in a deep saline aquifer disturbs dramatically the equilibrium of the medium, leading to highly coupled THMC (Thermal-Hydrodynamical-Mechanical-Chemical) behavior [4]. In particular, carbon dioxide injection is known to trigger strong salt precipitation in the zone close to the injection well, causing a strong decrease in the porosity and the permeability of the host rock [START_REF] Pruess | Formation dry-out from CO2 injection into saline aquifers. part 2, analytical model for salt precipitation[END_REF][START_REF] Pruess | Formation dry-out from CO2 injection into saline aquifers. part 1, effects of solids precipitation and their mitigation[END_REF]. Indeed, the addition of solid matter in the porosity clogs the percolation paths for carbon dioxide, lowering the storage capacity of the aquifer. Ultimately it has been shown that precipitation of salt during evaporation by a flow of supercritical carbon dioxide can completely clog the medium, stopping completely the injection operations and thus the carbon dioxide storage [START_REF] Peysson | Permeability alteration induced by drying of brines in porous media[END_REF]. However, precipitation of salt in a confined medium is also known to create strong stresses on the rock called crystallization pressure. Crystallization pressure has long been recognized as a major source of damage to concrete or masonry buildings [START_REF] Espinosa-Marzal | Advances in understading damage by salt crystallization[END_REF]. Numerous examples of degradation of construction rocks or concrete by salt can be easily found, be this latter brought by salt sprays of by capillary rise. Literature on crystallization pressure is more than one century rich of experimental and theoretical works [9][10][11] and even if several interpretations still co-

GENERAL INTRODUCTION

exist, the commonly accepted explanation for salt weathering relies on a thermodynamical modelling: stress on the rock is created by the constrained growth of a non-wetting crystal within the pore network. The presence of the pore wall, limiting the crystal's growth modifies the chemical equilibrium and induces the build up of a disjoining pressure that reestablishes the mechanical equilibrium. Crystallization pressure is then very likely to result from the precipitation of salt during CCS. The stresses created may be able to fracture the host rock enabling new percolation paths for the carbon dioxide, but also may raise security issues about leakage and fault reactivation [12].

This PhD work, supported by the Chaire CTSC (Captage Transport et Stockage du CO 2 ), is part of the SALTCO project (No. 1094C0013) supported by the BRGM (French Geological Survey), IFSTTAR (Institut Français des Sciences et Technologies des Transport, de l'Aménagement et des Réseaux), and by the French Environment and Energy Management Agency (ADEME) for a duration of 42 months. This project aims at determining the balance between the negative effects (clogging) and the positive effects (percolation path opening and matrix dissolution) of the crystallization of salt following the brine evaporation of carbon dioxide. As the clogging part was well developed in the literature, the accent has been put on the determination of the mechanical effect of salt crystallization and thus on the estimation of crystallization pressure. To support the modelling of this crystallization pressure and to observe with in situ conditions the drying of the porous medium, this work included an experimental part on an innovative prototype. The aim of these experiments was to quantify the permeability variation and to measure the deformation of a rock plug filled with brine and subjected to carbon dioxide injection and finally to compare these experimental results with the modelling of crystallization pressure. However, because of technical issues, this last part was not completely proceeded at the moment of the redaction of this thesis.

In a first part, we review the political and scientific contexts of CCS. The second part is dedicated to the definitions and the bibliography review. The first chapter describes the hydrodynamical behavior of the aquifer, the second the chemical behavior and the third the thermal behavior. Modelling of crystallization pressure and its estimation in the case of CCS is the subject of the third part. Three different models are presented, two macroscopic and one microscopic, based on different hypotheses and each of them highlighting characteristics and issues of the crystallization pressure estimation. Finally, the last part is the experimental part, presenting the prototype of reactive percolation and describing the different experiments carried-out, such as intrinsic/relative permeability measurement. Because of the novelty of the reactive percolation prototype, we faced several technical issues which are (along with the adopted solutions) described at the end of the experimental part. This increase of emissions is the main cause of an increase of the concentration of carbon dioxide in the atmosphere. It is now largely accepted that this increase of concentration is the cause of an augmentation of the intensity of the greenhouse effect (figure 1.2). Greenhouse effect is a consequence of the absorption of the thermal radiations from the Earth surface by the greenhouse gases: after absorption, the energy is reemitted in all directions and particularly back toward the earth and as a result, an important part of the thermal energy is kept in the atmosphere and is not lost in outer space. This effect is one of the phenomena which make Earth a suitable planet for life by maintaining a temperature allowing water to remain liquid. When the carbon dioxide molecule is hit by a thermal radiation from the Earth, it vibrates strongly because the frequency of this radiation is close to the resonance frequency of the molecules, and reemits this thermal radiation in all directions. This is the fact that the resonance of the carbon dioxide molecule is close to the frequency of the radiation emitted from the surface of the Earth which makes it a potent greenhouse gas. The increase of carbon dioxide content in the atmosphere causes then an augmentation of greenhouse effect and thus an augmentation of temperature, leading to what is called global warming: as thermal energy is more and more retained in the atmosphere, global temperature at the surface of the Earth is increasing as represented in figure 1.3. This figure shows the variation of temperature from 1880 to 2012 relative to the mean temperature from 1951-1980. The current warming is 1.1. CARBON DIOXIDE AND GREENHOUSE EFFECT about 0.6 ℃ but no slowing down is expected. The most pessimistic scenarios expect a warming up to 6.4 ℃ [START_REF]Special Report on Carbon Dioxide Capture and Storage[END_REF]. Consequences of this global warming are numerous and several scenarios foresee an augmentation of tornados, drought and floods along with a global augmentation of the sea level, drowning huge parts of lands particularly in Bangladesh, The Netherlands, or some Indian Ocean islands such as Seychelles or Maldives [START_REF] Meehl | Global climate projections[END_REF]. Anthropic carbon dioxide comes from different sources but as can be seen from figure 1.1, combustion of fossil energies and cement industry are the strongest emitters of this gas. Fossil energy combustion includes road and sea transport, thermal powerplants, heating. . . The cement industry produces carbon dioxide in two ways: the chemical reaction itself (CaCO 3 → CaO + CO 2 ) and the energy needed to produce the clinker (heating at 1450 ℃).

Part I Injection of supercritical carbon dioxide in deep saline aquifers

Political decisions and solutions 1.2.1 International conferences and political decisions

Because of the dramatic consequences of global warming, for both the environment and humankind, the United Nations decided to create a framework to cope with the increasing carbon dioxide emissions: the FCCC (Framework Convention on Climate Change), at the "Earth summit" in Rio de Janeiro during June 1992. The aim of this framework is to "stabilize greenhouse gas concentrations in the atmosphere at a level that would prevent dangerous anthropogenic interference with the climate system" [START_REF]The United Nation Framework Convention on Climate Change[END_REF]. Under this framework, several international conventions were held, the most famous is the Kyoto convention in which was signed the Kyoto protocol in December 11th 1997. The purpose of this protocol is to reduce the emissions of six greenhouse gases (carbon dioxide, methane, nitrous oxide and three substitutes of chlorofluorocarbons), by 5.2 % between 2008 and 2012 compared to the emissions of 1990. This protocol was signed by 192 countries but with several levels of commitment (for example the USA have signed the treaty but with no intention of ratifying, while the Canada has withdrawn from it in 2011).

The different ways to decrease carbon dioxide emissions

Besides the different policies of efficient energy use, one solution to decrease the emissions is to capture the carbon dioxide at the exit of the emitters and to store it underground where it will no longer contribute the greenhouse effect and global warming. This solution is called CCS (Carbon Capture and Storage) [START_REF] Pacala | Stabilisation wadges: solving the climate problem for the next 50 years with current technologies[END_REF]. The targets for storage as represented in figure 1.4 are mainly geological formations which can be depleted gas and oil fields (in this case CCS can also be combined with enhanced oil recovery (EOR) where carbon dioxide is used to increase the oilfield pressure and recover more oil, the carbon dioxide staying within the aquifer), unmined coal beds, salt formations. . . and finally deep saline aquifers. These latter geological formations are aquifers containing saline water whose salinity make it unsuitable for human consumption. These formations are often immense

Chapter 2

Problem at stake: injectivity and permeability evolution

One of the key parameters for an efficient and durable injection is the injectivity. This parameter describes the quantity of carbon dioxide it is possible to inject per unit of time in the aquifer. Injectivity depends on the petrophysical characteristics of the medium such as porosity and permeability (see section 3.2.2) and on the characteristics of the injected fluid (density, viscosity) which depend particularly on its state: supercritical, liquid or gaseous. A control of the injectivity is then particularly important in order to keep a high storage and injection capacity of the considered aquifer. Moreover, knowledge of the injectivity is mandatory in order to avoid any strong mechanical response of the aquifer such as microseismicity or fault reactivation. Injectivity during injection of carbon dioxide can evolve through three main mechanisms: relative permeability, matrix dissolution, and finally clogging and fracturation. Relative permeability has strong influence on the transport capacities of the medium. Indeed, depending on the relative occupancy of the medium, the fluids can move more or less efficiently thus impacting the injectivity. The second effect is matrix dissolution: the acidification of the brine by the dissolution of carbon dioxide disturbs the chemical equilibrium of the aquifer leading to strong dissolution of pH sensitive minerals like carbonates. This dissolution can locally enhance the permeability and thus increase injectivity. However, dissolved minerals can also reprecipitate further in the aquifer following the chemical conditions. Finally, the last effect is the precipitation of salt within the porosity. The brine originally filling the aquifer contains a high content of salt, and evaporation of this brine during the injection of dry carbon dioxide (particularly in a zone close to the injection well called the dry-out zone as described in chapter 6) triggers its precipitation. This precipitation adds solid matter within the porosity which has as main consequence the clogging of the percolation paths of carbon dioxide, reducing dramatically the permeability and thus the injectivity Problem at stake: injectivity and permeability evolution (a decrease of a few percent of porosity can lead to a decrease of several order of magnitude of intrinsic permeability [START_REF] Pruess | Formation dry-out from CO2 injection into saline aquifers. part 2, analytical model for salt precipitation[END_REF][START_REF] Pruess | Formation dry-out from CO2 injection into saline aquifers. part 1, effects of solids precipitation and their mitigation[END_REF]). However, precipitation of salt is also the cause of the creation of strong stresses on the solid matrix. These come from the so-called crystallization pressure which has been intensively studied in the case of civil engineering but not in the CCS context. They may be the cause of the fracturation of the host rock raising both interesting consequences on the injectivity and security issues. Indeed, the creation of a fracture network is known to enhance strongly the permeability and thus the injectivity. On the other hand this fracturation may cause or enhance unwanted fluid leakage toward the surface threatening the storage integrity. As a result an estimation of the effect of salt and chemical equilibrium on the injectivity seems particularly important in order to obtain a whole overview of the behavior of the aquifer during injection. The initial aim of this PhD work was then to evaluate the relative importance of the described processes on the evolution of the injectivity. However, as the porosity reduction and the chemical interaction is well grasped by the current reactive transport codes (principally TOUGHREACT, section 6.1) the accent is here put on the mechanical consequences of the precipitation through the crystallization pressure.

Part II

Thermo

Porosities and phase saturations

Definitions

In the following, we consider a porous medium, consisting in a solid matrix, m and a porous (void) network. Occluded porosity (empty spaces in the rock matrix which are not connected to the main porous network) is considered as belonging to the solid matrix.

Considering the possible deformation (we will consider in all the following that we are under the hypothesis of small deformations (see section 7.1.2)) of the considered element of porous solid, we can define two different porosities n and φ; n, called the Eulerian porosity, is defined as the ratio between the current elementary volume of the empty space dΩ pores and the current elementary volume of the porous medium dΩ :

n = dΩ pores dΩ (3.1)
φ, called Lagrangian porosity, is rather defined as the ratio between the current pore volume to the initial volume of porous medium:

φ = dΩ pores dΩ 0 (3.2)
and we have the following relation between φ and n:

φ = (1 + ǫ) n (3.3)
where ǫ is the volumetric deformation (see section 7.3) of the REV (Representative Elementary Volume, see section 7.1.2).

In absence of any deformation of the porous solid (e.g. ǫ = 0), Eulerian porosity and Lagrangian porosity are equal. Considering a deformation of the medium, we can define the change of porosity ϕ as:

ϕ = φ -φ 0 (3.4)
Expression of the sole porosity change is much simpler in the Lagrangian definition than with the Eulerian one [START_REF] Coussy | Mechanics and Physics of Porous Solids[END_REF]. Accordingly, in the following, because of equation (3.4), we will consider only the Lagrangian definition of porosities and saturations.

For a system with several phases, the porosity is occupied by the different phases J (for example liquid and supercritical phases). We can then write this partitioning as:

3.1. POROSITIES AND PHASE SATURATIONS φ = J φ J (3.5)
where φ J is the partial porosity (φ J dΩ 0 being the volume currently occupied by phase J).

The saturation degree (or simply saturation) of the phase J is referred to the reference configuration:

S J = φ J φ 0 ; J S J = 1 (3.6)
It corresponds to the ratio between the volume of pores occupied by the phase J and the total initial volume of pores.

We can then define the porosity variation associated with each phase J:

ϕ J = φ J -φ 0 S J ; J ϕ J = ϕ (3.7)
In figure 3.1 are represented the different porosity quantities and variations in the case of a porous material filled with two phases (here C and L).

Diagenesis and the formation of natural porous media

Sedimentary rocks, i.e. rocks present in the aquifers are created by the compaction of sediments deposited mainly in aqueous environments (rivers, lakes, seas. . . ). These sediments CHAPTER 3. HYDRODYNAMIC BEHAVIOR OF THE AQUIFER are slowly buried under other geological layers and then are submitted to compaction and diagenesis (set of processes leading to the modification of the rock structure under the effect of pressure and temperature during burial [START_REF] Zinszner | A Geoscientist's guide to Petrophysics[END_REF]). The porosity and also the mineralogical composition change during this event. Modification of porosity consists of the rearrangement of grains (compaction) and cementation by pressure solution (see 8.1.1). The latter is particularly important in limestones because of the chemical reactivity of calcite and dolomite (cementation by dolomitization).

Porosity is then a function of the initial distribution of the size of grains, the pressure and temperature conditions in which the rock has aged, and its mineralogy. Finally, the presence of fluids within the reservoir rock and the chemical reactions can alter the rock thus modifying its porosity. These processes are the cause of the interesting storage properties of the deep saline aquifer's rocks, which make them very interesting targets for carbon dioxide storage.

Mass balance and fluid velocity

Mass Balance

Evolution of the different saturations of the moving phases can be obtained by considering the mass balance equation. It is obtained by the calculation of inputs and outputs [START_REF] Pruess | Formation dry-out from CO2 injection into saline aquifers. part 1, effects of solids precipitation and their mitigation[END_REF][START_REF] Coussy | Poromechanics[END_REF][START_REF] Rutqvist | Thermohydromechanics of partially saturated geological media: governing equations and formulation for four finite element models[END_REF] in a small representative volume of porous medium (and under infinitesimal transformations, see 7.1.2):

∂m J ∂t = -∇. ω J - • m J→K (3.8)
m J is the in-pore mass of J per unit of initial volume dΩ 0 , vector ω is the flux vector which can be written as:

ω J = m J v J = ρ J φ 0 S J v J (3.9)
with ρ J the density, S J the saturation, φ J the partial porosity and v J the velocity of phase J. The last term,

• m J→K , is the creation (sink/source) term rate implying chemical transformation and phase partitioning (see 4.1.3).

MASS BALANCE AND FLUID VELOCITY

Fluid velocity: Darcy's law and relative permeabilities

The equation governing viscous fluid velocity in a porous medium is well known and based on Darcy's law:

v J = - k 0 k r J (S J ) η J ( ∇p J -ρ J g) (3.10)
k 0 is the intrinsic permeability, k r J is the relative permeability corresponding to the current saturation of the studied fluid [START_REF] Anderson | Wettability literature survey: Part 5. the effect of wettability on relative permeability[END_REF]; η J is the dynamic viscosity and ∇p J -ρ J g the driving force producing the flow (p J is the pressure and g the gravity). IS units of intrinsic permeability is m 2 . However, is it usual use the Darcy, D, which corresponds to 9.868233.10 -13 m 2 ( ≈ 1µm 2 ). Usual values of intrinsic permeabilities for aquifers rocks range from 0.1mD to 1D for the most permeable.

Relative permeability describes the capacity to transport a fluid compared to the other fluid present in the medium. It depends strongly on the wettability of the medium. There is no real theoretical way to calculate relative permeabilities and semi-empirical relations are used such as the famous van-Genuchten/Mualem1 equation [START_REF] Van Genuchten | A closed-form equation for predicting the hydraulic conductivity of unsaturated soils[END_REF], which is adopted here for liquids:

k r l = √ S * (1 -(1 -[S * ] (1/m) ) m ) 2 (3.11)
with S * = S l -S lr 1-S lr , S l the water saturation, S lr the residual water saturation (saturation at which water, being the wetting fluid, cannot be displaced any more, no matter the flow rate of carbon dioxide), m is an exponent depending on the porous medium only.

Gas relative permeability is usually represented by the Corey [START_REF] Corey | The interrelation between gas and oil relative permeabilities[END_REF] relation:

k r g = (1 -S h ) 2 (1 -S 2 h ) (3.12)
with S h = S l -S lr 1-S lr -Sgr , S gr is the residual gas saturation i.e. the remaining gas saturation after maximum imbibition with water (as brine is the wetting fluid, residual gas saturation is often very low). Chemical reactions occurring between the brine and the reservoir rock (see chapter 4.2) can lead to a modification of these two relations by modifying the pore structure of the system and thus the coefficients m, S lr and S gr .

A measure of relative permeabilities of reservoir rocks and a fit with the two presented laws are done in chapter 13.

Displacement of dissolved species: diffusion

Diffusion is the displacement of species because of a gradient of concentration in order to homogenize the solution. It is described by Fick's law which can be written in extended form as [START_REF] Coussy | Mechanics and Physics of Porous Solids[END_REF] 2 :

m i v i = -φτ D i ∇m i (3.13)
i is the considered species, m i is the molality of this species (concentration or mole fraction can also be used, see 4.9), D i (m 2 /s) is the diffusion coefficient relative to species i in free space, φ is the porosity and τ is the tortuosity. Tortuosity is a parameter characterizing the deviation of straightness and the connectivity of the porous network.

Klinkenberg effect

We have defined the intrinsic permeability as a constant coefficient in Darcy's law corresponding to a characteristic of the porous medium. However, if one carries out measures of intrinsic permeability with liquids and gas, results will show a bigger permeability for gas as opposed to liquids. This effect is called the Klinkenberg effect [START_REF] Coussy | Mechanics and Physics of Porous Solids[END_REF] and is a consequence of the slipping of gas molecules on the pore wall which then make them go faster than the liquid molecules, thus increasing comparatively the measured intrinsic permeability. This Klinkenberg effect occurs when the Knudsen number, measuring the ratio between the mean free path of the fluid molecules and a characteristic length (here the pore diameter), is close to unity. In this case, intrinsic permeability has to be modified as follows [START_REF] Coussy | Mechanics and Physics of Porous Solids[END_REF]:

k 0 klinkenberg = k 0 1 + ̟ p (3.14)
with ̟ is a characteristic pressure depending on the gas and the porous network geometry and p the gas pressure. In the case of injection of supercritical carbon dioxide, all fluids will be in dense phase (liquid and supercritical) and will thus not be concerned by the Klinkenberg effect.

SURFACE ENERGY AND CONSEQUENCES

Surface energy and consequences

Interfacial energy and wettability

Surface energy and surface stress

By definition, interfacial energy (written σ) is the excess of energy due to the presence of an interface. Indeed, molecules at an interface have an excess of energy because of the lack of bonding with respect to the bulk molecules [START_REF] De Gennes | Capillarity and wetting phenomena: drops, bubbles, pearls, waves[END_REF]. The energy W needed to create an interface of area A between two phases 1 and 2, as represented in figure 3.2, is:

W = σ 1,2 A (3.15)
Figure 3.2: Definition of the surface energy as the creation of new surfaces.

For a fluid, because of the mobility of the molecules, any stretch of the surface will bring new molecules from the bulk to this surface and is similar to the creation of a new surface. On the contrary, for a solid, a stretch of the surface will just stretch the bonds between the atoms (or molecules) constituting the surface and will thus be different from the surface energy. This stretching of a surface is characterised by the surface stress σ * as expressed by the Shuttleworth equation:

σ * = σ + dσ dǫ (3.16)
σ is the surface energy and ǫ is the strain of the considered surface, σ * = σ for a fluid-fluid interface.

Interfacial energy is really important in both transport and equilibrium situations. Indeed, it plays an important role on capillary pressure (see section 3.3.2) and on crystal CHAPTER 3. HYDRODYNAMIC BEHAVIOR OF THE AQUIFER equilibrium (see chapter 8). Values of interfacial energy between brine and carbon dioxide are tabulated3 but can evolve depending on the conditions (ionic strength [START_REF] Leroy | Predicting the surface tension of aqueous 1:1 electrolyte solutions at high salinity[END_REF], pressure, temperature). However, in order to simplify our modelling, we will not consider this effect on interfacial energies.

Wettability

Wettability is a consequence of the interfacial energy [START_REF] Anderson | Wettability literature survey: Part 1. rock-oil-brine interactions and the effect of core handling on wettability[END_REF] and defines the ability of a fluid to spread (wet) on a solid with respect to another fluid. Indeed because of the interfacial energies, contact between solid and fluid may be favorable or not favorable. If the contact is energetically favorable, the contact angle will be low (under 90°and even 0°for perfect wetting). On the contrary, if the contact is not favorable, the considered fluid will not spread on the surface and contact angle will be above 90°and can even reach 180°for complete non-wetting (see figure 3.3). The contact angle θ (here between the phase 2 represented in blue in figure 3.3 and the solid s) can be related to the interfacial energies according to the Young-Dupré law:

σ s,1 = σ 1,2 cos θ + σ 2,s (3.17)
In saline aquifers, rocks are water-wet [START_REF] Anderson | Wettability literature survey: Part 1. rock-oil-brine interactions and the effect of core handling on wettability[END_REF], and carbon dioxide injection is then called drainage. Wettability influences the irreducible saturations [START_REF] Al Mansoori | Measurement of non-wetting phase trapping applied to carbon dioxide storage[END_REF] and also the pattern of residual water, depending on the flow regime [START_REF] Lenormand | Capillary fingering -percolation and fractal dimension[END_REF] (see section 3.4). Pattern of residual water has huge implications on the crystallization in the dry-out zone because 3.3. SURFACE ENERGY AND CONSEQUENCES the distribution of crystal in the pores will be strongly different if residual water is trapped as clusters or as capillary film on the pore walls (see section 8.2).

Capillary pressure

Capillary pressure is commonly defined as the pressure difference between the non-wetting and the wetting fluid (see for example [START_REF] Anderson | Wettability literature survey: Part 4. effects of wettability on capillary pressure[END_REF]).

p cap = p nw -p w (3.18)
nw refers to the non-wetting fluid and w to the wetting one (in our case, carbon dioxide is the non-wetting fluid and brine the wetting one). Capillary pressure can have a microscopic interpretation related to the surface energies. Indeed, as represented in figure 3.4, two immiscible phases are separated by an interface whose curvature κ is linked to the variation of pressure across the interface by Laplace law:

p cap = σ nw,w κ (3.19)
Finally, capillary pressure can also be measured experimentally at a macroscopic scale4 (see chapter 13 for experimental capillary pressure measurement). An empiric expression obtained by van Genuchten [START_REF] Van Genuchten | A closed-form equation for predicting the hydraulic conductivity of unsaturated soils[END_REF] is often used to model capillary pressure as a function of water content. This relation follows the same hypothesis as the relative permeabilities relations (no chemical relations and no mechanical couplings):

p cap = p 0 (T ).([S * ] -1/m -1) 1-m (3.20)
CHAPTER 3. HYDRODYNAMIC BEHAVIOR OF THE AQUIFER with m the same exponent as in section 3.2.2 for relative permeability, p 0 is the pore entry pressure, i.e. the capillary pressure at which the non-wetting fluid can penetrate the considered porous medium for the initial drainage.

It is however important to recall that capillary pressure exhibits a strong hysteresis during imbibition (injection of wetting fluid) and drainage (injection of non-wetting fluid i.e. injection of carbon dioxide). Indeed, in the case of imbibition, the liquid fills spontaneously the pores, while for drainage, the non-wetting fluid is forced in the porous medium. Moreover, successive cycles of imbibition/drainage lead to different capillary curves as represented in figure 3.5. However, in the case of CO 2 injection, there is only one cycle of drainage (if we do not consider the resaturation following a stop of the injection), and we can then neglect the hysteresis of the capillary pressure curve. Figure 3.5: Evolution of capillary pressure with water saturation S w in drainage and imbibition conditions (from [START_REF] Md | Bubble free resine for infusion process[END_REF]).

Capillary Currents According to figure 3.5, we can see that as water saturation decreases, capillary pressure increases. Definition of capillary pressure given in equation (3.18), shows that an increase of capillary pressure leads to a decrease of water pressure (if carbon dioxide pressure is considered as constant) and thus a decrease of water saturation leads to a decrease in water pressure. This means that if the system considered presents a gradient in water saturation, there will be a gradient in water pressure which triggers movement of water within the system, bringing back water toward the drier zones [START_REF] Mahadevan | Flow-through drying of porous media[END_REF][START_REF] Mahadevan | Water removal from porous media by gas injection: experiments and simulation[END_REF][START_REF] Peysson | Diffusive and convective drying on a shaly sandstone using NMR measurements[END_REF]. This effect which plays a role both at sample and field scales, has 3.4. HYDRODYNAMIC REGIMES a strong influence on the displacement of salts [START_REF] Le | Impact of capillary-driven liquid films on salt crystallization[END_REF] during carbon dioxide injection (see section 6.3) Another manifestation of these capillary currents occurs when the injection of CO 2 is stopped for any reason (technical or other. . . ). As a result, there is no more carbon dioxide displacing brine or evaporating it. There is then a strong flow of brine which comes from the rest of the aquifer to resaturate the zones close to the injection well. This phenomenon may have a strong impact on the salt quantity and on the behavior of the aquifer and is very likely to happen as stops in the injection process are highly probable during the exploitation.

However, the existence of these capillary currents is conditioned by the mobility of the aqueous phase and thus by the current relative permeability of the considered medium. As a consequence, once the residual water saturation is reached, capillary currents will be strongly reduced and the considered medium can in first approximation be considered as isolated from the rest of the aquifer.

Hydrodynamic regimes

For most of the reservoir rocks, relative permeability of water reaches 0 while water saturation is still high (20 to 40 %). This is explained by the fact that water being the wetting fluid, carbon dioxide can only percolate and cannot displace completely the remaining water, trapped by capillary forces. The work of Lenormand [START_REF] Lenormand | Capillary fingering -percolation and fractal dimension[END_REF][START_REF] Lenormand | Numerical models and experiments in immiscible displacements in porous media[END_REF][START_REF] Lenormand | Liquids in porous media[END_REF] 5 highlights the different regimes and the behavior of fluids during multiphase transport. Indeed, depending on the capillary number (measuring the ratio between capillary forces and viscous forces), on the viscosity ratio, and the type of flow (imbibition or drainage), the hydrodynamical behavior of the fluids in the porous medium can vary strongly.

Classical capillary number for supercritical carbon dioxide injection is:

Ca = vη sc σ ≈ 1.5.10 -6 (3.21)
with v= 1 mm/s the flow velocity of the injected fluid (here CO 2 ), η sc = 4.5.10 -5 Pa.s its dynamic viscosity, and σ = 30 mJ.m -2 [START_REF] Chiquet | CO2/water interfacial tensions under pressure and temperature conditions of CO2 geological storage[END_REF] its surface tension at 80℃ and 20MPa. Capillary number is scaled on the supercritical carbon dioxide viscosity because it is the injected fluid.

Viscosity ratio is defined as: CHAPTER 3. HYDRODYNAMIC BEHAVIOR OF THE AQUIFER

M = η sc η l ≈ 0.126.10 -3 (3.22)
η sc is the viscosity of carbon dioxide in the supercritical state and η l = 0.3596 mPa.s is the dynamic viscosity of liquid pure water (at 20 MPa and 80 ℃). According to the sketch from [START_REF] Lenormand | Capillary fingering -percolation and fractal dimension[END_REF] presented in figure 3.6, we can see that the hydrodynamic regime during injection of carbon dioxide will be viscous fingering6 (see also [START_REF] Maloy | Viscous fingering fractals in porous media[END_REF]). The behavior is then controlled by viscous forces and carbon dioxide will follow the paths of less hydrodynamic resistance R h . For a cylindrical pore, we have using the Poiseuille law [START_REF] Sahimi | Flow and transport in porous media and fractured rock[END_REF]:

R h = 8η sc L πr 4 p in kg m 4 .s -1 (3.23)
with r p , the radius of the pore, L its length.

HYDRODYNAMIC REGIMES

We see easily that carbon dioxide will follow the biggest pores and create a pattern of fingers (cluster), separating blobs of water confined in the smaller pores as represented in figure 3.7. The consequence of this behaviour is that residual water will have three different forms: capillary trapped water in the corner of the biggest pore space, thin film (around 1nm [START_REF] Sgaier | On ions transport during drying in a porous medium[END_REF]) of water coating the surface of the pore walls (water is the wetting fluid), and blobs of water in the smallest pores acting as confined bulk water and separated by the percolating paths of carbon dioxide.

This residual water will be subjected to evaporation by the carbon dioxide. Evaporation will obviously happen at the interface between the supercritical phase and the brine. Consequently, the strongest evaporation will occur when the carbon dioxide flow rate is the highest and when the surface between the two phases is the largest. The largest surface of contact is the one exposed by the physisorbed wetting fluid, namely the thin film of water. However, because of interaction forces, evaporation in this film is negligible. There will then be two different kinds of drying in the rock: the drying of the capillary trapped water and the drying of the blob water. As we will see in part III, these two kinds of drying will lead to two different kinds of crystallization pressure generation.

Carbon dioxide/brine partitioning

Before carbon dioxide injection, the reservoir rock is in chemical equilibrium with the filling brine. When carbon dioxide is injected, the chemical equilibrium is disturbed and several reactions are triggered. Indeed, when carbon dioxide is in contact with brine, partitioning occurs, leading to a dissolution of CO 2 in the brine and the evaporation of water in the carbon dioxide phase (section 4.1.3). The dissolution of CO 2 in brine modifies the pH and the molalities of the different dissolved species displacing the different equilibriums. In the following chapter, we review the chemical behavior of the aquifer subjected to CO 2 injection. The first section describes the basis of geochemistry and the water/CO 2 partitioning, and the second section describes the chemical reactions between dissolved species and minerals induced by the dissolution of CO 2 in the brine.

Thermodynamics of mixtures

Gibbs free energy and chemical potentials

Gibbs free energy (or free enthalpy) is the most useful thermodynamic potential for chemistry and chemical equilibria1 .

It is defined as:

G = H -T S (4.1)
H is the enthalpy, T the temperature and S the entropy.

The partial molar quantity related to G is the chemical potential:

µ i = ∂G ∂n i p,T,n j,j =i (4.2)
with n j the molar quantities of the species present in the system considered, µ i is called the "chemical potential of the species i". As we will see later, chemical potentials determine the evolution of a system and its equilibrium condition.

The exact differential of G can be written

dG = V dp -S dT + i µ i dn i (4.3)
V is the volume of the system considered. Combination of (4.3) and the differential form of (4.4) gives the Gibbs-Duhem equation:

i n i dµ i = V dp -S dT (4.5)
Let us now consider a system with only one component i divided in two phases J and K. At equilibrium, T and p are constant and equation (4.3) writes:

µ J i dn J i = -µ K i dn K i (4.6)
Considering a closed system, we have dn i = 0 = dn J i +dn K i . As a result for an infinitesimal variation of n J i , at equilibrium:

dn J i µ J i -µ K i = 0 (4.7)
That is to say:

µ J i = µ K i (4.8)
The evolution of a chemical system is controlled by the evolution of the chemical potentials. Indeed, the condition of equilibrium is that the chemical potential of a species in the system is homogeneous whatever the phase. Any variation of chemical potential will result in a transformation leading to a homogeneization of the chemical potential of the considered species.

Chemical Potential Expression

Demonstration of the expressions of chemical potentials can be found for example in [START_REF] Coussy | Mechanics and Physics of Porous Solids[END_REF] or [START_REF] Gaskell | Introduction to the thermodynamic of materials[END_REF] 2 . Chemical potentials are commonly written as the sum of a reference chemical potential (corresponding to a reference state called the standard state) and a term of CHAPTER 4. CHEMICAL BEHAVIOR OF THE AQUIFER mixing:

µ i = µ 0 i + RT ln a i (4.9)
a i is the activity of the component i in the mixture.

• For a liquid mixture (or an aqueous solution), activity of a compound i is described by an activity coefficient γ i (see section 4.1.4.1 for activity coefficients model and experimental data) and a function of the mixture composition: mole fraction (

x i = n i /n tot ), molar concentration (c i = n i /V solution ) or molality (m i = n i /m solvant ).
For example:

a i = γ i m i (4.10)
Following the geochemistry convention as in [START_REF] Hc | Theoretical prediction of the thermodynamic behavior of aqueous-electrolytes at high-presures and temperatures .4. calculation of activity-coefficients, osmotic coefficients, and apparent molal and standard and relative partial molal properties to 600-degrees-c and 5kb[END_REF], the standard state of compound (i.e. a i = 1 and µ i = µ 0 i ) is defined as 1 mole of the considered aqueous species in a hypothetical infinitely diluted solution (ideal behavior) containing 1kg of water, for any temperature and pressure.

• For a gas, activity is described by the fugacity (see Appendix A.2 for models of fugacities and experimental data). Similarly to the activity of aqueous solution, fugacity is expressed as the product of the pressure and a fugacity coefficient. One has:

a i = f i f 0 i with f i = Φ i p i (4.11)
f i is the fugacity of component i in the gas mixture, Φ i its fugacity coefficient and p i its partial pressure; f 0 i = p 0 = 1bar is the reference fugacity of the gas in its standard state, that is when it is pure and perfect, whatever the temperature T of the system.

• Solids behave similarly to gas and liquids. However, in the case of precipitation/dissolution from solutions, there is no solid solution ("mixture" between two solids). As a result, all solids considered in the following will be taken as pure i.e.

a s = 1 and µ s = µ 0 s .

Mixing Quantities

Gibbs free energy of mixing The Gibbs free energy of mixing is defined as the Gibbs free energy between the mixture of components i and the components in their standard

∆G = G(p, T, n i ) -G 0 (p, T, n i ) (4.12a) ∆G = i n i µ i - i n i µ 0 i (4.12b) ∆G = i n i RT ln(a i ) (4.12c)
We can divide the free enthalpy of mixing into two terms. The first one, noted G comp , contains the composition of the system, x i while the second one, noted G ex contains the activity coefficients γ i .

∆G = RT i n i ln(x i ) G comp + RT i n i ln(γ i ) G ex (4.13)
We will see in next paragraph that G ex characterises the non-ideality of the mixture.

Entropy and enthalpy of mixing

The exact differential of the Gibbs free energy (4.3), allows us to express the entropy of mixing as:

∆S = - ∂∆G ∂T p,n i = - ∂G comp ∂T p,n i - ∂G ex ∂T p,n i (4.14a) ∆S = -R i n i ln(x i ) S comp -   R i n i ln(γ i ) + RT i n i ∂ ln γ i ∂T p,x i   S ex (4.14b)
The entropy of composition variation, S comp , is positive (the argument of the logarithm is smaller than one). Mixing creates disorder because of the augmentation of possible configurations induced by the addition of new molecules: each molecule can have neighbours of different kinds increasing the lack of information on the system. Excess entropy S ex corresponds to the fact that the interactions between the species in the mixture are not the same and thus influence the repartition of the molecules.

Finally, because ∆G = -T ∆S + ∆H we can write the enthalpy of mixing corresponding to the energy difference between compounds in the mixture and in their standard state: Ideal mixtures and consequences There are several equivalent definitions of ideal mixtures. The most intuitive one is the excess definition: an ideal mixture is a mixture in which all the activity coefficients are constant and equal to one (i.e. activities are equal to concentrations) and thus all the excess quantities are equal to 0. From equation (4.16), this condition leads to ∆H = 0 which means that the chemical compounds in the ideal mixture still behave as if they were in their standard state. Ideal mixtures are then composed of similar molecules (like benzene and toluene for example) 3 .

The second consequence of the ideal behavior of a mixture is that the entropy of mixing ∆S reduces to S comp = -R i n i ln(x i ), according to equation (4.14b).

Finally, in an ideal solution, molar volume of species in their standard state and partial molar volumes of these species in the mixture at the same T and p are equal (demonstration with Gibbs-Duhem relation).

However, the ideal behavior of a mixture is non-realistic and exists only for very dilute solutions, or mixtures of the same kind of molecules. In the case of CO 2 injection, the huge amount of salt in the brine creates a strong non-ideal behavior (cf 4.1.4.1).

Molar quantities of reaction

Using molar quantities is the most convenient and most usual way to describe the thermodynamic properties of chemically reactive systems. If we consider the chemical reaction ,(A,B . . . being chemical species: ions, gas . . . ), we define the molar quantity of reaction for an extensive thermodynamic property as:

A ν A B ν B . . . --⇀ ↽--ν A A + ν B B + . . .
∆ r X = i ν i ∂X ∂n i p,T,n k =i (4.17)
with i the species in the mixture and ν i the stoichiometric coefficient of species i. Usually, X is the volume, the free enthalpy, the entropy. . . The ν i are considered as algebraical 4.1. CARBON DIOXIDE/BRINE PARTITIONING coefficients: ν i < 0 for the reactants (left-hand side) and ν i > 0 for the products (righthand side) of the reaction. In the particular case of the Gibbs free energy we have:

∆ r G = i ν i ∂G ∂n i P,T,n k =i = i ν i µ i (4.18)
Using the expressions of chemical potentials given previously, we have:

∆ r G(p, T ) = i ν i µ 0 i (p, T ) + RT i ν i ln(a i ) (4.19a) ∆ r G(p, T ) = ∆ r G 0 (p, T ) + RT ln Q R (4.19b) Q R = i a ν i
i is called the reaction quotient (or ion activity product) and ∆ r G 0 (p, T ) is the standard Gibbs free energy of the reaction at pressure p and temperature T . ∆ r G is a decreasing function (in absolute value) of the reaction progress and equilibrium is reached when ∆ r G = 0. As the standard state of a gas is defined only at 1 bar (see section 4.1.1.3), if a gas is implied in the reaction considered, the standard Gibbs free energy of reaction is only defined at p = 1bar. The effect of a different pressure is then only grasped in the ion activity product.

It is then possible to define the property K s as follows:

∆ r G 0 (p, T ) = -RT ln K s (p, T ) (4.20)
Then equation (4.19b) reads:

∆ r G = RT ln(Q R /K s ) (4.21)
K s is the equilibrium constant of the reaction. It corresponds to the value taken by Q R at the equilibrium, when ∆ r G = 0. Data for K s for different chemical reactions are available in handbooks [START_REF]Handbook of Chemistry and Physics[END_REF] and thermodynamic databases (like [START_REF] Ph | THERMODDEM: a geochemical database focused on low temperature water/rock interactions and waste materials[END_REF]).

With this expression we find a demonstration of the Guldberg and Waage relation for equilibria of chemical reactions (commonly know as the Law of Mass Action), giving the activities of the concerned species at chemical equilibrium:

K s = i a ν i i eq (4.22)
The condition ∆ r G = 0 expresses a thermodynamical equilibrium situation. However, this situation may not be reached for all reactions in a given chemical system because of kinetics constraints. The system is then in a metastable state. The nucleation phenomenon (section 4.2.3.1) for example is one cause of a system to be in a metastable state. To measure the departure from equilibrium, we use the saturation ratio which is defined as:

SR = Q R K s = i a ν i i K s (p, T ) (4.23)
Its logarithmic version is called the saturation index (SI):

SI = log (SR) (4.24)
Saturation index and saturation ratio are key parameters for geochemistry because they inform on the direction of evolution of a chemical reaction. SR < 1 (or SI < 0) means that the mineral, if present in the system, has the tendency to dissolve. SR > 1 (or SI > 0) means the mineral is thermodynamically able to precipitate. In this case, the saturation ratio is also called supersaturation and is written S. Supersaturation is particularly useful for the study of crystallization and crystallization pressure (see section 4.2.3 and chapter 8).

Laws for CO 2 /H 2 O partitioning

Numerous authors have published experimental data and models of solubility of carbon dioxide in saline water (see for example [START_REF] Driesner | The system H 2 O-N aCl. part I: Correlation formulae for phase relations in temperature-pressure-composition space from 0 to 1000C, 0 to 5000 bar, and 0 to 1X N aCl[END_REF] and [START_REF] Driesner | The system H 2 O-N aCl. part II: Correlations for molar volume, enthalpy and isobaric heat capacity from 0 to 1000C, 0 to 5000 bar, and 0 to 1 X N aCl[END_REF] for calculation of thermodynamical properties of the system H 2 O, NaCl, CO 2 , see also [START_REF] Lassin | Geochemistry of unsaturated soil systems: Aqueous speciation and solubility of minerals and gases in capillary solutions[END_REF] for the description of Thermo-ZNS calculation for geochemistry). They are all based on the same ground i.e. thermodynamic equilibria, but they become more and more complex and more and more accurate, incorporating Pitzer's equation [START_REF] Pitzer | Thermodynamics of electrolytes. 1. theory basis and general equations[END_REF][START_REF] Pitzer | Thermodynamics of electrolytes. 2. activity and osmotic coefficients for strong electrolytes with one or both ions univalent[END_REF][START_REF] Steiger | An improved model calculating Pitzer's equations for calculation of thermodynamic properties of pore solutions implemented into an efficient program code[END_REF] for activity corrections in saline solutions; corrections for the supercritical state . . . In this section we will consider Henry's law which is the most used relation for carbon dioxide solubility, and the Spycher et al. model [START_REF] Spycher | CO2-H2O mixtures in the geological sequestration of CO2. I. assessment and calculation of mutual solubilities from 12 to 100C and up to 600 bar[END_REF][START_REF] Spycher | CO2-H2O mixtures in the geological sequestration of CO2. II. partitioning in chloride brines at 12-100C and up to 600 bar[END_REF] because it is implemented in most of the releases of a numerical code called TOUGH which is extensively used in CCS simulations(see section 6.1).

Henry's law

Equilibrium of partitioning is calculated differently depending on the author. Extended Henry's law (here the formulation from [START_REF] Diamond | Solubility of CO2 in water from -1.5 to 100C and from 0.1 to 100MPa: evaluation of literature data and thermodynamic modelling[END_REF]) is the most used formula (for a thermodynamic proof of this law see for example [START_REF] Coussy | Mechanics and Physics of Porous Solids[END_REF]): The same reasoning can be applied to the equilibrium between liquid water and vapor in the CO 2 rich phase. As temperature is quite low and pressure quite high, the quantity of water vapor in the CO 2 phase is quite small and is reached very quickly (see section 4.1.5).

Calculation with equilibrium constants

Another method to calculate CO 2 /brine partitioning is to follow Spycher's procedure, using directly the free enthalpy of separation [START_REF] Spycher | CO2-H2O mixtures in the geological sequestration of CO2. I. assessment and calculation of mutual solubilities from 12 to 100C and up to 600 bar[END_REF][START_REF] Spycher | CO2-H2O mixtures in the geological sequestration of CO2. II. partitioning in chloride brines at 12-100C and up to 600 bar[END_REF][START_REF] Spycher | A phase partitioning model for CO2-brine mixtures at elevated temperatures and pressures: Application to CO2-enhanced geothermal systems[END_REF]. They start from the liquid/gas equilibria for carbon dioxide and water:

H 2 O (l) --⇀ ↽--H 2 O (g) K H 2 O = f H 2 O /a H 2 O (4.26) CO 2(l) --⇀ ↽--CO 2(g) K CO 2 = f CO 2 /a CO 2 (4.27)
K are the thermodynamic equilibrium constants (tabulated for example in [START_REF] Ph | THERMODDEM: a geochemical database focused on low temperature water/rock interactions and waste materials[END_REF] and not to be mistaken with the Henry's constant which is an empirical parameter), f the fugacities in the vapor and a the activities in solution. Using the pressure variation of K and the definition of fugacity (f = Φp tot y), one can write the mole fractions of the components (y is the mole fraction in the gas phase, and x the one in the aqueous phase):

y H 2 O = K 0 H 2 O (1 -x CO 2 -x salt ) Φ H 2 O p tot exp (p -p 0 )ν H 2 O RT (4.28)
In first approximation, we consider the solvent as ideal (see 4.1.1.4): activity of water is taken equal to its mole fraction:

x CO 2 = Φ CO 2 (1 -y CO 2 )p tot 55.508γ ′ x K 0 CO 2 exp - (p -p 0 )ν CO 2 RT (4.29)
CHAPTER 4. CHEMICAL BEHAVIOR OF THE AQUIFER considering that a CO 2 = 55.508x CO 2 (computed from the molalities); γ ′ x is the activity coefficient of aqueous CO 2 on a mole fraction scale and Φ CO 2 the fugacity coefficient of the supercritical carbon dioxide.

Setting A = K 0 H 2 O Φ H 2 O ptot exp (p-p 0 )ν H 2 O RT and B = Φ CO 2 ptot 55.508γ ′ x K 0 CO 2 exp - (p-p 0 )ν CO 2

RT

, we can rewrite these equations with the molality of the dissolved salt: [START_REF] Ph | THERMODDEM: a geochemical database focused on low temperature water/rock interactions and waste materials[END_REF].508 (1/A -B)(νm salt + 55.508) + νm salt B) (4.30)

y H 2 O = (1 -B)
x CO 2 = B(1 -y CO 2 ) (4.31)
ν is the stoichiometric number of ions in the dissolved salt (e.g. ν = 2 for NaCl, and ν = 3 for CaCl 2 ) The kinetics of this partitioning is really fast and can be considered in first approximation as instantaneous [START_REF] Spycher | CO2-H2O mixtures in the geological sequestration of CO2. I. assessment and calculation of mutual solubilities from 12 to 100C and up to 600 bar[END_REF][START_REF] Penner | On the kinetics of evaporation[END_REF]. The two approaches are similar and both give accurate results.

Activity and fugacity coefficients

In an ideal mixture (see section 4.1.1.4), activities of the different species can be approximated by the molar fraction (or molality, molar concentration . . . ). Real mixtures are not as simple and activities can depart significantly from these quantities. Activity and fugacity coefficients quantify this gap between real behavior and ideality. They can be obtained through several semi-empirical relations.

Activity coefficients of dissolved species: extended-Debye-Hückel law

The most simple law for activity coefficients is the Debye-Hückel law and its extensions.

In particular, in Thermo-ZNS code (see section 9.1.1) is used an extension of the Debye-Hückel theory called B-dot [START_REF] Hc | Theoretical prediction of the thermodynamic behavior of aqueous-electrolytes at high-presures and temperatures .4. calculation of activity-coefficients, osmotic coefficients, and apparent molal and standard and relative partial molal properties to 600-degrees-c and 5kb[END_REF]. The code calculates (among other things) the B-dot parameters which can then be used in chemical codes like PhreeqC. In this model, the activity coefficient of species i is calculated from the ionic strength

I = 1 2 i z 2 i m i as: log(γ i ) = - A γ z 2 i √ I 1 + • a i B γ √ I + ḂI (4.32)
i refers to the considered ion, z i is its electric charge,

• a i represents its mean effective electrostatic radius. Ḃ is the characteristic B-dot parameter. A γ and B γ are the temperature 4.1. CARBON DIOXIDE/BRINE PARTITIONING and pressure dependent Debye-Hückel parameters for aqueous solutions:

A γ = √ 2πN a ρ ln 10 √ 1000 e 2 4πε 0 εk B T b (4.33a) B γ = 2e 2 N a ρ 1000ε 0 εk B T (4.33b)
where ρ is the density of pure water, ε 0 is permittivity of vacuum, ε is the relative permittivity of water, e is the elementary charge and k B is the Boltzmann constant. At 200 bar and 75℃ A γ = 0.5568 kg 1/2 mol -1/2 and B γ = 0.3367.10 10 kg 1/2 mol -1/2 m -1 . At 45 ℃ we have for NaCl, Ḃ = 0.043.

Unfortunately, this model is not valid for high ionic strength (I > 1) in aqueous solution where NaCl is not the major dissolved salt. In order to obtain more accurate evaluation of the activity coefficients, Pitzer's model has to be used (see appendix A).

Activity coefficient of water

Relation (4.32) gives the activity coefficient for solutes. Activity of water which is important in the calculation of equilibrium for hydrated minerals (like gypsum, see section 4.2), has to be calculated differently. Indeed, it cannot be calculated on a molality basis (which would give a constant activity whatever the composition) and has thus to be calculated on a molar fraction basis: a w = λ w x w (λ w is the activity coefficient for water on a molar fraction basis). For non ideal solutions, water activity is obtained through the following relation:

ln a w = - Σm m w + 1 RT ∂G ex ∂n w (4.34)
with Σm = i m i , the sum of molalities of all solutes (without water), G ex the excess Gibbs free energy (as defined in equation (4.13)), m w , the number of moles of water in 1kg of pure water ( ≈ 55.508), and n w the number of moles of water in the considered solution. The expression of water activity is then dependent on the expression of the excess Gibbs free energy. A precise calculation of the water activity with Pitzer model is made in appendix A. In the B-dot model, water activity has a much simpler expression:

ln a w = 1 m w -Σm + 1 3 A γ I 3/2

Fugacity coefficients

Calculation of fugacity coefficients is based on the same ideas as activity coefficients but is strongly dependent on the equation of state used in the modelling of the phase. Calculation of fugacity coefficients with a virial development and a Redlich-Kwong equation are presented in appendix A.2.

Kinetics of evaporation and simplifications

Kinetics of evaporation can be expressed with the Hertz-Knudsen relation [START_REF] Penner | On the kinetics of evaporation[END_REF] which is calculated with the Kinetic Theory of Gases (KTG).

v evap = - p w -p sat w 2πM H 2 O RT (4.36)
v evap is the speed of evaporation (mol/m 2 /s). This speed characterises the quantity of evaporated water per unit of contact area between the gas and the liquid and per unit of time. p w is the current vapor pressure of water, p sat w is the saturation vapour pressure, M H 2 O is the molar mass of water, R the perfect gases constant and T the current temperature. However, this equation is only valid for evaporation of water with water vapor as a perfect gas and does not describe the evaporation of water in a supercritical phase. Nevertheless, we can consider that the physics behind Hertz-Knudsen relation can be applied to supercritical carbon dioxide without strong modifications. Numerical estimations of the evaporation rate in the case of CCS show that the kinetics of evaporation of water and thus the saturation of the CO 2 -rich phase can be considered as instantaneous. As a result, drying of the aquifer will only be controlled by the carbon dioxide flow rate.

CHEMICAL REACTIONS INDUCED BY CARBON DIOXIDE DISSOLUTION (MINERAL TRAPPING)

The TOUGH releases (see section 6.1) for example consider an instantaneous equilibrium between liquid water and vapor.

Chemical reactions induced by carbon dioxide dissolution (mineral trapping)

Chemical reactions are likely to happen during CO 2 injection4 . Indeed, injection of carbon dioxide in the aquifer disrupts dramatically the chemical equilibria leading to several reactions as precipitation, dissolution, and mineral transformation. In the following we will only consider aqueous chemical reactions and brine/rock reactions. Gas/rock chemical reactions will not be considered.

Carbon dioxide dissolution and pH evolution

The dissolution of CO 2 in water disrupts the initial chemical equilibrium of the aqueous phase. Indeed, CO 2 in water decreases pH through the reactions [4,[START_REF] Johnson | Reactive transport modelling of CO2 storage in saline aquifers to elucidate fundamental processes, trapping mechanisms, and sequestration partitioning[END_REF][START_REF] André | Numerical modeling of fluid rock chemical interactions at the supercritical CO2-liquid interface during CO2 injection into a carbonate reservoir, the Dogger aquifer (Paris Basin, France)[END_REF][START_REF] Izgec | CO2 injection into saline carbonate aquifer formations II: Comparison of numerical simulations to experiments[END_REF]:

CO 2 (g) + H 2 O(aq) --⇀ ↽--CO 2 (aq) + H 2 O --⇀ ↽--HCO 3 -+ H + --⇀ ↽--CO 3 2-+ 2 H + (4.37)
The acidity constants for the two reactions are (under ambient conditions of temperature and pressure: 25℃ and 1 bar):

pKa(CO 2 /HCO 3 -) = 6.37 (4.38a) pKa(HCO 3 -/CO 3 2-) = 10.3 (4.38b)
In figure 4.1 is plotted the evolution of the proportion of the different carbon dioxide aqueous species with pH. At low pH, CO 2 (aq) is predominant and the other species are negligible. At moderate pH, between pKa(CO 2 /HCO 3 -) and pKa(HCO 3 -/CO 3 2-), the hydrogenocarbonate ion HCO 3 -is predominant, while at high pH, it is the carbonate ion CO 3 2-which is predominant.

Mineral reactions and carbon dioxide trapping

Mineral trapping (CO 2 reacting with the rock and becoming part of the mineral) is a really long process [START_REF] Bachu | Sequestration of CO2 in geological media in response to climate change: Capacity of deep saline aquifers to sequester CO2 in solution[END_REF][START_REF] Johnson | Reactive transport modelling of CO2 storage in saline aquifers to elucidate fundamental processes, trapping mechanisms, and sequestration partitioning[END_REF][START_REF] Lagneau | Reactive transport modelling of CO2 sequestration in deep saline aquifers[END_REF], but is the final aim of CO 2 storage. Reactions are numerous and depend of the initial composition of the reservoir rock and of the brine (see for example the fluid composition of the Montmiral site in [START_REF] Pauwels | Chemistry of fluids from a natural analogue for a geological CO2 storage site (Montmiral, France): Lessons for CO2-water-rock interaction assessment and monitoring[END_REF]). They include K-feldspar transformation into dawsonite:

KAlSi 3 O 8 K-feldspar + Na + + CO 2 (aq) --N aAlCO 3 (OH) 2 dawsonite + 3SiO2 qtz/chal/crist + K + (4.39) (qtz = quartz, chal = chalcedony, crist = cristobalite)
Calcite-group precipitation:

M 2+ + CO 2 (aq) + H 2 O --MCO 3 + 2 H + (4.40)
Clays can also react with dissolved CO 2 as magnesite cementation (Fe-Mg-rich clays One must remember that these are only examples. There exist numerous reactions which can obviously not be written here. In a complex mineralogical system, reactions can be ordered through the thermodynamic constants by identifying the most favorable reactions. Numerical methods using software like PhreeqC allow to calculate easily the complete behavior of a complex system. Solubility of minerals depends on the pH and on the amphoteric properties of the dissolved species. Equation of reaction for a considered mineral is then written according to the existing species at the considered pH. Below are written the equations for dissolution with an excess of protons. The last reaction is independent of the pH. Dissolution of anhydrite and gypsum are written for pH > 2 (this condition is never always fulfilled the classic conditions of CCS). Obviously the presence of these equilibria is subjected to the initial conditions i.e. type of reservoir rock and brine composition.

Calcite + H + = HCO 3 -+ Ca 2+ Dolomite + 2 H + = 2 HCO 3 -+ Ca 2+ + Mg 2+ Magnesite + H + = Mg 2+ + HCO 3 - Siderite + H + = HCO 3 -+ Fe 2+ Dawsonite + 3 H + = HCO 3 -+ Na + + Al 3+ + 2 H 2 O K-feldspar + 4 H + = K + + Al 3+ + 3 SiO 2 + 2 H 2 O Kaolinite + 6 H + = 2 Al +3 + 2 SiO 2 + 5 H 2 O Albite + 4 H + = Na + + Al 3+ + 3 SiO 2 + 2 H 2 O Anhydrite = SO 4 2-+ Ca 2+ Gypsum = SO 4 2-+ Ca 2+ + 2H 2 O Halite = Na + + Cl -
All these equilibria are of concern in the field (depending on the initial conditions) but their roles depend on the local environment and are strongly dependent on the comparison between chemical kinetics, flow rates and diffusion processes. Moreover, one mineral can dissolve at one place and re-precipitate in another place. Initial composition CHAPTER 4. CHEMICAL BEHAVIOR OF THE AQUIFER of the brine is then really important as are temperature and flow rate. The last three minerals which are the most susceptible to precipitate will be studied more precisely in the modelling (chapter 9).

Temperature acts in three ways: a reduction of the temperature increases CO 2 dissolution, modifies minerals solubility and acts on the kinetics of reaction through the Arrhenius law (4.55). The flow rate acts on the reaction in parallel to the kinetics and defines the zones of dissolution and precipitation.

Kinetics of mineral reactions: nucleation and crystal growth 4.2.3.1 Nucleation

Nucleation is the process by which a new solid phase is created, from a solution or a melt. Homogeneous nucleation, which is the nucleation of a new solid phase from the bulk is often referred to as a "difficult process". Indeed, the creation of a solid phase within a bulk solution is statistically and thermodynamically unfavorable until strong supersaturations (for example, for halite, homogeneous nucleation occurs at supersaturations as high as 1.8 -1.9). This homogeneous nucleation is then a rare phenomenon, and in natural media, the presence of impurities or other crystalline phases decrease the energy barrier for nucleation and thus leads to nucleation at a lower supersaturation (this is called heterogeneous nucleation and is described at the end of this section). The reason of the difficulty of nucleation arises from the fact that, in a comparison to a big crystal, for a small crystal, the surface energy contribution is much higher relative to the volumetric energy contribution. Indeed, let us consider the Gibbs free energy difference between ions in solution and the same ions in a small crystal of radius r in a supersaturated solution:

∆G(r) = - 4πr 3 3υ s ∆ r G + 4πr 2 σ (4.42)
with υ s the molar volume of the crystal, σ its surface energy.

From equations (4.21) and (4.24), we substitute ∆ r G by the supersaturation S:

∆G(r) = - 4πr 3 3υ s RT ln S + 4πr 2 σ (4.43)
If the Gibbs energy difference ∆G(r) is positive, the nucleation is thermodynamically unfavorable and will not happen. Let us now consider the variation of ∆G(r), with the radius of the crystal:

CHEMICAL REACTIONS INDUCED BY CARBON DIOXIDE DISSOLUTION (MINERAL TRAPPING)

∂∆G ∂r = - 4π υ s r 2 RT ln S + 8πrσ (4.44)
The variation of energy has then a maximum for a critical value of r: We can then rewrite the energy cost in terms of critical radius:

r c = 2συ s RT ln S (4.45)
∆G(r) = 2πr 3 c 3υ s RT ln S ∆Gc 3 r r c 2 -2 r r c 3 (4.46)
Figure 4.2 shows the evolution of this energy with the critical radius. We clearly see the metastable characteristic of the nucleation considering ∆G c as the energy barrier to the nucleation. ∆G c corresponds to the value of ∆G when r = r c . Nucleation is a statistical process: nuclei of different sizes appear in the bulk because of density fluctuations [START_REF] Christian | The theory of transformations in metals and alloys[END_REF][START_REF] Markov | Crystal Growth for beginners: fundamental of nucleation, crystal growth, and epitaxy[END_REF]. However, the stability of these nuclei depends on the energy cost. If the radius of the nucleus is smaller than the critical radius, the nucleus is not stable because any growth following the nucleation of this nucleus will cost some energy which is unfavorable. Only nuclei whose radii are bigger than the critical radius are stable, because their following growth is accompanied with a decrease in the global energy. Considering the critical radius we can see that the more the solution is supersaturated the smaller the critical radius will be and thus the easier the nucleation will be. Similarly, the smaller the interfacial energy is, the smaller the critical nuclei will be because the expense of energy to create the surface will be lower.

Nucleation of crystals is controlled by the nucleation rate. Theory of nucleation kinetics involves several models, considering the statistical fluctuations of density, the energy barrier of attachment of atoms on the formed embryo. . . Nucleation rate can be expressed as follows:

I = K 0 (T ) exp - ∆G c + ∆G m k B T (4.47)
with k B the Boltzmann constant, K 0 (T ) is the prefactor taking into account the mobility of the ions in the solution and their rate of attachment to the embryo and thus make it grow above the critical radius. The exponential term is an Arrhenius-like term considering the activation energy for the nucleation as ∆G c + ∆G m (with ∆G c defined in equation (4.46). The term ∆G m corresponds to the energy barrier for an atom to join the embryo.

Accordingly we can write using the expression of ∆G c obtained previously:

I = K 0 (T ) exp -∆G m k B T exp - 16πσ 3 υ 2 s 3k B T (RT ln S) 2 (4.48a) I = K 0 (T ) exp -∆G m k B T exp - 16πσ 3 α 2 3 (k B T ) 3 (ln S) 2 (4.48b)
with α the volume of a formula unit (α = υ s /N a , N a is the Avogadro number). This relation is valid only for spherical crystals. In order to be more general, we can also express the energy barrier with a shape factor λ such as:

λ = A n 2/3 (4.49)
with A the surface of the nucleus, and n the number of moles in it. Energy cost (4.43) and critical parameter (4.45) can then be rewritten: with n c the number of moles in the critical embryo.

∆G(n) = -nRT ln S + λn 2/3 σ = ∆G c 3 n n c 2/3 -2 n n c (4.50a) n c = 8 27
σ 3 λ 3 (RT ln S)
The energy barrier becomes: ∆G c =

4λ 3 σ 3
27(RT ln S) 2 which allows to rewrite the nucleation rate as:

I = K 0 (T ) exp -∆G m k B T exp   - 4 λ 3 υ 2 s σ 3 α 2 27 (k B T ) 3 (ln S) 2   (4.51) remarking that λ 3 υ 2 s = A 3 V 2 ,
with A the area of the nucleus and V its volume. In table 4.2 are listed the shape factors for several classical volumes.

The prefactor K 0 (T ) is difficult to estimate. Discussion on its expression can be found in [START_REF] Christian | The theory of transformations in metals and alloys[END_REF] or [START_REF] Markov | Crystal Growth for beginners: fundamental of nucleation, crystal growth, and epitaxy[END_REF]. It is a function of the mean free path of ions in the solution, of the concentration, of the critical radius . . . The calculation resulting in equation (4.51) has been made with an homogeneous nucleation situation. As said in the introduction of this section, homogeneous nucleation is a difficult process and does not happen often in reality. Indeed, the presence of substrates or impurities acts as help for the nucleation by reducing the cost of creation of the interfacial energy. The reduction of the cost can be calculated by considering an angle of contact θ between the crystal nucleated and the substrate5 : [START_REF] Espinosa-Marzal | Advances in understading damage by salt crystallization[END_REF] we can find a practical expression for the nucleation rate of crystallization as: ). Nucleation increases rapidly with supersaturation but this increase slows down at high supersaturations because of the evolution of the solution viscosity η acting as an obstacle: the increase of viscosity following the increase of concentration limits the diffusion of ions and thus the possibility to nucleate embryos. As we can see from this figure, the increase of nucleation rate is very steep and can reach really high values. This allows to define the nucleation threshold: if the supersaturation is below this threshold, the nucleation rate is negligible. However, if the supersaturation is above this threshold, nucleation is strong. The threshold is maximum for homogeneous nucleation and decreases with decreasing contact angle for heterogeneous nucleation. Nucleation threshold is clearly seen in the numerical simulations in chapter 9.

∆G hetero c = ∆G homo × f (θ) (4.52) with f (θ) = 1 2 -3 4 cos θ + 1 4 cos 3 θ. Finally, in
I = k B T πα 5/3 η exp - 256σ 3 α 2 27(k B T ) 3 (ln S) 2 ) f (θ) (4.53)

Crystal growth

The aim of chemical kinetics is the knowledge of the rate of reaction for the different species involved. If we consider that there is no physical limitation (diffusion. . . ), and that the growth is not constrained (e.g. unconfined crystallization), the rate of chemical reaction between dissolved species and minerals is obtained from the comparison of the rates of attachment and detachment of ions from the crystal phase and can be described as [START_REF] Lasaga | Chemical kinetics of water-rocks interactions[END_REF][START_REF] Espinosa | Phase changes of salts in porous materials: crystallization, hydration and deliquescence[END_REF][START_REF] Derluyn | Deformation and damage due to dryinginduced salt crystallization in porous limestone[END_REF]

r cryst = k cryst A S θ -1 η if S θ > S start ≥ 1 (4.54a) r sol = k sol A 1 -S θ η ′ if S θ < 1 (4.54b)
r cryst and r sol are respectively the reaction rates (mol/s) for precipitation (and dissolution); according to relations (4.54), r cryst and r sol are positive. We then have r cryst = dn/dt and r sol = -dn/dt with n the mole quantity of the crystal; k cryst and k sol are the kinetic constants (mol/m 2 /s) for precipitation and dissolution (usually we consider that k cryst = k sol ), A is the reactive surface area of the mineral (m 2 ), S is the supersaturation (or saturation ratio SR) of the mineral as defined in equation (4.23). S start represents the nucleation threshold. If there are already crystals in the solution, S start = 1. If not, the nucleation metastability (see section 4.2.3.1), forces the supersaturation to be higher than one for crystals to appear (S start > 1). Finally, θ, η, η ′ are empirical parameters which allow improving the description of the kinetics behavior of the reaction close to the equilibrium. Their values depend on the crystal considered and also on the environment (porous medium or bulk. . . ). For example, for mirabilite (Na

2 SO 4 • 10 H 2 O), η = η ′ = 1.9;
while for potassium nitrate (KNO 3 ), η = η ′ = 1; θ is almost always equal to unity [START_REF] Espinosa | Phase changes of salts in porous materials: crystallization, hydration and deliquescence[END_REF][START_REF] Espinosa | Model for the mechanical stress due to salt crystallization in porous materials[END_REF]. Equation (4.54) is however valid only if crystals can be at equilibrium with supersaturation equal to unity. In the case of confined crystals, equilibrium supersaturation can be higher than unity (see section 8.1.1): a modification of this kinetics relation for small crystallites is presented in section 9.2.1.1.

The kinetics constant k for a single mechanism, is a function of the chemical mechanism considered and temperature. It obeys the Arrhenius law:

k = k 0 exp -E a R 1 T - 1 T 0 (4.55)
with k 0 the intrinsic constant at the reference temperature, E a the activation energy (J/mol).

As said previously (section 4.2.2), chemical speciation controls the main active dissolution/reaction mechanism of any mineral. This also affect reaction kinetics and one has to distinguish the corresponding mechanisms in the rate laws. In general, acid, neutral and alkaline (or basic) mechanisms are considered for alumino-silicate minerals and acid, neutral and carbonate mechanisms can be considered for carbonate minerals:

• For alumino-silicates:

k = k nu 0 exp -E nu a R 1 T - 1 T 0 + k H 0 exp -E H a R 1 T - 1 T 0 a n H H + k OH 0 exp -E OH a R 1 T - 1 T 0 a n OH OH (4.56)
• For carbonates

k = k nu 0 exp -E nu a R 1 T - 1 T 0 + k H 0 exp -E H a R 1 T - 1 T 0 a n H H + k CO 2 0 exp   -E CO 2 a R 1 T - 1 T 0 a n CO 2 CO 2   (4.57)
Superscripts nu, H, OH, and CO 2 correspond respectively to neutral, acid, alkaline and carbonate mechanisms, a represents the activity of the concerned species. Values of the constants can be found in [4,[START_REF] Palandri | A compilation of rate parameters of water-mineral interaction kinetics for application to geochemical modeling[END_REF]. However, these values of kinetics constants are obtained experimentally in a batch and their extension to a porous model is questionable.

Porosity and permeability variations induced from dissolution/precipitation

Precipitation/dissolution of the minerals has a strong impact on transport through variation of (local) porosity and thus intrinsic permeability. Indeed, precipitation can lead to clogging of the material while dissolution opens new paths for the fluids. A simple model allowing to calculate intrinsic permeability reduction linked to porosity reduction is the tubes-in-series model [START_REF] Verma | Thermohydrological conditions and silica redistribution near high-level nuclear wastes emplaced in saturated geological formations[END_REF] as represented in figure 4.4:

K K 0 = θ 2 1 -Γ + Γ/ω 2 1 -Γ + Γ [θ/(θ + ω -1)] 2 with θ = 1 -S s -φ r 1 -φ r and ω = 1 + 1/Γ 1/φ r -1 (4.58)
K is the current permeability, K 0 is the initial permeability, Γ is the fractional length of pore bodies (see figure 4.4) and φ r a parameter denoting the fraction of initial porosity at which permeability is reduced to 0, S s is the volume fraction of the solid. Usually Γ and φ r are taken equal to 0.8. All the phenomena described in the previous chapters have an influence on the thermal behavior of the aquifer. Among the effects of temperature we can list the modification of thermodynamic constants of chemical reactions, of viscosity, of density, phase changes. . . In the following chapter, we review the different thermal effects and their relative importance in order to determine the evolution of the temperature in the aquifer during the injection. At first we consider thermo-hydrodynamic effects (Joule-Thomson expansion), then we consider latent heat of phase transition and heat of reaction, and finally we study the diffusion of heat through the Fourier law.

CHEMICAL REACTIONS INDUCED BY CARBON DIOXIDE DISSOLUTION (MINERAL TRAPPING)

Joule-Thomson expansion

The Joule-Thomson (JT) effect characterises the change of temperature of a fluid submitted to an adiabatic expansion (change of pressure without any exchange of heat i.e. at constant enthalpy).

The Joule-Thomson coefficient µ JT can be defined as [START_REF] Oldenburg | Joule-Thomson cooling due to CO2 injection into natural gas reservoirs[END_REF]:

µ JT = ∂T ∂P H ≈ ∆T ∆P (5.1)
Its values are obtained experimentally. Because of the small values of the JT coefficient (between 0 and 10 K/MPa), we can assume that this effect will be rather negligible [START_REF] Oldenburg | Joule-Thomson cooling due to CO2 injection into natural gas reservoirs[END_REF]. The principal consequence would be a reduction of a few ℃ close to the injection well.

Heat of reaction

For a reaction, the first law of thermodynamics gives the variation of the enthalpy of the considered system as:

dH = δQ p + V dp (5.2)
with δQ p , the heat exchange during the transformation. Considering the definition of heat, we can develop it as:

δQ p = C p dT + hdp + Q c dξ (5.3)

THERMAL DIFFUSION

The first term in the right-hand side of equation ( 5.3) corresponds to the modification of heat due to a change of temperature at constant pressure (with C p = i n i c i p the molar calorific capacity of the system, c i p being the molar calorific capacity of the species i and n i its mole quantity), the term second corresponds to the change of heat following a change of pressure at constant temperature, and the last to the heat algebraically released by the considered chemical reaction at constant pressure and temperature, ξ is the extent of reaction defined as dn i = ν i dξ.

As a conclusion, we can write that:

dH = C p dT + (h + V )dp + Q c dξ (5.4)
and finally:

Q c = ∂H ∂ξ T,p = i ν i ∂H ∂n i p,T,n k =i = ∆ r H (5.5)
As expected the heat exchanged during the isothermal and isobaric transformation is equal to the enthalpy of reaction at constant pressure and temperature. As a conclusion, the total heat variation following a chemical reaction at constant T and p is finally:

∆H = 2 1 ∆ r Hdξ = ∆ r H (ξ 2 -ξ 1 ) (5.6)
with ξ 2 the final extent of the considered reaction. In the case of a phase change, the enthalpy of reaction reduces to the difference of chemical potential in the phases 1 and 2 and is called latent heat. Values of latent heat and enthalpy of reaction are tabulated for numerous chemical reactions [START_REF]Handbook of Chemistry and Physics[END_REF].

Thermal diffusion

Like diffusion of species in solution, heat can diffuse through a medium following a gradient of temperature. The heat flux q (in W.m -2 ) through a surface is obtained by Fourier's law [START_REF] Coussy | Poromechanics[END_REF]:

q = -λ ∇T (5.7)
With λ the heat diffusion coefficient, depending on the thermal properties of the rock (a combination of the thermal properties that are constitutive of the rock), on the pore water content, temperature and pressure.

Importance of the different phenomena

Finally, in order to obtain the temperature variation, we need to consider all the phenomena together. Let us consider a REV in the aquifer. The evolution of enthalpy between a state 1 and a state 2 is:

∆H = 2 1 k n k c k p dT + R ∆ r H R ξ R 2 -ξ R 1 + 2 1 (h + V )dp = t 2 t 1 A q(M )dA(M )dt (5.8)
with h the enthalpy of the system (as a sum of the enthalpy of the different constituents), R a reaction (chemical, partionning...), M a point belonging to the surface A of the REV and t 1 and t 2 the times of the states 1 and 2. Knowing the characteristics of the considered REV, the different chemical reactions, and the calorific capacities, we can then determine the evolution of the temperature of the REV.

In [4], a simulation calculating the evolution of the temperature in the aquifer during the injection, is proceeded with TOUGH (see section 6.1). It appears that the latent heats and Joule-Thomson effect have a small impact, modifying the temperature by a few ℃ while the heat exchange with the rest of the aquifer and the other geological layers has much more importance. Moreover, close to the injection well, the temperature is set by the temperature of the injected CO 2 as the heat exchange from the rest of the aquifer is not sufficient to compensate the thermal disturbance induced by the injected carbon dioxide. As a result, in our calculations, as we only study the zone close to the well, we will consider a constant temperature which is set by the temperature of the injected carbon dioxide.

Previous paragraphs gave us the physical background to understand the hydrodynamic and chemical behavior of an aquifer subjected to carbon dioxide injection. Analytical solutions in the literature ( [START_REF] Burton | CO2 injectivity into brine aquifers: why relative permeability matters as much as absolute permeability[END_REF][START_REF] Zeidouni | Analytical solution to evaluate salt precipitation during CO2 injection in saline aquifers[END_REF] for example, solving the transport through Buckley-Leverett resolution [START_REF] Buckley | Mechanism of fluid displacement in sands[END_REF]; or [START_REF] Pruess | Formation dry-out from CO2 injection into saline aquifers. part 2, analytical model for salt precipitation[END_REF] for a resolution by mass balance) cannot take into account all the chemical and hydrodynamic processes and thus make strong simplifications. Then, numerical simulations are necessary in order to reproduce accurately the behavior of the aquifer.

TOUGH, an intensively used family of codes

The development of simulations of carbon injection into saline aquifers is closely related to the development of the TOUGH ("Transport Of Unsaturated Groundwater and Heat") code family. This code developed by Karsten Pruess [START_REF] Pruess | The TOUGH codes: A family of simulation tools for multiphase flow and transport processes in permeable media[END_REF] at the Lawrence Berkeley National Laboratory in the 1990s is an evolution of MULKOM developed in the early 1980s. MULKOM was dedicated to solve multicomponent multiphase fluid flows. TOUGH became more and more efficient including several modules and packages to become TOUGHREACT [START_REF] Xu | Toughreact: A simulation program for non-isothermal multiphase reactive geochemical transport in variably saturated geologic media: application to geothermal injectivity and CO2 geological sequestration[END_REF]. The first idea was to simulate movement of groundwater for geothermal and nuclear waste applications but the CO 2 storage evolution led to a module named ECO2N [START_REF] Pruess | ECO2N -a fluid property module for the TOUGH2 code for studies of CO2 storage in saline aquifers[END_REF] which is able to take into account CO 2 reactivity with water. The code is based on the resolution of mass conservation equation in the discretized space. The resolution is made by integral finite difference which allows more flexibility and efficiency in the resolution, particularly the possibility to use unstructured meshes not linked to a global system of coordinates.

TOUGH has several drawbacks, most of them in the chemistry field. The most important in the field of CO 2 injection in deep saline aquifers, is the fact that activity coefficients are computed using extended Debye-Hückel law (see equation (4.32)) which is questionable or even wrong for high ionic strengths (even THOUGHREACT includes the B-dot activity model [START_REF] Hc | Theoretical prediction of the thermodynamic behavior of aqueous-electrolytes at high-presures and temperatures .4. calculation of activity-coefficients, osmotic coefficients, and apparent molal and standard and relative partial molal properties to 600-degrees-c and 5kb[END_REF] which allows to model activity coefficients up to relatively high ionic strengths but only for brines where NaCl is the major electrolyte). Evaporation/precipitation process following CO 2 injection leads to very high ionic strengths which causes strong deviation of the chemical behavior of aqueous species compared to ideality. On of the most robust model for highly saline systems is the Pitzer model but the related thermodynamic database is not extended enough, which limits its application. A new version of the TOUGH code including the Pitzer formalism is currently developed at the Lawrence Berkeley National Laboratory but it has not been released yet.

Presentation of different simulations

One of the first articles on CO 2 storage simulation is "Multiphase flow dynamics during CO 2 disposal into saline aquifers" [START_REF] Pruess | Multiphase flow dynamics during CO2 disposal into saline aquifers[END_REF]. This article explains the simulation of CO 2 injection, using TOUGH2 and the EWASG module (the predecessor of ECO2N). In this article we find the basis of the CO 2 injection simulations: a cylindrical infinite reservoir and permeability evolution through salt precipitation. The 2D problem is solved as a 1D one, thanks to the similarity variable [START_REF] O'sullivan | A similarity method for geothermal well test analysis[END_REF]: ξ = r 2 t , r being the distance to the injection well and t the time. The CO 2 considered in this simulation is below the critical point. The behavior is supposed isothermal and with no gravity effect and no chemical reaction (except the precipitation in the dry-out zone). The authors highlight three zones in the aquifer: inner region with complete dry-out and precipitation, intermediate region with two-phase flow and finally a non impacted region with only saline water. Other similar simulations under the same hypotheses (isothermal, no gravity and no chemistry apart the precipitation in the dry-out zone) have been run using other codes with the same conclusions: [START_REF] Prévost | Numerical modeling of carbon dioixide injection and transport in deep saline aquifers[END_REF] with Dynaflow or [START_REF] Hurter | Impact of mutual solubility of H2O and CO2 on injection operations for geological storage of CO2[END_REF] with Eclipse and [START_REF] Müller | CO2 injection impairment due to halite precipitation[END_REF] for a stream-line based simulator. Halite precipitation has also been studied in the case of carbon dioxide injection within depleted gas reservoirs [START_REF] Giorgis | 2D modeling of salt precipitation during the injection of dry CO2 in a depleted gas reservoir[END_REF] leading to similar conclusions.

Finally, [START_REF] Alkan | Impact of capillary pressure, salinity and in situ conditions on CO2 injection into saline aquifers[END_REF] studied the impact of different parameters on the simulation (TOUGH2 with ECO2N) with the same hypotheses as before, but including gravity. They studied the effect of capillary pressure, salinity and in situ thermodynamic conditions. Capillary pressure, viscous forces and gravity have an effect on the radius of the dry-out zone (buoyancy for example makes the CO 2 move up and the radius of the zone decreases when depth increases). Among the most complete simulations are those from [START_REF] André | Numerical modeling of fluid rock chemical interactions at the supercritical CO2-liquid interface during CO2 injection into a carbonate reservoir, the Dogger aquifer (Paris Basin, France)[END_REF] and [4] including chemical reactions and thermal effects for the latter. The effect of salt precipitation is actually partially balanced by mineral dissolution.

The results of these simulations and the inferred behavior of the aquifer, will then serve as a basis for our modelling of the crystallization of salt during CO 2 injection.

Summary of the simulations: the THC Behavior

First of all, carbon dioxide displaces water and partitioning occurs: CO 2 dissolves in the brine and water evaporates into the carbon dioxide rich phase. During the transport, chemical equilibria are disturbed and several reactions can occur: dissolution of the minerals constituting the rock matrix and re-precipitation, or incorporation of carbon dioxide within newly formed minerals in the rock. While fluids are flowing, water saturation and thus water relative permeability decreases. Once the residual saturation is reached, water flow stops. Because the carbon dioxide injected is dry, evaporation continues, and water saturation decreases below the residual saturation and finally reaches zero. In this zone, called the dry-out zone (zone 5 with S g = 1 on figure 6.1), there is no water left, only rock and carbon dioxide. Before this drying-out, evaporation increases the concentration of the dissolved salts in the brine, leading to very high ionic strengths. Brine is then supersaturated with respect to the salts which begin to precipitate, leading to a huge amount of crystallization.

The region undergoing evaporation (zones 3 and 4 on figure 6.1), which are the object of study of this PhD work, has a very low water saturation compared to the rest of the aquifer because of evaporation. Strong capillary currents will then take place between this zone and the rest of the aquifer, bringing back brine to the dessicated zone. These counter-currents have two main effects. First they limit the size of the dry-out zone: the zone cannot expand infinitely and is restricted by the equilibrium between the evaporation rate, water displacement rate induced by carbon dioxide, and the counter-currents, and thus by the carbon dioxide flow rate. The other effect of these counter-currents has huge implications. Indeed, the brine brought back by capillary pressure is still full of salts, increasing the total amount of salt in the dry-out zone and leading to a really localized precipitation at the border of the dry-out zone [START_REF] André | Impact of porous medium dessication during anhydrous CO2 injection in deep saline aquifer: up scaling from experimental results at laboratory scale to near-well region[END_REF]. Clogging and permeability reduction is then dramatically worsened. Finally, it is important to remember that a field is in 3D and that depth also plays an important role through gravity. Carbon dioxide is lighter than brine, and the plume has a tendency to rise while water sinks creating a convection cell as represented in figure 6.2. Finally a cross-section of the aquifer from above is represented in figure 6.1.

Fundamental hypothesis 7.1.1 Presentation and definitions

Poromechanics is the branch of the mechanics of deformable solids applied to porous media saturated with fluids. The presence of fluids inside the porous network modifies the classical equilibrium of deformable solids. The first work on this field has been done by Biot and extended by O. Coussy [START_REF] Coussy | Mechanics and Physics of Porous Solids[END_REF][START_REF] Coussy | Unsaturated poroelasticity for crystallization in pores[END_REF]. This is the approach adopted here. We will consider an elastic model meaning no energy dissipation during the deformation. The material will then be considered as brittle, leading to a rupture without any viscoplastic effect. All the elastic energy stored will be released when fracture appears.

Let us consider a porous system consisting of four constituents: the solid matrix, the precipitated salt(s), and the two fluids, brine and carbon dioxide. We will use the following definitions: the porous medium consists of all the phases; the skeleton (sk) is the porous medium from which the in-pore bulk phases (brine, CO 2 , and crystallized salts) have been removed. It is important to notice that the interfaces belong to the skeleton. It is then composed of the matrix (empty material) (m) and the interfaces.

Frame of the Study and Constitutive Hypothesis

The object we study is an elementary volume Ω of porous medium. It is named Representative Elementary Volume (REV) because its size is small enough to be elementary, but big enough to be representative of the global porous medium: if l h is the characteristic length of the heterogeneities, l Ω the characteristic length of the REV, and L a characteristic length of the structure containing Ω, we must have:

l h ≪ l Ω ≪ L (7.1)
Moreover, we consider the problem under the approximation of small perturbations, which means:

∇ξ ≪ 1 ; ξ l Ω ≪ 1 ; |φ -φ 0 | ≪ 1 (7.2)
ξ is the displacement field of the skeleton.

Under this approximation, the strain is linearised as:

= ε = 1 2 ∇ξ + ∇ T ξ (7.3)

CONSTITUTIVE EQUATIONS OF UNSATURATED POROELASTICITY

As argued in chapter 5, temperature can be assumed constant and homogeneous.

Constitutive equations of unsaturated poroelasticity

In this section, we review the fundamental equations of poromechanics. Two postulates of local state are mandatory in order to derive the equations of poromechanics: the first one considers that the thermodynamic state of an homogeneous system is characterized by the same variables as for equilibrium. For a reversible evolution, equations of state are sufficient while for a non reversible transformation, complementary behavior laws has to be added, depending on the speed of evolution of the state variables. The second one says that if we consider a volume Ω constituted of elementary volumes dΩ able to exchange work, heat and matter, the thermodynamic state of Ω is equal to the sum of the thermodynamic states of the elementary volumes.

In order to calculate the evolution of the porous volume, we will consider the Helmholtz free energy state function Ψ = U -T S . The (isothermal) infinitesimal evolution of Ψ is expressed as follows:

dψ = dw + i µ i dn i (7.4)
with ψ = dΨ dΩ 0 the volumetric free energy, dw, the infinitesimal strain work density to the domain dΩ 0 , µ i the chemical potential of species i. Expressing the strain work in terms of stress and strain, we obtain:

dψ = = Σ : = dε + i µ i dn i (7.5)
This relation expresses the variation of Helmholtz free energy for the whole porous solid. If we only consider the skeleton, its free energy is:

ψ sk = ψ -ψ in-pore (7.6)
with dψ in-pore is the Helmholtz free energy of the in-pore bulk phases. Considering equations (4.1) and (4.4), we can then express the in-pore bulk free energy variation as dψ in-pore = i µ i dn i -J p j dφ J and thus obtain the infinitesimal variation of the free energy of the skeleton (under isothermal condition):

CHAPTER 7. BASICS OF POROMECHANICS

dψ sk = = Σ : = dε + J p J dφ J (7.7)
if we suppose a homogeneous pressure p J for each phase J (even for crystallized salt as argued in [START_REF] Coussy | Mechanics and Physics of Porous Solids[END_REF] and [START_REF] Osselin | Dependence on injection temperature and on aquifer's petrophysical properties of the local stress applying on the pore wall of a crystallized pore in the context of co 2 storage in deep saline aquifers[END_REF]).

Considering that dφ J = dϕ J + φ 0 dS J (3.7), we finally obtain:

dψ sk = = Σ : d = ε + J p J dϕ J + φ 0 J p J dS J (7.8)
As equation (7.8) is an exact differential of the free energy of the solid skeleton, we can obtain the following equations of state:

= Σ = ∂ψ sk ∂ = ε ; p J = ∂ψ sk ∂ϕ J ; φ 0 p J = ∂ψ sk ∂S J (7.9)
However, as recalled in section 3.3.1, interfacial energies act at the interfaces of the different phases. To express these interfacial energies, we will consider that we can decompose the free energy of the skeleton as the sum of the solid matrix free energy W and the energy of the interfaces U , considering that U is independent of the skeleton strain (such hypothesis has been validated in the case of freezing [START_REF] Scherer | Poromechanics of frost damage[END_REF][START_REF] Sun | Effect of air voids on salt scaling and internal freezing[END_REF]):

ψ sk = ε, ϕ J , S J = W = ε, ϕ J , S J + U (S J , ϕ J ) (7.10)
If we consider infinitesimal transformations and linear behavior, we can neglect the evolution of the solid matrix free energy with saturation in front of the evolution of the interfacial energies with saturation (∂W/∂S J ≪ ∂U/∂S J ). Using now equation (7.10) in (7.9), we obtain the following equations:

= Σ = ∂W ∂ = ε ; p * J = ∂W ∂ϕ J ; φ 0 p J = ∂U ∂S J (7.11)
where:

p * J = p J - ∂U ∂ϕ J (7.12)
is the pressure effectively transmitted to the rock matrix through the internal solid walls delimiting the porous volume occupied by the constituent J.

Finally, we can use the Legendre transform W = W -J p * J ϕ J to obtain the following 7.3. APPLICATION OF THE EQUATIONS TO THE CCS CASE equations of state:

= Σ = ∂ W ∂ = ε ; ϕ J = - ∂ W ∂p * J (7.13)
Considering only unsaturated isotropic linear poroelastic solids, the energy W will be a quadratic function of the pore pressures p J and of the two invariants of the strain tensor.

As a result we can write:

W = 1 2   Kǫ 2 -2 J b J p * j -p * J,0 ǫ - I,J p * I -p * I,0 p * J -p * J,0 N I,J + 2µ = e : = e   + = Σ 0 : = ε (7.14)
where the subscript 0 refers to the reference state condition, ǫ = tr If we combine the state laws (7.13) and the expression (7.14), we finally obtain the constitutive equations for unsaturated poroelasticity:

Σ -Σ 0 = Kǫ - J b J (p * J -p * J,0 ) (7.15a) = s - = s 0 = 2µ = e (7.15b) ϕ J = b J ǫ + I 1 N I,J (p * I -p * I,0 ) (7.15c)

Application of the equations to the CCS case

Let us consider that we are in a 1D problem (the carbon dioxide is flowing in only one direction) as it is usual in the simulations (section 6.2). This allows us to write the strain tensor as (similar to an oedometer test):

= ε =     ε xx 0 0 0 0 0 0 0 0     (7.16)
which allows to rewrite equation (7.15a) and (7.15b) as:

Σ xx -Σ xx,0 = K + 4 3 µ ε xx - J b J p * J -p * J,0 (7.17)
The other diagonal terms of the stress tensor (Σ yy and Σ zz ) are non-zero (deviatoric stress).

Finally, the mechanical equilibrium imposes ∇ • = Σ = 0 (neglecting gravity), which leads to:

Σ xx = cste (7.18)
Setting as limit conditions Σ xx = Σ xx,0 at x = 0, we have that Σ xx = Σ xx,0 , ∀x. We finally obtain the strain as:

ε xx = J b J p * J -p * J,0 K + 4 3 µ (7.19)
Considering the relations b J = S J b and 1 N IJ = S I S J N [START_REF] Fen-Chong | Poroelastic analysis of partial freezing in cohesive porous materials[END_REF], with b = 1 -K K M the Biot coefficient (K M is the bulk modulus of the matrix), and 1/N = (b -φ 0 )/K M the Biot modulus, we can express the strain of equation (7.19) as:

ε xx = b K + 4/3µ J S J p * J -p * J,0 (7.20)
The elastic energy stored in the matrix is finally:

W = 1 2 b 2 K + (4/3)µ + 1 N J S J p * J -p * J,0 2 (7.21)
As defined in [START_REF] Coussy | Deformation and stress from in-pore drying-induced crystallization of salt[END_REF], the equivalent macroscopic tensile stress is defined as the hypothetical tensile stress needed to reach the elastic energy calculated in equation (7.21):

APPLICATION OF THE EQUATIONS TO THE CCS CASE

̟ = 2 K + 4 3 µ W = b 2 + K + (4/3)µ N J S J p * J -p * J,0 (7.22)
Let us now consider the three phases in CCS: s = crystal, l=brine and G=supercritical CO 2 . Following [START_REF] Coussy | Mechanics and Physics of Porous Solids[END_REF][START_REF] Fabbri | Physico-Mécanique des matériaux cimentaires soumis au gel-dégel[END_REF][START_REF] Fabbri | Indirect measurement of the ice content curve of partially frozen cement based materials[END_REF], we will neglect the derivatives ∂U ∂ϕ J for the liquid phase and the gas phase. Accordingly, considering an initial state at the beginning of the evaporation (no crystal), we can rewrite the equivalent tensile stress as:

̟ = b 2 + K + (4/3)µ N S s p * s + S l (p l -p 0 l ) + S G (p G -p 0 G (7.23)
Determination of p * s is linked to the crystallization pressure of the growing crystal and is detailed in the following chapter.

Chapter 8

Thermodynamics of in-pore crystallization In the following, we will use the subscript j for the species concerned by the equation of reaction (aqueous species, crystal and water), and the subscript i for the aqueous species only (i.e. the j without water and crystal). The dissolution equation of the crystal is then:

A ν A B ν B . . . • ν w H 2 O ⇋ ν A A + ν B B + . . . + ν w H 2 O (8.1)
where the stoichiometric coefficients satisfy both the conservation of matter and charge.

Let us now consider the thermochemistry associated with this reaction. The Gibbs free energy of this reaction can be written with equation (4.18):

∆ r G = j ν j µ j (8.2)
where µ j is the chemical potential of the species j as defined in section 4.1.1.3.

At equilibrium, ∆ r G = 0 which leads to:

i ν i µ i + ν w µ w = ν s µ s (8.3)
Let us now consider a crystal which does not have the same pressure as the surrounding solution (this can come from its size or from external interactions as we will see in the following). Its chemical potential is then modified through its molar volume:

∂µ s ∂p T,ns = υ s (p, T ) (8.4)
with υ s the molar volume of the crystal at the pressure and temperature considered.

In the following we will neglect the variation of the molar volume of the crystal with temperature and pressure in the range considered. As a consequence of this linearisation, we can write:

∆µ s = µ s (p s , T ) -µ s (p l , T ) = υ s (p s -p l ) (8.5)
8.1. OSTWALD-FREUNDLICH EQUATION AND WULFF THEOREM ∆µ s is the difference of chemical potential between a crystal at the pressure p s and a crystal at a pressure p l . The Gibbs free energy of the reaction of dissolution (8.2) can then be expanded as:

∆ r G = ∆ r G 0 (p l , T ) -υ s (p s -p l ) + RT ln Q R (8.6)
This leads to the modification of the law of mass action:

p s -p l = RT υ s ln S(p l , T ) (8.7)
with the supersaturation S classically defined as S = j a ν j j

Ks(p l ,T ) (see equation (4.23)), K s being the isobar equilibrium constant at the liquid pressure. Equation (8.7) describes the equilibrium between a crystal at pressure p s surrounded by a solution at pressure p l and supersaturation S. This relation is strongly dependent on the kind of crystal considered (through the factor 1/υ s ), but the trend is the same for all crystals: the more the crystal is under stress, the more soluble it is at fixed p l , T , and ions concentrations. This phenomenon is well known in geosciences and is called pressure solution [START_REF] Weyl | Pressure solution and the force of crystallization -a phenomenological theory[END_REF]: an increase in crystal pressure will lead to a dissolution of the crystal; if a crystal (or parts of the crystal) is more stressed than another because of the distribution of mechanical stresses (or external forces), it will dissolve in favor to the less stressed crystal (or other crystal parts) because of the diffusion of dissolved ions from high stressed parts to low stressed parts. Finally for negative overpressure (corresponding to supersaturations below 1), the crystal is not stable and dissolves completely.

Ostwald-Freundlich equation for small crystallite

Growth of a crystal in a supersaturated solution results from the competition of an over concentration which tends to make it grow and the surface energy which tends to make it shrink. If we consider a crystal with curved interfaces, we can link its inner pressure to the curvature κ s of its interface as:

p s -p l = σ dA dV = σκ s (8.8)
with dA and dV the variation of the surface and the volume of added solid (during growth) [START_REF] Steiger | Crystal growth in porous materials -II: Influence of crystal size on the crystallization pressure[END_REF], and σ the surface energy of the crystal considered. This surface energy should not be mistaken with the surface tension of a liquid implied in Laplace equation. Indeed, as developed in [START_REF] Scherer | Crystallization in pores[END_REF][START_REF] Scherer | Stress from crystallization of salt[END_REF][START_REF] Vandamme | Adsorption and strain: the CO2-induced swelling of coal[END_REF] (and in recalled in section 3.3.1), surface tension character-izes stretching of the interface while surface energy characterises the addition of matter (growth) to the interface. Combining equations (8.7) and (8.8), we obtain the Ostwald-Freundlich equation which links the curvature of the crystal to the current supersaturation of the solution at equilibrium:

σκ s = RT υ s ln S. (8.9)
For a spherical crystallite, κ s = 2 r , which leads to:

2σ r = RT υ s ln S (8.10)
An interesting but unsurprising remark is that the equilibrium radius given by equation (8.10) corresponds to the critical radius for nucleation (see equation (4.45)). If a crystal is in equilibrium in a solution, it must have a curvature corresponding to the current supersaturation. The less concentrated the solution, the bigger the crystal has to be, in order to be in equilibrium with the solution, the limit being a crystal of infinite size (i.e. macroscopic size) for supersaturation equal to 1 (leading also to p s = p l ). However, this equilibrium is unstable as we will see below (section 8.3).

Wulff theorem and equilibrium shape

The Ostwald-Freundlich equation considers that the crystal is a sphere in the conditions of homogeneous surrounding supersaturation. However, it is known that crystals which are allowed to reach their equilibrium shape do not shape like spheres, but rather like polyhedron composed of planes, which underlines the microscopic structure of the crystalline arrangement of atoms in the crystal [START_REF] Sekerka | Equilibrium and growth shapes of crystals: how they differ and why should we care[END_REF]. In this equilibrium shape, each face i is at a precise distance h i from the geometric center of the polyhedron called the Wulff length (see figure 8.1). These lengths are linked to the surface energy of the considered phase through the following relation:

σ 1 h 1 = . . . = σ N h N = constant (8.11)
with σ i the surface energy of the face i. The proportionality constant in equation (8.11) represents the work needed to displace the interface of a length dh i and thus depends on the difference of pressure across the considered face. The Gibbs-Curie-Wulff theorem (demonstration in appendix B) gives this constant as:

8.1. OSTWALD-FREUNDLICH EQUATION AND WULFF THEOREM

∆µ s = υ s (p s -p l ) = 2σ i υ s h i (8.12)
Finally, combining equation (8.12), with equation (8.7), we obtain the Ostwald-Freundlich equation for crystals in their Wulff equilibrium shape:

p s -p l = RT υ s ln S = 2σ i h i (8.13)
Figure 8.1: Illustration of the Wulff shape for a cubic halite crystal.

Non-flat surfaces

In the previous section, crystals were considered in their equilibrium shape, which supposes that crystals have sufficient time to grow and thus their growth is controlled by thermodynamic and not kinetics. However, in reality crystals will grow more or less rapidly depending on the supersaturation. For fast growth, it is likely that the surface will not be smoothly flat but will rather present a step/terrace shape as in figure 8.2 [START_REF] Markov | Crystal Growth for beginners: fundamental of nucleation, crystal growth, and epitaxy[END_REF][START_REF] Sekerka | Equilibrium and growth shapes of crystals: how they differ and why should we care[END_REF]. We consider that the surface energy of the surface can be calculated from the surface energy of the steps and that of the terraces. As a consequence, surface energy of a surface will depend on the point considered on this surface. If we consider a point on a surface determined by the spherical coordinates (r, θ 1 , θ 2 ), the surface energy associated to this point is σ(r, θ 1 , θ 2 ), and the Gibbs-Curie-Wulff theorem is modified to become Herring's formula as follows [START_REF] Markov | Crystal Growth for beginners: fundamental of nucleation, crystal growth, and epitaxy[END_REF]:

∆µ s = υ s (p s -p l ) = υ s R 1 σ + d 2 σ dθ 2 1 + υ s R 2 σ + d 2 σ dθ 2 2 (8.14)
with R 1 and R 2 , the two principal radii of curvature. The quantity σ = σ + d 2 σ dθ 2 is called the surface stiffness. It corresponds to the surface resistance to bending at a constant number of atoms. If we consider a flat surface, this means that the two radii of curvature are infinite. This is only possible when the surface stiffness is also infinite. A flat surface (which is also called a singular face) is a surface with an infinite resistance to bending.

It is interesting to consider the relation between Ostwald-Freundlich equation (8.9) and Herring formula (8.14). Ostwald-Freundlich equation is a particular case of Herring formula in which the surface stiffness is taken equal to the surface energy ( d 2 σ dθ 2 ). In other terms it means that the crystal surface has no resistance to bending during growth (much like a liquid) and can thus grow to any shape corresponding to the solicitations: we will call them liquid-like crystals. On the other end of the spectrum the Wulff theorem considers that all surfaces have infinite resistance to bending and thus are all flat. Finally, the surface stiffness described here is a consequence of growth and should not be mistaken with the surface stress as in Shuttleworth equation (3.16).

Interaction energy between two solid surfaces 8.2.1 Equilibrium of a crystal surface close to the pore wall

Let us now consider two infinite flat crystal surfaces separated by a thin film of solution containing ions of the same substance. If the two solid faces are close enough, their energy profiles will overlap and the two surfaces will interact (see figure 8.3). Let us call W (e) the interaction energy between the two surfaces, e being the thickness of the film of solution. We will consider that this film is thick enough to allow ions to diffuse within and thus maintain a constant supersaturation of the solution1 . The interaction energy is bound by the two following conditions [START_REF] De Gennes | Capillarity and wetting phenomena: drops, bubbles, pearls, waves[END_REF]:

W (∞) = 0 ; W (0) = S = σ 1,2 -σ 1 -σ 2 (8.15)

INTERACTION ENERGY BETWEEN TWO SOLID SURFACES

with σ 1 (respectively σ 2 ) the interfacial energy between the solid 1 (respectively 2) and the solution; σ 1,2 is the interfacial energy between the two solids in contact, S is called the spreading coefficient. The existence of the film is linked to the sign of this spreading coefficient. If S is negative, the contact between the crystal and the pore wall is energetically favorable and no film of solution will exist. On the contrary, is S is positive, it means that the contact between the crystal and the pore wall is energetically unfavorable and a repulsive force will arise when the two surfaces are pushed together, leaving the film of solution to separate both. The presence of this repulsive force will modify the equilibrium of the solids. To analyse the case of solid 1 as pore wall and solid 2 as crystal growing toward this pore, we extend the thermodynamic treatment of the premelted film and crystallization used in [START_REF] Coussy | Mechanics and Physics of Porous Solids[END_REF]. Let us consider the free energy cost per unit surface associated to the change of thickness of the film from e to e + de by dissolution or precipitation of the crystal. Similarly to the nucleation cost in equation (4.43), this cost can be written as:

dG(e) = i µ i (p l , T )dn i + µ w (p l , T )dn w + µ s (p l , T )dn s + dW (e) (8.16)
considering that the change in thickness does not modify the liquid free energy. Note that here the mole quantity n are here mole quantity per unit surface (mol/m 2 ) for homogeneity of equation (8.16) .

Considering the stoichiometry of the reaction, we have dn i ν i = dnw νw = -dn s = de υs , an increase of the thickness e meaning a dissolution of the crystal. This leads to:

CHAPTER 8. THERMODYNAMICS OF IN-PORE CRYSTALLIZATION dG(e) = de υ s i ν i µ i (p l , T ) + ν w µ w (p l , T ) -µ s (p l , T ) + dW (e) (8.17)
The term between the brackets corresponds to the energy cost of the reaction for a crystal being at the same pressure as the liquid. The second term corresponds to the influence of the interaction energy on the equilibrium. Following [START_REF] Coussy | Mechanics and Physics of Porous Solids[END_REF][START_REF] Markov | Crystal Growth for beginners: fundamental of nucleation, crystal growth, and epitaxy[END_REF], we define the degree of metastability also called supersaturation (not to be mistaken with the supersaturation S) ∆µ = i ν i µ i (p l , T )+ν w µ w -µ s (p l , T ). This quantity measures the difference between the (molar) chemical potential of the crystal in solution and the (molar) chemical potential of the precipitated crystal at the same pressure. The cost associated to the change of thickness is then:

dG(e) = de υ s ∆µ + dW (e) (8.18)
Equation ( 8.16) can finally be rewritten:

dG(e) de = ∆µ υ s -ω p (8.19)
The term ω p = -dW (e) de is akin to a pressure and is called disjoining pressure. As W (e) decreases from W (0) = S > 0 to W (∞) = 0, dW/de is negative and the disjoining pressure is in this case always positive (in the case of a negative spreading parameter, the two surfaces would be attracted to each other and the interaction energy would be decreasing with decreasing thickness. In this case the quantity ω p would be called joining pressure because the interaction is attractive).

Let us now consider the sign of the quantity dG(e)/de. If the degree of metastability is negative, dG(e)/de is always negative. As a result the system will be more and more stable as the thickness of the film e increases. The presence of the pore wall and the degree of metastability act both toward a dissolution of the crystal, the thickness of the film will increase indefinitely. On the other hand, if ∆µ is positive, there are two cases. For small values of ∆µ, dG(e)/de is negative (dissolution), and for high values of ∆µ, dG(e)/de is positive (growth). There is then an equilibrium situation corresponding to dG(e)/de = 0. This equilibrium is reached for a value of the thickness e = δ. We can then write:

ω p (δ) = ∆µ υ s (8.20)
In absence of any external forces (i.e. ω p = 0), the degree of metastability reduces

INTERACTION ENERGY BETWEEN TWO SOLID SURFACES

to the free enthalpy of reaction ∆ r G and equation (8.20) reduces to ∆ r G = 0. The disjoining pressure, by adding an external pressure on the crystal surface displaces the equilibrium from the equilibrium of the same crystal in a bulk solution. The degree of metastability scales this displacement from the bulk equilibrium: any change in the degree of metastability will modify the equilibrium thickness.

Considering the equilibrium, we have as usual:

∆ r G = 0 = i ν i µ i (p l , T ) + ν w µ w (p l , T ) -µ s (p s , T ) (8.21)
Injecting the degree of metastability in equation (8.21), we obtain:

0 = µ s (p l , T ) + ∆µ -µ s (p s , T ) (8.22)
which leads to:

∆µ = µ s (p s , T ) -µ s (p l , T ) (8.23)
The degree of metastability can be interpreted as the difference of chemical potential between a crystal at a pressure p s and the same crystal at the pressure p l as expressed in equation (8.5). Once again if we consider the crystal out of interaction from the pore wall, ∆µ = ∆ r G = 0 which leads to p s = p l (infinite surface).

Finally, using equations (8.5), (8.20) and (8.23), we obtain:

p s -p l = ω p (δ) (8.24)
The disjoining pressure modifies the pressure of the crystal which displaces the equilibrium much like the size of the crystal modifies its inner pressure and changes the equilibrium. Finally, the crystal pressure is then p s = p l + ω p (δ).

Considering now equation (8.7), we obtain:

ω p (δ) = RT υ s ln S (8.25)
The equilibrium disjoining pressure (or the equilibrium thickness δ) is then linked to the current supersaturation of the solution. The higher the supersaturation is, the closer the crystal will be to the pore wall. Equation (8.25) gives the link between the disjoining pressure and the current supersaturation of the solution and has been calculated in the case of two infinite flat surface in interaction through a liquid film of solution. However, in the case of porous media, surfaces are often non-flat because of the curvature of the pores. As a result, equation (8.25) has to be modified to take into account the effect of curvature on the crystal pressure. This calculation is made in section 8.3.

Expression of the interaction energy

Interactions between two surfaces can be divided into two types: DLVO forces (from the names of Derjaguin, Landau, Verwey and Overbeek) and the non DLVO forces. DLVO forces are the classical long range forces that include Van der Waals interactions and electrostatic forces.

DLVO forces

Van der Waals forces: These forces exist between molecules or atoms. They are responsible of the cohesion of solids and liquids (along with covalent bonding and hydrogen bonding). Van der Waals forces can be classified into three categories depending on the polarity of the considered molecules: Keesom interaction, Debye interaction and London interaction (the intensity of the forces is decreasing in this order).

• Keesom interaction exists between polar molecules and the forces is created by the electric field surrounding each molecule (this force is also called orientation effect).

• Debye interaction appears between a polar molecule and a non-polar one. The electric field of the polar molecule induces a polarization in the non-polar one, which creates the interaction.

• London force which exists between two non-polar molecules is created because of the variation of the position of the electrons surrounding the atoms. Their fluctuation of position create sometimes an instantaneous polarization of the molecule and finally an electric field and an interaction. Because of the statistics implied by the London interaction, its intensity is much lower than the two previous, but as this interaction exists for all entities, it is the most important of the Van der Waals interaction.

Electrostatic interaction Contrary to the Van der Waals forces, this interaction acts between two charged surfaces through a ionic solution. This interaction is screened by the ions in the solution and its range can be calculated with the Debye length (here for a 1:1 electrolyte such as NaCl):

8.2. INTERACTION ENERGY BETWEEN TWO SOLID SURFACES κ -1 = ε r ε 0 RT 2F 2 C 0 (8.26)
with ε r the dielectric constant, ε 0 , the permittivity of free space, R the perfect gas constant, T the temperature, F the Faraday constant, and C 0 the concentration of the solution. In our case, we are dealing with highly concentrated electrolytes (S > 1) leading to a Debye length magnitude smaller than 1

• A.
Considering that the film of solution between the two surfaces is usually of 1 nm (section 8.3.4), we can consider that the electrostatic interaction can be neglected.

Non-DLVO forces

At high concentration (i.e. small Debye length), it is observed experimentally that the interaction between two surfaces departs significantly from the DLVO form. The reason of this departure are the so-called non-DLVO forces usually called solvation forces or hydration forces. These forces are repulsive and arise from the steric constraint at small thickness of the solvent molecules. These forces are very complicated to model theoretically, but can be measured experimentally [START_REF] Israelachvili | Intermolecular and Surface forces[END_REF]. It appears that they have an exponential form with two parameters:

W (e) = W 0 exp - e λ 0 (8.27)
If we consider the limit condition W (0) = S, we can identify the prefactor W 0 as the spreading parameter. The disjoining pressure becomes then:

ω p (e) = S λ 0 exp - e λ 0 (8.28)
According to [START_REF] Israelachvili | Intermolecular and Surface forces[END_REF], the parameter λ 0 ranges between 0.6 and 1.1 nm for 1:1 electrolytes. This parameter decreases with increasing concentration of the solution. As we deal with highly concentrated solutions, we will take λ = 0.6 nm in the following. Surface energy of halite is found in the literature to be σ = 0.038 J.m -2 . We can also consider different kinds of pore wall according to the type of reservoir rock considered. For limestones, the pore wall will mainly be composed of calcite CaCO 3 , while for sandstones, the pore wall will be made of silica SiO 2 . Literature reports a large dispersion of values for calcite surface energy (depending on the crystalline surface considered), from 0.039 J/m -2 to 0.15 J.m -2 . Moreover, to our knowledge, there is no experimental value for the interfacial energy between calcite and halite (or silica). As a conclusion, there is no experimental value of the spreading parameter S. Following [START_REF] Israelachvili | Intermolecular and Surface forces[END_REF], we will take for the spreading parameter a CHAPTER 8. THERMODYNAMICS OF IN-PORE CRYSTALLIZATION mean value between the two classic boundaries (3 -30 mJ.m -2 ). A recent measurement of the interaction forces with conventional surface force apparatus [START_REF] Espinosa-Marzal | Confinement during in-pore crystallization[END_REF] during confined crystallization confirms the effect of non-DLVO forces between the growing and confining surfaces. Far from the pore wall, there is no interaction between the two solid surfaces: the crystal (with isotropic surface energy) has a spherical shape. When the crystal is big enough, the presence of the pore wall forces it to adopt a cigar-shape as shown in figure 8.4. It exhibits then two different curvatures κ side = 1 rp-δ and κ tips , with r p the radius of the pore and δ the thickness of the film. The constrained side of the crystal being cylindrical and in equilibrium with the liquid film at ω p + p l (section 8.2), equation (8.8) becomes:

p s -(p l + ω p ) = σκ side = σ r p -δ (8.29)
Similarly, the tips being in equilibrium with the sole liquid, we have from (8.8):

IN-PORE GROWTH OF A SINGLE CRYSTAL AND CRYSTALLIZATION PRESSURE

p s -p l = σκ tips (8.30) As the crystal is under homogeneous pressure, we can obtain the expression for the disjoining pressure:

ω p = σ (κ tips -κ side ) (8.31)
In this equation, we have two unknowns: the curvature of the tips and the disjoining pressure. To determine the latter, we will use the fact that the side of the crystal does not grow anymore because of the presence of the pore wall. This results in a forced equilibrium situation with the thin film of solution at supersaturation S similarly to section 8.2. Injection of equation (8.7), p s -p l = RT υs ln S(p l , T ), in (8.29) leads to:

ω p = RT υ s ln S(p l , T ) - σ r p -δ (8.32)
which is the classical expression of the so-called crystallization pressure [START_REF] Scherer | Crystallization in pores[END_REF][START_REF] Scherer | Stress from crystallization of salt[END_REF][START_REF] Flatt | Salt damage in porous material: How high supersaturations are generated[END_REF][START_REF] Steiger | Crystal growth in porous materials -I: The crystallization pressure of large crystals[END_REF] (among others). It is important to notice that this relation has been obtained by considering only the equilibrium of the side of the crystal. This relation which is always valid as soon as the crystal has a cigar-shape does not rely on any equilibrium of the tips of the crystals. We can now obtain the curvature of the tips with equations (8.30) and (8.7):

σκ tips = RT υ s ln S (8.33)
Considering that the tips are hemispherical, we finally obtain:

2σ r tip = RT υ s ln S(p l , T ) (8.34)

Second case: crystallization pressure in a spherical pore with small entry channels

As we have seen in previously, calculation of crystallization pressure is strongly geometrydependent. Indeed, let us consider the geometry presented in figure 8.5: a spherical pore with two small cylindrical entry channels in which is growing a liquid-like crystal. When the crystal arrives at interaction distance from the pore wall, it has at first a spherical shape of radius r p -δ (figure 8.5(a)). Then it begins to invade the side channels with a curvature depending on the current supersaturation as described in equation (8.10).

CHAPTER 8. THERMODYNAMICS OF IN-PORE CRYSTALLIZATION

This is depicted in figures 8.5(b) and 8.5(c), the last one representing the classical case of curvature corresponding to the radius of the entry channel. Crystallization pressure is created as soon as the crystal reaches the pore wall i.e. as soon as the configuration described in figure 8.5(a) is obtained. Equations (8.29) and (8.30) can then be written:

p s -p l = σκ tips p s -(p l + ω p ) = 2σ r p -δ (8.35)
κ tips , the curvature of the tips is coloured in red in figure 8.5. Once again, we do not make any hypothesis on the shape of the crystal growing in the entry channel. The crystallization pressure becomes then:

ω p = σ κ tips - 2 r p -δ (8.36)
Considering now the impeded growth condition of the spherical part, as previously, we have the value of the crystallization pressure:

ω p = RT υ s ln S - 2σ r p -δ (8.37)
If figure 8.5(a) is an equilibrium situation (i.e. RT υs ln S = 2σ rp-δ ), there is no crystallization pressure created and the crystal remains perfectly spherical. This is obviously a very peculiar situation. In every other situations, the supersaturation is high enough to allow a curvature smaller than this given by the pore and thus create a crystallization pressure. Considering that the tips are spherical and in equilibrium with the solution we obtain: The penetration of the crystal within the side channel is however not automatic because it forces the crystal to adopt the cylindrical curvature κ channel side = σ rs-δ . As a result, in order to penetrate this channel, the supersaturation has to follow this condition:

RT υ s ln S(p l , T ) > σ r s -δ (8.39)
This condition is automatically filled if the supersaturation is high enough to force the crystallization pressure on the spherical pore, and if the ratio of radii between the channel and the spherical pore is bigger than 1/2. In the other cases, it is possible to have a crystal creating a crystallization pressure on the pore wall while not being able to penetrate the side channels. However, once the crystal has begun to penetrate the channel, it becomes confined exactly as the crystal was confined in the cylindrical pore in section 8.3.1. As a result, the crystal exerts a crystallization pressure on the wall of the side channel following relation (8.32):

ω p = RT υ s ln S(p l , T ) - σ r s -δ (8.40)
Comparison with the literature As we can see, the results obtained here are slightly different from the pure cigar-shape and the spherical pore situation analysed by [START_REF] Scherer | Crystallization in pores[END_REF][START_REF] Scherer | Stress from crystallization of salt[END_REF][START_REF] Flatt | Salt damage in porous material: How high supersaturations are generated[END_REF][START_REF] Steiger | Crystal growth in porous materials -I: The crystallization pressure of large crystals[END_REF]. The equations for crystallization pressure are completely similar, however, the phenomena described here are broader than theirs. In their calculation they have considered that the crystal tips for both the cigar-shape and the spherical pore were hemispherical with a radius corresponding to the pore radius (as presented in figure 8.6 or in figure 8.5(c)). This representation comes from the analogy with bubbles or drops but is not necessarily happening in the case of liquid-like crystals as we have seen previously.

In [START_REF] Scherer | Crystallization in pores[END_REF] the radius of curvature is linked to the wetting properties of the crystal and to the contact angle θ. However, for crystallization pressure to occur, there must be a film of liquid between the crystal and the pore wall as presented in section 8.2. Because of this film, the crystal is separated from the pore wall, which forbids to define a contact angle and particularly, which does not imply that the contact angle is 180°. In the development presented here, the tips of the crystal can have any curvature corresponding to the supersaturation of the solution, independently from the radius of the pore in which the crystal is growing. The case with hemispherical curvature for the tips with radius matching the pore radius is actually a particular case happening only for precise conditions. Forcing the crystal tips curvatures to match the pore radius leads to the following condition on the tips radius: r tips = r p -δ (respectively r s ). This new condition forces a relationship between the pore radius r p (or r s ) and the supersaturation S: RT υs ln S(p l , T ) = 2σ rp-δ (respectively r s ). Considering an unsaturated porous volume, there is at first approximation a unique supersaturation but a continuous distribution of pore radii. The relation between the supersaturation and the pore radius implies that the crystal precipitate in only one class of pores which is obviously not the case (the corresponding class of pore may not even exist in the considered material). Figure 8.6: Particular case of a crystal growing in a cylindrical pore, with tips radius of curvature corresponding to the pore radius (adapted from [START_REF] Scherer | Crystallization in pores[END_REF]).

Finally, it is important to remember that these two sections are based on the hypothesis that the crystal can grow with any shape (liquid-like crystals). This is wrong in reality. In the following we consider crystals in their Wulff shape and see how the reasoning is changed.

Crystallization pressure and Wulff shape

As we have seen in section 8.1.3 crystals allowed to grow with their equilibrium shape will not have smooth curved surfaces but will rather present a polyhedron shape corresponding to the crystalline arrangement of atoms in the crystal. In the following, we consider a halite crystal (cubic) for simplicity but the development can be adapted for any crystal. As we can see in Appendix B, the Wulff theorem is still valid if a crystal is submitted to different pressures on its faces.

Let us consider the growth of a halite crystal within a pore of square section as presented in figure 8.7. At first, the crystal is out of interaction range of the pore and thus can adopt its cubic shape. When the crystal enters the interaction range, use of Gibbs-Curie-Wulff theorem (8.12) for the crystal at p s in contact with a liquid at p l + ω p , and of chemical equilibrium (8.7) provides with the expression of the disjoining pressure:

ω p = RT υ s ln S - 2σ r p -δ (8.41)
This equation has a form similar to the previous relations in the case of liquid-like crystals.

For the free face, use of equation (8.13) gives:

2σ h * = RT υ s ln S (8.42)
Consequently we can write the crystallization pressure in terms of the geometrical parameters:

ω p = 2σ h * - 2σ r p -δ (8.43)
The positive condition for the disjoining pressure leads to:

h * ≤ r p -δ (8.44)
In contrary to the previous case of the liquid-like crystal, as soon as the crystal reaches the pore wall, its growth in all directions is hindered to maintain stability. This means that the filling of the pore by the crystal will not occur by the growth of a single crystal invading progressively the porous medium, but rather by the nucleation of several small crystals which will grow until they reach the pore wall.

Figure 8.8: Tilted halite crystal at homogeneous pressure confined in a pore with a square section.

The case presented above considers that the crystal grows with faces parallel to the pore wall. If this is not the case, the confined crystal will present the geometry as shown in figure 8.8. This creates a crystalline surface along the pore wall which does not belong to the Wulff shape. As a conclusion, the crystallization pressure equation has to be modified as:

ω p = RT υ s ln S - 2σ ′ r p -δ (8.45)
with σ ′ the value of the surface energy for the new surface.

We have seen that whatever the geometry used to model the pore, the crystallization pressure can be written as a term of supersaturation and a term corresponding to the size of the crystal. Let remark that Steiger [START_REF] Steiger | Crystal growth in porous materials -II: Influence of crystal size on the crystallization pressure[END_REF][START_REF] Steiger | Crystal growth in porous materials -I: The crystallization pressure of large crystals[END_REF] summarized both liquid-like and Wulff shape approaches using the dA/dV term which corresponds to the variation of the surface of the crystal dA with a small addition of matter creating a change in volume dV :

ω p = RT υ s ln S -

Thickness δ and interaction forces

In the section 8.2 we have seen that the disjoining pressure can be expressed as a function of the thickness of the film:

ω p (e) = S λ 0 exp - e λ 0 (8.47)
with e the thickness of the film. At equilibrium, the thickness takes the value δ. Using the previous equation and equation (8.41), we have the following relation:

RT υ s ln S - 2σ r p -δ = S λ 0 exp - δ λ 0 (8.48)
The thickness of the film is then completely determined by the supersaturation of the solution. It is usual to neglect the value of δ before r p . As a conclusion, the equilibrium thickness is:

δ = -λ 0 ln RT υs ln S -2σ rp S λ 0 (8.49)
In figure 8.9 is plotted the evolution of δ with supersaturation for a pore size of 1 µm, λ 0 = 0.6nm, S = 10 mJ.m 2 .

It is interesting to notice that there is a critical supersaturation above which the film disappears creating contact between the crystal and the pore wall. This happens when δ = 0 in equation (8.48):

RT υ s ln S c = S λ 0 + 2σ r p (8.50)
The maximum disjoining pressure reachable before the film disappears is finally:

ω max p = S λ 0 (8.51)
From the data recalled in section 8.2.2.2, the value of the maximum disjoining pressure ranges from 2.7MPa (λ 0 = 1.1nm and S = 3 mJ.m 2 ) to 50MPa (λ 0 = 0.6nm and S = 30 mJ.m 2 ). As a conclusion, crystallization pressure may not be harmful if the interaction is weak between the crystal and the pore, no matter how high the supersaturation is. The cementation of the film with secondary precipitation is experimentally shown in [START_REF] Bouzid | In-pore tensile stress by drying-induced capillary bridges inside porous materials[END_REF].

Effectively transmitted stress

As we have seen in section 8.2, the disjoining pressure is added to the liquid pressure to obtain the stress exerted on the pore wall. In equation (7.23) giving the macroscopic tensile stress, the important parameter is p * s i.e. the effectively transmitted normal stress exerted by the crystals in the REV. Let us consider a porosity deformation ϕ s at constant saturations. The energy balance for the interface between the crystal and the adsorbed water [START_REF] Fabbri | Stress from NaCl crystallization by co 2 injection in aquifers[END_REF] is then:

∂U s,M ∂ϕ s = p s -(p l + ω p ) (8.52)
where U s,M is the energy of the crystal-pore wall interface. At first order, we can assume that ∂U s,M ∂ϕs = ∂U ∂ϕs , thus the use of equation (7.12) allows to obtain:

p * s = p l + ω p (8.53)

Final remark

Microfluidics offers an interesting way to study experimentally the mechanisms of confined crystallization (see for example [START_REF] Scherer | Pressure form crystallization in pore channels[END_REF]). These experiments allow to study single crystals in microfluidics channels simulating the in-pore crystallization, and giving interesting insights in the nucleation/growth processes of confined phase transition. Moreover, the polymer used in the fabrication of these microchannels (PDMS -Polydimethylsiloxane) is deformable and birefringent, which may allow to determine the intensity of the crystallization pressure by the measure of the deformation of the channel subjected to crystallization [START_REF] Song | Elastomer based tunable optofluidic devices[END_REF]. These experiments could unfortunately not be developed in the frame of this PhD work but several teams are currently working on this kind of projects, hopefully leading to interesting results in the field of confined crystallization.

Chapter 9

Estimation of crystallization pressure in the case of CCS In-pore crystallization is the source of the creation of stresses on the porous medium. The knowledge of the crystallization pressure during CCS is then mandatory to determine the stress state of the aquifer. In the following chapter, we study a REV subjected to carbon dioxide injection and thus to salt precipitation. The first section presents two macroscopic modellings and the associated poromechanical calculations. The second section presents a microscopic calculation of the nucleation and growth of crystals followed by an upscaling in order to obtain the poromechanical behavior on a REV. As we have defined in section 3.4, residual water is divided into three different forms: wetting water, blob water and capillary trapped water. As we discuss in the following, the first two modellings apply to blob water, while the microscopic modelling applies to capillary trapped water. The complete behavior of the REV should then be a combination of both effects.

Macroscopic behavior: poromechanics at the

REV scale

Modelling at constant concentration

Modelling of crystal invasion at equilibrium

Let us consider a REV close to the injection well. Considering the scenario described in section 6.3, this REV will be subjected at first to the multiphase flow of carbon dioxide displacing brine, then to the evaporation process leading to the precipitation of salt and the creation of the crystallization pressure, the temperature remaining constant (see section 5.4). Let us model the porous medium as in figure 9.1. The porous medium is composed of spherical pores of radius r p linked together by small cylindrical channels of radius r s . We will consider liquid-like crystals (see section 8.1.4).

We will also consider that the ratio between the radius of a pore and the adjacent channel is constant and equal to a parameter λ (shape ratio of the porous medium):

λ = r s1 r p1 = . . . = r sx r px (9.1)
Supersaturation will increase with evaporation until it reaches the nucleation threshold, i.e. when the nucleation rate becomes non negligible, as described in section 4.2.3. In the model we will consider a continuous crystalline phase. This means that all new crystals created through nucleation dissolve or coalesce within this single crystalline phase (Ostwald ripening). We have then only one crystal growing in the porous medium, invading 9.1. MACROSCOPIC BEHAVIOR: POROMECHANICS AT THE REV SCALE progressively the pores of decreasing size. The initial pore where nucleation occurs can be determined by the knowledge of both the pore size distribution and the brine saturation at the beginning of the evaporation. Considering that crystals can only precipitate where brine is present, it means that the first pore subjected to nucleation and crystal growth will be the biggest pore filled with brine. Indeed, because of the stability of crystals as defined in section 8.1.2, crystals will try to have the biggest curvature possible and thus to fill the biggest pores they can reach. The knowledge of the radius of the pore in which precipitation occurs first can be then obtained from the Mercury Intrusion Porosimetry (see section 12.3) of the considered rock. Indeed, as described in appendix E the maximal pore radius filled with brine can be determined as a function of the water saturation. As a result, if we know the saturation S prec l at which precipitation occurs, we know the radius r p1 in which this precipitation occurs. There is then a link between the nucleation supersaturation S prec , the pore in which it occurs (r p1 ) and S prec l . Figure 9.1: Modelling of the porous medium as a succession of spherical pores of decreasing size linked by cylindrical channels of decreasing radius.

Evaporation and precipitation

Let us consider the porous medium at the end of the displacement regime. We will consider that after the percolation, the REV is still saturated with S 0 l =0.45 of brine. Consequently, there is an initial quantity n 0 of salt in the REV corresponding to the quantity of the limiting ion α (the ion belonging to the considered crystal composition having the smallest concentration in the initial solution):

n 0 = t α ν α M s ρ 0 w S 0 l φdΩ (9.2)
with t α , the mass proportion (g/kg) of the limiting ion, ν α its stoichiometric coefficient, M s the molar mass of the considered salt, ρ 0 w the initial density of the brine, φ the porosity and dΩ the volume of the REV considered.

As expressed in section 3.3.2 capillary currents can only occur if the brine mobility is high enough. As we consider a REV at its residual saturation and that evaporation begins only at this moment, we can neglect the effect of capillary current and consider that the concentration of the solution at the beginning of evaporation is the initial concentration of the aquifer. Supersaturation of the solution can then be written from equation (4.23):

S = i a ν i i a νw w K s (p l , T ) (9.3a) S = γ N +-m N α β N K s (p l , T ) a νw w (9.3b) with m i the molality of ion i, N = i ν i , γ N +-= i γ ν i i , β N = i β ν i i , (β i = m i /m α
with m α the molality of the limiting ion, β α = 1).

Writing that m α = ν α n 0 /(V w × ρ w ) we have:

S = γ N +- tα ναMs S 0 l N β N S N l K s (p l , T ) ρ 0 w ρ w N a νw w (9.4)
where S l is the current brine saturation.

Precipitation begins when the supersaturation has reached the nucleation threshold S prec . We can then calculate the value of the brine saturation at which the precipitation begins:

S prec l = γ +- tα ναMs S 0 l β (S prec K s (p l , T )) 1/N ρ w ρ 0 w a νw/N w (9.5)
For a simple salt as halite NaCl or anhydrite CaSO 4 (see chapter 4.2), the equation reduces to:

S prec l = γ +- tα Ms S 0 l β S prec K s (p l , T ) ρ w ρ 0 w (9.6)
This equation is more complex than it looks like. Indeed, because of capillarity (p cap = f (S l )), K s (p l , T ) is a function of S l , as is S prec . Indeed, the crystal must adopt a curvature corresponding at least to the radius of the pore in which it precipitates (i.e. RT υs ln S prec = 2σ r p1 -δ ). The smaller the pore in which precipitation occurs (i.e. the smaller S prec l is), the higher the supersaturation has to be for the nucleation to occur. Resolution of equation 9.1. MACROSCOPIC BEHAVIOR: POROMECHANICS AT THE REV SCALE (9.6) has been proceeded with Maple software.

In order to simplify the model, we will in the following, consider that the concentration of salt in the brine is constant after the beginning of the precipitation. It is similar to say that precipitation rate matches the concentrating of the solution by evaporation. In reality, precipitation should be slower than evaporation, leading to an increasing concentration and thus leading to higher crystallization pressure through the increase of the supersaturation. As a conclusion, values of stress calculated under constant concentration are a lower bound of the real stresses.

The constant concentration hypothesis leads to

m α = ν α n 0 S l φdΩρ w - n cr S l φdΩρ w = m α (S prec l ) (9.7)
with m α (S prec l ), the molality of the salt at the moment of the precipitation and n cr the mole quantity of crystal precipitated. This condition leads to the following crystal saturation S s :

S s = ρ 0 w υ s t α ν α M s S 0 l 1 - ρ w ρ prec w S l S prec l (9.8)
where ρ prec w is the density of water at the beginning of precipitation.

We can finally neglect the variation of water density with the concentration leading to:

S s = ρ w υ s t α αM s S 0 l 1 - S l S prec l (9.9)

Invasion of the crystal phase and crystallization pressure

After the first nucleation, the crystal phase grows in the pore of radius r p1 at a rate following the previous relation (9.9). As soon as the crystal reaches the pore wall, as presented in section 8.3.2, crystallization pressure is created. The volume of crystal needed to reach the pore wall is small because we are considering precipitation in a unique pore. We can thus assume that confinement happens almost simultaneously with the beginning of precipitation. As a result, the crystallization pressure applied on the pore wall of the spherical pore follows equation (8.37): As presented in section 8.3.2, as soon as the supersaturation is high enough, the crystal penetrates the entry channel and creates another crystallization pressure which corresponds to the cylindrical pore case :

ω p = RT υ s ln S(p l , T ) - 2σ r p1 -δ (9.
ω s = RT υ s ln S(p l , T ) - σ r s1 -δ (9.12a
)

ω s = RT υ s ln S(p l , T ) - σ λr p1 -δ (9.12b)
using the definition of the shape ratio λ (equation (9.1)). The supersaturation at which the crystal begins to invade the entry channel is:

RT υ s ln S rs = σ λr p1 -δ (9.13)
Finally, if the crystal continues to fill the pore volume, it will reach the next pore r p2 . Once this second pore is filled, the growing crystal will create also a crystallization pressure ω p2 = RT υs ln S(p l , T )-2σ r p2 -δ . These subsequent crystallization pressures occurring after the creation of the crystallization pressure on the pore of radius r p1 are of smaller magnitude than the first one because of the higher value of the size term. Moreover, as we have seen in section 7.3 with equation (7.23), the influence of the stress depends on the saturation of the considered phase. As the quantity of crystal penetrating the side channels and the subsequent pores is much lower than the quantity of crystal in pore of radius r p1 we will neglect the influence of these crystallization pressures and consider that the whole crystal phase applies a crystallization pressure as in equation (9.11) on the rock matrix.

In order to determine the value of the crystallization pressure, we need now to express the supersaturation after the beginning of the precipitation. Using the constant concentration hypothesis (Q R = cste), the supersaturation can be expressed with the current ion activity product and the equilibrium constant

S S prec = Q R Ks(p l ,T ) × Ks(p prec l ,T ) Q R = Ks(p prec l ,T )
Ks(p l ,T ) , with K s (p prec l , T ) the equilibrium constant at the beginning of the precipitation, p prec l being the brine pressure at the moment of the precipitation. We can then write 9.1: Relative chemical composition of the brine from the Dogger aquifer in the region of Fontainebleau (in g/kg) from [START_REF] André | Numerical modeling of fluid rock chemical interactions at the supercritical CO2-liquid interface during CO2 injection into a carbonate reservoir, the Dogger aquifer (Paris Basin, France)[END_REF].

ω p = RT υ s ln K s (p prec l , T ) K s (p l , T ) (9.14)
Crystallization pressure depends on the crystal through the molar volume and its solubility (K s ), but it also depends on the characteristics of the medium through the evolution of p l and p prec l . Evolution of crystallization pressure after the beginning of the precipitation will then depend on the evolution of capillary pressure with evaporation.

Poromechanical calculation

The evolution of crystallization pressure during the evaporation depends only on the evolution of the equilibrium constant. This equilibrium constant depends on the brine pressure and on the temperature. Brine pressure will evolve through capillary pressure as defined in section 3.3.2. In the following we will use the capillary pressure curve from the MDV and BBL rocks (see chapter 12) as calculated in appendix D. In the calculation we will consider that the carbon dioxide pressure is constant, neglecting the pressure drop on the REV. We will take P CO2 = 22MPa which is a mean value for the different simulations in the literature (see for example [4]). Finally, capillary pressure depends on the different saturation degrees of gas and liquid phase. As crystals are growing in the liquid phase, the gas saturation becomes S g = 1 -(S s + S l ). As a result, capillary pressure will be calculated as

p cap = f (1 -S g ).
The brine concentration will be that of the Dogger aquifer of the Paris basin (5g/kg) and the different proportions for the different ions are presented in the table 9.1.

As presented in the section 8.3.1 the stress transmitted to the pore wall is the addition of the crystallization pressure ω p and the brine pressure. In figure 9.2 is presented the different contribution of the transmitted stress (p * s ) to the pore wall at 40 ℃ for halite for the MDV capillary curve. Because of the decrease of the liquid pressure, the crystallization pressure increases through the diminution of the equilibrium constant which leads to a further precipitation. The addition of this crystallization pressure and the liquid OF CCS pressure is finally decreasing. As we know from thermodynamics, liquid pressure is not the only parameter to affect the crystallization pressure. Temperature also modifies the equilibrium through the equilibrium constant. As we have seen in the section 5.4, close to the injection well, the temperature of the aquifer is set mostly by the temperature of the injected carbon dioxide.

Figure 9.3 shows the influence of temperature on the transmitted stress for halite/hydrohalite (NaCl/NaCl • H 2 O) and for the couple anhydrite/gypsum (CaSO 4 /CaSO 4 • 2 H 2 O) for the two rocks considered. The thermodynamics of calcium sulfate shows indeed, that below 40 ℃, gypsum is the stable form, while above this temperature, anhydrite is the stable form. Both forms have been represented at 40 ℃ to show the incertitude on the existing form. For halite, the stable form below 0 ℃ is hydrohalite

(NaCl • H 2 O).
In these figures we can see the strong influence of temperature on the transmitted stress. First of all, it plays a role on the moment when the precipitation begins by modifying the equilibrium constant and the solubility. It is interesting to notice that, if for halite an increase in temperature increases the solubility and thus delays the precipitation, it is the opposite for anhydrite: solubility decreases with temperature. Another remark is that the solubility of anhydrite/gypsum is much lower than the solubility of halite/hydrohalite. Indeed, anhydrite/gypsum begins to create a stress at high water saturation (between 0.15 and 0.20) when halite/hydrohalite needs a much stronger evaporation to precipitate (S prec l ≈ 0.003). Finally, because of the strong difference in capillary pressure, the values of the stress change strongly between the MDV and BBL rock (one order of magnitude).

If we consider the derivative of the transmitted stress with temperature we have:

∂p * s ∂T = 1 υ s ∆ r S 0 (p l , T ) -∆ r S 0 (p 0 l , T ) (9.15)
The effect of temperature is then completely grasped by the entropy of reaction. The distance between two curves on figure 9.3 depends then on the effect of liquid pressure on entropy. Molar entropy of reaction's variations with liquid pressure can be accurately approximated with a linear fit. As a result, variation of stress with temperature writes:

∂p * s ∂T = p l -p 0 l υ s ∂∆ r S 0 (T ) ∂p l (9.16)
As a conclusion, the relative position of two curves depends on the sign and the value of the derivative of the entropy with liquid pressure at the temperatures considered and the distance between them increases with capillary pressure.

Another effect of the temperature is to modify the slope of the evolution of the stress. Indeed, let us consider the derivative of the transmitted stress with water saturation:

∂p * s ∂S l = -1 + ∆ r υ 0 (T ) υ s ∂p cap ∂S l (9.17)
with ∂p l ∂S l = -∂pcap ∂S l . As capillary pressure is monotonously decreasing with increasing water saturation, , the evolution of the transmitted local stress will depend on the sign of 1 + ∆rυ 0 (T ) υs and thus to the value of the molar volume of reaction. For halite and gypsum, the local stress is always decreasing , while for anhydrite the stress is always increasing, with decreasing water saturation. The slope is determined by the evolution of the molar volume of reaction with T . For gypsum, the decrease is much stronger than for halite.

The poromechanical calculation presented in section 7.3 gives us the equivalent tensile stress on the REV for crystallization processes:

̟ = b 2 + K + (4/3)µ N (S l + S s ) p l -p 0 l + S s ω p (9.18)
with p 0 l and p 0 cap the liquid pressure and the capillary pressure at the beginning of evaporation. Considering that evaporation begins at S 0 l = 0.45, the corresponding capillary 9.1. MACROSCOPIC BEHAVIOR: POROMECHANICS AT THE REV SCALE pressure is 0.01879 MPa. As a result p 0 l = 21.98MPa.

Equation (9.18) shows the evolution of the equivalent tensile stress for the considered REV following the precipitation. The first term in the bracket accounts for the evolution of the liquid pressure and the second term accounts for the evolution of the crystallization pressure. An important point is that the term in factor of the liquid pressure is S l + S s . Indeed, as the crystal grows in the solution, it pushes brine in pores where it had been displaced by the CO 2 increasing then the surface of contact between brine and the porous matrix. In other words, crystallization increases the impact of capillarity. The second term is the effect of crystallization pressure and depends only on the quantity of salt precipitated.

In figure 9.4 are plotted the equivalent tensile stress on the REV without crystallization i.e. with pure water and thus only capillary forces for MDV and BBL. We see at first a decrease in stress because it follows the decrease of liquid pressure. Then at the end the stress increases toward zero because of the decrease in water saturation which decreases in absolute value the equivalent tensile stress.

In figure 9.5 is plotted the evolution of the difference between the equivalent tensile stress for halite and anhydrite with the equivalent tensile stress without any crystallization (plotted in figure 9.4) for BBL and MDV. Behavior of the stress results from a competition between the expansion effect which increases the compressive effect of capillary pressure and the crystallization pressure which opposes to this effect. For halite and gypsum, the crystallization pressure magnitude is too small to oppose to the expansion effect. As a result and paradoxically, crystallization increases the magnitude of the compressive effect of the negative liquid pressure. Behavior of anhydrite is less straightforward, because of the strong influence of liquid pressure on the equilibrium constant. In anhydrite equivalent tensile stress difference, we clearly see the influence of capillary pressure curve on the behavior in figures 9.5(d) and 9.5(b).

Discussion and conclusion

The model presented here allows to calculate the intensity of the crystallization pressure within an REV close to the injection well. Thanks to this simple modelling, we have been able to highlight the importance of capillary pressure and pore network in the creation of crystallization pressure. The constant concentration hypothesis is interesting by giving a lower bound of the crystallization pressure and thus of the obtained stress. Indeed, precipitation should be slower than evaporation leading to concentrations and thus supersaturations higher than those calculated here. As a conclusion, the values obtained here clearly show an almost non significative effect of crystallization pressure on the pore OF CCS structure. This comes from the constant concentration hypothesis but also on the pore structure considered. Indeed, the rocks considered here have large pores which leads to a quite small supersaturation at the precipitation. A rock with a finer structure would lead to higher crystallization pressures. Note also that this simulation is based on an equilibrium hypothesis considering a single crystal.

The very low temperatures presented here (0 ℃ and -10 ℃) are not unrealistic. Indeed, some studies (see for example [START_REF] Vilarrasa | Liquid CO2 injection for geological storage in deep saline aquifers[END_REF]) suggest to use liquid carbon dioxide instead of supercritical leading to very low injection temperatures mainly because of Joule-Thomson effect. It is finally interesting to underline the potentially dramatic effect of the anhydrite/gypsum transformation if the temperature is close to the transition temperature and one form becomes unstable. The behavior will then be close to the behavior of mirabilite/thenardite [START_REF] Espinosa-Marzal | Crystallization of sodium sulfate salts in limestones[END_REF][START_REF] Shahidzadeh-Bonn | Damage in porous media due to salt crystallization[END_REF][START_REF] Steiger | Crystallization of sodium sulfate phases in porous materials: The phase diagram of Na2SO4-H2O and the generation of stress[END_REF]: a change in temperature will make one of the form unstable and will force its dissolution. The reprecipitation of the other form will then occur at a high supersaturation (and likely at a smaller brine saturation), leading to high transient crystallization pressure which cannot be grasped with the equilibrium modelling.

Simulation of a REV submitted to carbon dioxide evaporation

In order to determine correctly the moment of evaporation and to go beyond the previous simulation we need to determine the complete evaporation/precipitation process. To do In this section we will keep the classic hypothesis of poromechanics (REV composed by the solid matrix m, the porous network filled with supercritical CO 2 -rich phase also called rich-gas phase (index G), salt crystals (s) and pore aqueous solution (l)). We will keep the same notations as in chapter 8, i.e. the index i for the ions in solution, w for liquid water, v for water vapor (in phase G) and CO 2 for the supercritical carbon dioxide. We will consider here only halite as the possible salt to precipitate (S = NaCl, i = Na + , Cl -). Following the scenario described in section 6.3, the evaporation regime and thus the precipitation begins after the end of the displacement regime i.e. when the relative permeability of water has reached zero. We consider that the porosity still occupied by the brine is 45% of the initial porosity at the beginning of the evaporation.

Brine evaporation

As recalled before, the injected carbon dioxide is dry at the entrance of the aquifer. Evaporation then occurs at the interfaces between carbon dioxide and the remaining brine. The liquid-vapor equilibrium for water is pressure dependent and can be estimated with Spycher's work as presented in section 4.1.3.2.

The molar fraction of vapor in the supercritical phase is then obtained by the following relation:

y = K 0 w β v p G exp υ w RT (p l -p ref ) (9.19)
with p ref the reference pressure taken as 1 bar, K 0 w is the equilibrium constant at p l = p ref and T , υ w is the molar volume of water taken as constant in the range of pressure considered. We will consider as usual that this equilibrium is reached instantaneously (evaporation kinetics is instantaneous), which means that the decrease of brine saturation will only depend on the CO 2 flow rate.

The mass balance for water can then be written using the relation (3.8):

∂m H 2 O ∂t = ∇ • ρ G D∇ ρ v ρ G -ρ v V G (9.20)
where m H 2 O = φ G ρ v + φ l ρ w is the total mass of water per unit of initial volume, ρ w is the liquid brine density taken as constant and equal to 1 kg/L, ρ G and ρ v are the density of the rich-gas phase and the mass of vapor per unit of rich-gas phase volume. They are linked to the molar fraction of vapor in the carbon dioxide phase as:

9.1. MACROSCOPIC BEHAVIOR: POROMECHANICS AT THE REV SCALE ρ G = 1 υ G M CO 2 (1 -y) + M H 2 O y ; ρ v = M H 2 O υ G y (9.21)
where υ G is the molar volume of the rich-gas phase while M H 2 O and M CO 2 are respectively the molar mass of water and carbon dioxide. Following [START_REF] Spycher | CO2-H2O mixtures in the geological sequestration of CO2. I. assessment and calculation of mutual solubilities from 12 to 100C and up to 600 bar[END_REF], we can neglect the influence of water on the equation of state of G. Then υ G can be calculated by the Redlich-Kwong's model as presented in appendix A.2.

Finally, D stands for the effective diffusion coefficient of the vapor within G, and V G is the advective velocity which is due to the pressure gradient of G and is given by a generalized Darcy's law as equation (3.10):

V G = - k 0 k G r (S G ) η G ∇p G (9.22)
where k 0 is the intrinsic permeability, k G r the gas relative permeability as defined in section 3.2.2, η G is the dynamic viscosity of G.

Salt crystallization and crystallization pressure

The macroscopic behavior of evaporation in a REV is well known. The difficulty in modelling the behavior is to model the crystallization process. Indeed, as expressed in paragraph 8.3, crystallization will proceed both by nucleation and crystal growth. However, if we continue to keep the equilibrium hypothesis as in the previous section (and the single crystal hypothesis), we can obtain a consistent thermo and poromechanical model. We consider the ideal behavior of a 1:1 electrolyte (like halite or anhydrite) in stoichiometric proportions. The temporal derivation of the supersaturation S is then:

∂S ∂t = - 2S φ l ∂φ l ∂t + ρ s ρ i ∂φ s ∂t - S K s ∂K s ∂t (9.23)
with φ l the porosity filled by the liquid, φ s the porosity filled by the crystals, K s is as usual the equilibrium constant, ρ i is the density of the ion i and ρ s is the density of the considered salt.

The first term of rhs of equation (9.23) accounts for the supersaturation evolution due to the evaporation process, the second one describes the species consumption by the crystallization and the last one describes the evolution of the equilibrium constant with the liquid pressure (i.e. through capillary effects). Following the invasion process described in the previous section we will consider that the crystal consists of one unique phase which fills the biggest pores possible. Indeed, as OF CCS water is the wetting fluid, the capillary effects will reduce the carbon dioxide presence only in the biggest pores. As a conclusion, considering that crystals can only grow in brine filled pores, the radius of the pores filled by crystals will be the same as the biggest radius of the pores filled by brine. We suppose then that at equilibrium we have r s = r h , with r s the curvature radius of the brine/crystal interface, and r h the radius of the CO 2 /brine interface. As a conclusion, using the definition of capillary pressure and equation (8.9), we can write;

S eq = exp συ s (p G -p l ) σ lG RT (9.24)
If the supersaturation is bigger than the equilibrium supersaturation of equation (9.24), the crystal will grow in order to consume this excess supersaturation and will penetrate new pores. On the contrary, if the supersaturation is smaller than the equilibrium value, the crystal will shrink. This behavior obviously supposes the presence of a crystal. We then have the following conditions:

m s ≥ 0 ; S ≤ S eq ; m s (S -S eq ) = 0 (9.25)
To describe the crystallization pressure, we will once again consider the classic equation (8.37):

ω p = RT υ s ln S eq -σκ s (9.26)
with κ s the curvature of the side of the crystal interacting with the pore wall. As expressed in the previous sections, crystallization pressure depends on the curvature of the pore in which the crystal is precipitating. Considering a distribution of crystals in different pores, and supposing as usual that the supersaturation is constant in the considered REV, we can average the curvature term in the crystallization pressure equation using the parameter c = rs rm , where r s is as defined previously the curvature between the crystal and the brine and r m , the averaged radius of the pores in which the crystal are precipitating. This parameter should be a function of S l and S s , however, we will consider in first approximation that is is constant. The crystallization pressure equation can then finally be written:

ω p = (1 -c)
RT υ s ln S eq (9.27)

9.1. MACROSCOPIC BEHAVIOR: POROMECHANICS AT THE REV SCALE

Model simulation

In order to simulate the drying-out of the REV we consider a one-dimensional cubic structure made of an isotropic medium of length L = 1m. It is supposed insulated on its lateral surfaces and is submitted to a constant flow of carbon dioxide along the x-axis. Transport of fluid and species occurs only in the x-direction.

The initial state is just after the end of the transport process, when brine has reached its residual saturation. As a consequence, brine relative permeability is supposed constant and equal to zero whereas the relative permeability of CO 2 is assumed to be equal to 1 D (≈ 1.10 -12 m 2 ). Because of the high permeability of the medium, we are able to neglect the pressure gradient of CO 2 (it is equal to 10kPa/m considering a density of 600 kg/m -3 , a dynamic viscosity of 4.10 -5 Pa.s and an advective flow of 1 kg/s/m -2 ). The initial condition for the carbon dioxide phase will then be isobaric and the resulting initial hydraulic conditions are:

S l (x, t = 0) = S 0 l ; p 0 l (x, t = 0) = p 0 l = p 0 G -f (S 0 l ) (9.28)
where S 0 l is the residual saturation, f (S l ) being the capillary pressure function. Finally, we assume that there is no flow at the boundary of the medium for brine and salts. Evaporation will proceed until a minimal liquid saturation is reached. This minimal saturation corresponds to the quantity of water adsorbed on the porous matrix and the crystals and is considered as non evaporable.

The system of equation to solve is then:

B ∂p l ∂t = ∇ • ρ v D υ w RT ∇p l -V G (9.29)
where:

B = ρ v φ G υ w RT + m i 2ρ s 1 K s dK s dp l - σ sl σ lG υ s RT + φ 0 ρ w dS l dp l if S = S eq (9.30) B = ρ v φ G υ w RT + φ 0 ρ w dS l dp l if m s = 0 and S < S eq (9.31)
Finally, as in the previous section, the stress is obtained by:

̟ = b 2 + K + (4/3)µ N (S l + S s ) p l -p 0 l + S s ω p (9.32) OF CCS (a) (b) (c) (d)

Results and discussion

The input data for the numerical resolution are reported in table 9.3. The effective diffusion coefficient of vapor in supercritical CO 2 is supposed to be equal to the diffusion of water within liquid CO 2 i.e. D = 1.10 -7 m 2 /s.

Calculations where made for a velocity of 5.10 -3 m.s -1 which approximatively refers to a carbon dioxide flow rate of 0.1MT per year through a surface of 1m 2 of porous medium. Figure 9.6 show the profiles of liquid pressure, equivalent tensile stress, liquid saturation and volume ratio of salt crystal for the two retention curves described in Appendix D. The rocks considered are described in chapter 12 and are the Mondeville (MDV) and Bois-Brûlé (BBL) cores. In order to ease the comparison, the profiles are not drawn at the same time, but at the same length of penetration of the drying-out front. The difference in this penetration is a consequence of the boundary condition at x = 0, which imposes for MDV a drier CO 2 income than for the BBL (p l is lower at S min l for MDV than for The first significant difference between the two simulations is the amount of salt that precipitates. Indeed, the capillary pressure curve influences strongly this amount of crystal. This is not surprising as the saturation profiles are the same, while the equilibrium supersaturation which increases with decreasing p l is much higher for the MDV sample. MDV uses less salt than BBL to reach the equilibrium supersaturation.

Let us now consider the tensile stress undergone by a REV of rock as the drying-out progresses within it. For the two simulations, it first decreases as expected by the equivalent compressive stress of the decreasing liquid pressure. As for the previous modelling when the amount and the equivalent pressure of salt becomes sufficient to counterbalance this effect, ̟ rises rapidly and reaches a maximal value of 0.3MPa for the BBL and a significantly higher value of 1.9MPa for the MDV rock, then for this latter a small decrease of 0.1MPa is observed, corresponding to the variation of φ G due to the salt precipitation which locally increases the mass concentration of vapor ρ G . This variation induces a new equilibrium that forces, as shown in figure 9.6(d) the dissolution of a small amount of salt crystal, and thus a small reduction of ̟. This higher value of ̟ for MDV can be explained by the higher supersaturation at equilibrium leading to a higher crystallization pressure, while the reduction of the volume of salt that precipitates is tempered by the reduction of the solubility constant caused by the extremely low in-pore brine pressure.

In conclusion, and as underlined by the previous model, the capillary pressure curve and thus the pore size distribution of the rock appears to have a great impact on the two main consequences of the drying-out which are the amount of salt S s formed responsible of the clogging of the porosity and the equivalent tensile stress ̟ which can lead to a damage of the host rock: a finer porosity will induce a higher ̟ but a lower S s .

For both simulations, a sharp drying front is observed, which is not surprising as a high advective velocity of CO 2 is chosen for these simulations. Indeed, the characteristic time of vapor transport by advection through a cell of length l = 0.01m if

τ adv = φ G l/V G ≈1s, while the characteristic time for the diffusion process is τ dif f = (φ G l) 2 /D ≈500s.
The vapor diffusion is much longer than the advection and is then overshadowed by this latter. In order to observe the effect of diffusion during this evaporation, the same simulation is proceeded with a velocity of 5 µm/s. As expected, the result in figure 9.6 underline that due to the diffusion process, the profiles are smoother. A noticeable difference in S s and ̟ can also be observed. Indeed, the constant income of water by diffusion tends to keep a higher saturation degree that leads to a larger amount of salt that precipitates and a lower equivalent tensile stress. This effect however vanishes near x = 0 because of the boundary condition. To summarize, the flow rate appears to significantly impact the shape of the ̟ and S s profiles: a lower flow rate leads to a lower 9.1. MACROSCOPIC BEHAVIOR: POROMECHANICS AT THE REV SCALE equivalent tensile stress. This conclusion comes however from by the assumption of a discontinuous water phase at the beginning of the evaporation. In consequence, capillary currents and matter transport through the liquid phase are prevented.

Comparison of the two modellings

The two models presented here are based on different hypotheses. The first one neglects all kinetics and considers that the concentration in the REV is constant. It also presents a invasion model based on a single crystal growth. The second model takes into account the evaporation kinetics but not the precipitation kinetics and takes into account transport phenomena of water and carbon dioxide.

In figure 9.7 are presented the equivalent tensile stress obtained for these models. For the second model, the plot has been made at a distance of 0.6m (same time as for figure 9.6) from the entrance of the simulated volume. The behavior of the two plots are similar: at first there is a decrease in the macroscopic stress because of the decrease of liquid pressure. Then the equivalent tensile stress departs from the behavior of the pure liquid as crystallization pressure begins to apply stress on the REV. However, on figure 9.7(a) precipitation decreases the equivalent tensile stress while for 9.7(b) crystallization increases the equivalent tensile stress. Moreover, there are two orders of magnitude between the stress of the two simulations. This is due to the fact that the concentration of the solution in the second simulations is taken as S=1. There is then much more salt in the REV than for the first simulation which modifies strongly the behavior. This also explains why the precipitation occurs much earlier in the second simulation.

Another strong difference between the first and the second macroscopic model is the calculation of the crystallization pressure and the supersaturation. In the first model we have considered a constant concentration which allowed to know the supersaturation at every water saturation. In the second model, the condition has been replaced with the condition r s = r h considering that the crystal penetrates the pores as fast as the carbon dioxide.

Finally, it is important to note that the crystallization pressure described here can reach values as high as 50-60MPa. However, we have shown in section 8.2 that the crystallization pressure has an upper bound corresponding to the moment when the liquid film vanishes and the two solid surfaces make contact. Accordingly, depending on the kind of reservoir rock and on the salt considered, it is likely that in some condition, the real stress would be lower than this calculated here as all stresses above the critical stress vanish. OF CCS

Microscopic behavior: nucleation and stress creation at the pore level

The previous models have highlighted the importance of the crystallization pressure and have allowed to estimate its value in macroscopic cases. However, the strong drawbacks of these models is the equilibrium hypothesis. Indeed, for both models the crystals were at equilibrium with the solution and no kinetics of nucleation and growth were considered. This modelling is well fitted for the blob water where evaporation is slow because of the small interfacial area between carbon dioxide and brine. However, for capillary trapped brine which undergoes a strong evaporation, the equilibrium hypothesis is no longer valid. In the following section, we present a microscopic model [START_REF] Osselin | Poromechanics of salt nucleation within an unsaturated reservoir rock[END_REF] allowing to take into account the nucleation and crystal growth on the estimation of crystallization pressure. The numerical calculation is nevertheless very time consuming and makes strong assumption on the pore geometry, which forbids a macroscopical use as in the previous calculation. However, an uspcaling considering a homogeneous behavior of the REV is presented at the end of the model.

Following the hypothesis of the previous calculation, we suppose that the brine has reached its residual saturation and presents then a discontinuous pattern. In order to model the microscopic behavior, we need to model the pore geometry. To do so, we consider the work of Tuller and Or [START_REF] Tuller | Unsaturated hydraulic conductivity of structured porous media: A review of liquid configuration-based models[END_REF][START_REF] Tuller | Hydraulic conductivity of variably saturated porous media: Film and corner flow in angular pore space[END_REF][START_REF] Or | Hydraulic conductivity of partially saturated fractured porous media: flow in a cross-section[END_REF][START_REF] Or | Flow in unsaturated fractured porous media: Hydraulic conductivity of rough surfaces[END_REF] which have underlined the drawbacks of the cylindrical/spherical pore model. Indeed, in cylindrical pores, the saturation state is either full of empty, not allowing any bi-occupancy. However, as carbon dioxide flows in the biggest pores, evaporation will mainly occur in the corners of the big pores as shown in section 3.4. Let us now consider a corner of a big pore in which carbon dioxide is flowing as described in figure 9.8. We will call pore size the side length of the corner L. Capillarity will be taken into account by considering the meniscus of radius r on the figure (we will suppose that brine wets perfectly the pore wall). Evaporation of the residual brine in this corner will lead as usual to an increase in concentration of the solution and thus to a nucleation and a crystal growth. In this calculation capillarity effect on the equilibrium constant has been neglected.

Evolution of the salt quantity in the corner

Let us consider an ideal solution of a 1:1 electrolyte (NaCl, CaSO 4 ). Supersaturation is defined as usual as:

S = n 2 s S 2 l v 2 pore K s (9.33)
with n s the number of mole of salt in solution (considered in stoichiometric proportions), v pore the volume of the pore and K s the equilibrium constant. The evolution of dissolved salt quantity in the corner will then be determined by nucleation and crystal growth of the crystals present in the solution.

Kinetics of crystal growth for small crystals

In section 4.2.3.1 we have recalled the macroscopic equations (4.54a) and (4.54b) for the growth/dissolution of crystals in solution. The discriminating value between growth or dissolution (considering that nucleation has occured) is S = 1. As we have seen in chapter 8, equilibrium supersaturation of a small crystal in solution is bigger than unity.

Let us now consider a small spherical crystal of radius r. Ostwald-Freundlich equation (8.10) gives us the equilibrium supersaturation for this crystal:

S * = exp 2συ s RT r (9.34)
If the supersaturation is bigger than this value, the crystal will grow. On the contrary, if the supersaturation is smaller than this value the crystal will shrink and disappear. We

v = kA S S * θ -1 η if S S * > 1 (9.35a) v = kA 1 - S S * θ η ′ if S S * < 1 (9.35b) (9.35c)
considering that nucleation has already occurred. For simple salts like KNO 3 , NaCl, CaSO 4 . . . we can consider that θ, η and η ′ are equal to unity. As a result, if nucleation has occurred, the previous equations can be written in a single form:

v = kA S S * -1 (9.36)
As we can see, if S = S * , there is no growth. If S < S * , v is negative and the crystal shrinks, and finally, if S > S * , the crystal grows. If we use small spherical crystals, we can express v and A with the radius r of the crystal:

v = dn s dt = d dt 4 3 πr 3 υ s = 4πr 2 k S S * -1 (9.37)
This leads to the following differential equation on r:

dr dt = υ s k S S * -1 (9.38)
with S * = exp rc r , where r c = 2συs RT . If we consider a constant supersaturation (neglecting the consumption of salt by the growing crystal, i.e. small crystal in a big bulk), we can see in figure 9.9 that once the crystal has begun to go one way or another it will not stop. Indeed, if the supersaturation is high enough to make the crystal grow, the crystal radius will grow which will lower its equilibrium supersaturation S * and make it grow more and more following equation 9.38 (the same process happens when the crystal is shrinking). As a consequence, the equilibrium radius given by equation (8.10) is an unstable equilibrium (this is the same radius and unstable equilibrium as the critical radius for nucleation, see section 4.2.3.1). OF CCS Figure 9.9: Speed of growth for a small crystal as a fonction of the radius for constant supersaturation.

Global differential equation on the salt quantity

Crystal growth equation has been obtain in the previous section. The rate of nucleation for an homogeneous nucleation is:

I = k B T πα 5/3 η w exp - 256σ 3 α 2 27 (k B T ) 3 (ln S) 2 (9.39)
with α the volume of a formula unit (i.e. υ s /N A with N a the Avogadro Number). We consider that the nucleii created have a radius corresponding to the critical radius of nucleation as presented in section 4.2.3.1:

r cr = 2υ s σ RT ln S (9.40)
Once created, the crystal will evolve considering its size and the current supersaturation as expressed by equation (9.36). As a consequence, the evolution of the crystal quantity in solution containing N nuclei is then: RT . The first term accounts for the growth/dissolution of all the existing crystals, while the second term accounts for the addition of new crystals by nucleation. This equation is particularly complicated to solve because of the temporal dependence on the history of each crystal which will evolve accordingly to its size. However, we can compare both terms of the evolution. It appears that the second term plays a role only when the supersaturation is high enough because of the nucleation threshold (see section 4.2.3.1). As we will see with the simulations, the homogeneous nucleation threshold is around 1.9 for halite.

dn s dt = N n=1 k ′ r 2 (n, t) S(t) exp - r c r(n, t) -1 + I(S, t)f (S l ) 4π 3υ s
Let us now consider the temporal evolution of the supersaturation. From equation (9.33), we obtain:

dS dt = 2n s S 2 l v 2 pore K s dns dt + 2n 2 s K s d dt 1 S 2 l v pore = • S cr + • S ev (9.42)
where

• S cr and

• S ev respectively stands for the crystal growth/nucleation and the evaporation contributions to the variation of S. Using the crystal growth expression (9.37) we can define:

Γ = 2S Vw • ξ 4N kπr 2 s S S * -1 = A    • ξ N r 2 s    √ S S S * -1 (9.43) with
• ξ = -dVw dt the evaporation rate in m 3 /s, N the total number of nuclei, and A = √ Ks 4kπ . Γ measures the capacity of crystal growth to consume the supersaturation created by evaporation. If Γ is bigger than 1, evaporation is faster and the supersaturation will increase until the nucleation threshold is reached. On the contrary, if Γ is smaller than 1, crystal growth will be faster than evaporation and the supersaturation will decrease. The triggering of nucleation depends then on the current supersaturation, but also on the number and on the size of the nuclei already present in solution. If there are enough nuclei, crystal growth will consume supersaturation faster and the nucleation will not occur.

As a conclusion, the corner will have the behavior described in figure 9.10. Nucleation and crystal growth will fill the corner with crystals of different sizes, while drying reduces the available space. At the end, crystals becomes confined and the crystallization pressure is created. OF CCS Figure 9.10: Description of the evolution of the corner during evaporation and precipitation.

Algorithm

In order to solve equation (9.41), we will use the algorithm presented in figure 9.11. In this algorithm, we assume that we can calculate separately the nucleation of new crystals within the bulk pore solution, and the growth of existing crystals. At each time step, the new supersaturation after evaporation is calculated. Then, for all existing crystals, we calculate their growth during this time step considering equation (9.38). The new supersaturation after this growth is then computed and the corresponding nucleation rate is obtained. After nucleation we obtain the final supersaturation for this time step. The last step is to determine the occupation of the volume by the crystals. If the crystals are confined, then the crystallization pressure is created and the stress is determined with the classical equation (8.25) (zero curvature at the contact point).

Drying is calculated by simply considering a constant evaporation rate i.e. a linear decrease of the brine quantity in the corner (we neglect the influence of liquid pressure on the evaporation rate).

The stress is created when the crystals are confined. We consider that the crystals are hard spheres which are randomly packed (packing coefficient V crystals /V occupied = 0.7). The confinement appears then when the volume of the crystal packing is higher than the volume of the remaining brine.

Details of the calculation and inputs

In order to facilitate the calculation, we consider a pore of equilateral shape having three identical corners of length l/2 as shown in figure 9.12. The calculation is made in 2D and at 60 ℃ . The "pore volume" is then: The equation governing the evolution of the water saturation is:

v pore = √ 3 4 l 2 .
dS l dt = q CO 2 y v pore υ w υ CO 2 (9.44)
with q CO 2 , the flow rate of carbon dioxide in the pore, y the mole fraction of vapor in the CO 2 (considering y ≪ 1), υ w , the molar volume of water and υ CO 2 the molar volume of supercritical carbon dioxide. In order to avoid any problem of infinite values, the brine saturation cannot go below a minimal water saturation defined as the water saturation

CHAPTER 9. ESTIMATION OF CRYSTALLIZATION PRESSURE IN THE CASE OF CCS

corresponding to a wetting of the pore wall with a film of water of 1 nm (similarly to section 9.1.2).

At each time step, a new family of crystals is created, with a number of crystals (which have all the same size) corresponding to the nucleation rate. There is then as many families as time steps in the calculation, even if some are empty. Each family evolves then following the evolution of the supersaturation and its current size.

Finally, considering a flat pore wall, the curvature of the crystal in contact with the pore is null as represented in figure 9.13, and the crystallization pressure is then: The flow rate of carbon dioxide in the single pore is calculated by considering the values for a whole aquifer. If we consider the classical description of a radial aquifer of h=1m height and a CO 2 flow rate of 1kg/s, with a porosity φ, the flow rate in a pore at a distance R from the injection well is then:

ω p = RT υ s ln S ( 9
q CO 2 = Q CO 2 l 2 ρ CO 2 φR4h = 3.8.10 -4 l 2 φR (9.46)
Values of the flow rate will then depend on the porosity, the pore size and the distance to the injection well. With a pore of 10µm, a porosity of 7.6 %, and a distance of 10cm, we obtain q CO 2 = 5.10 -14 m 3 /s.

Finally, the calculation is very time consuming because of the numerous crystals present in the solution. As a result, the time steps have been chosen in order to minimize the error, while being large enough to allow moderate times of calculation.

MICROSCOPIC BEHAVIOR: NUCLEATION AND STRESS CREATION AT THE PORE LEVEL

Results of the simulation

The simulation has been processed for several carbon dioxide flow rates and several pore sizes. In figure 9.14 is plotted the evolution of the supersaturation (in black) and the water saturation (in red) in the corner (l = 10 µm, q CO 2 = 5.10 -14 m 3 /s and initial supersaturation S 0 l = 0.45).

As stated above, the evolution of brine saturation is linear because of the constant evaporation rate in the pore. It reaches its minimal value in approximatively 0.63s. The evolution of the supersaturation can be divided into three different phases. During the first phase, there is no nuclei in the corner. As a result, the supersaturation increases with the brine saturation decrease. A fit of the beginning shows well the evolution as 1/S 2 l . The first nucleation occurs at the end of this phase, when the supersaturation has reached the (homogeneous) nucleation threshold value. During the beginning of phase 2, nucleation is the main process by which supersaturation is consumed. Then there are enough crystals for the crystal growth to take the lead (Γ becomes smaller than one) and the number of nuclei becomes constant. The supersaturation decreases then, consumed by the crystal growth, tempered by the evaporation. Finally, when the brine saturation reaches its minimal value, evaporation stops. As the supersaturation is not maintained by evaporation, the only process remaining is crystal growth which leads to a fast decrease of the supersaturation to a value close to 1.1 (figure 9.14(b)). This final value corresponds to the equilibrium supersaturation of most crystals. During phase 3, crystals will evolve slowly toward a homogeneization of the crystals radii following Ostwald ripening.

In figure 9.15 is plotted the corresponding crystallization pressure. As expected, it is only created when the crystals become confined, which happens for a value close to 0.5s. Evolution of crystallization pressure can be easily related to the evolution of the supersaturation, we see clearly the end of phase 2 with the drop of supersaturation and then the phase 3. Interestingly, crystallization pressure during phase 3 can be fit with and exponential decay form:

ω p = ω max p -ω ∞ p exp - t τ + ω ∞ p (9.47)

Influence of the pore size and carbon dioxide flow rate

Two main parameters have a strong importance on the behavior of the corner: the CO 2 flow rate (i.e. the drying rate) and the pore size.

In figures 9.16(a) and 9.16(b) is represented the evolution of the maximum crystallization pressure (i.e. the crystallization pressure at the moment of its creation) with THE PORE LEVEL carbon dioxide flow rate and pore size (respectively for a corner size of 10 µm and a flow rate of 5.10 -14 m 3 /s). Evolution of maximum crystallization pressure with the flow rate is rather straightforward: the faster the evaporation is, the higher the crystallization pressure will be. Indeed, the faster the drying is, the less time the crystals have to consume supersaturation which means that at the moment of its creation, crystallization pressure will be higher than for a lower evaporation rate. The log/log scale allows to obtain a linear evolution of the maximum crystallization pressure with the flow rate.

The evolution with the corner size is less straight but still easy to understand. The smaller the corner is, the higher the crystallization pressure will be, because once again if a corner is smaller, there is less brine and the drying will happen faster.

Nucleation and crystal growth

As explained in the introduction, the calculation considers all the families of crystals of same size and their evolution with time as the evaporation proceeds. In this section we study these families and their evolution during the simulation. Let us consider figure 9.17. In this figure is plotted the critical radius of nucleation (i.e. a direct function of the supersaturation as shown by equation (9.40)), and the nucleation rate. Nucleation rate is significant only during a short period of time, between the moment the nucleation threshold is reached and the moment when there are enough crystals for the crystal growth to become preponderant (Γ > 1). As a result, the crystals are all created with similar radii as shown in figure 9.18: the initial radius of all crystals has a maximum variation of OF CCS What we can see from these figures, is that the biggest crystals will grow faster, while the smallest crystals will dissolve if the supersaturation decreases below their equilibrium value. Ultimately if given sufficient time, there would be only one size of crystals in the corner.

However, it is important to compare these sizes with the nucleation rate. Indeed, a bigger crystal means a smaller supersaturation at the creation and thus fewer big crystals than small crystals. There will then be a competition between the few big crystals which will grow and the numerous small crystals which will ultimately dissolve. As we can see in figure 9.20 where we have plotted the normalized number of nuclei N (t)/N max and the normalized nucleation crystal radius (r s (t)/r max s ), there are much fewer big crystals than small crystals. The maximum nucleation corresponds then to families which will dissolve. The behavior is thus dominated by the few crystals created at the beginning of the nucleation. In order to precise the relative importance of each family, let us now consider the figure 9.21. In this figure is plotted the volume occupied by each size of nuclei with time. We can indeed see what we expected, which is that there is a large distribution of the occupied volume, dominated by a mean radius, which corresponds to an average between the number of nuclei created and their size. It is interesting to see that this mean radius increases with time and that the distribution increases in height (which follows the crystallization).

As a result we see that the distribution of crystals during the simulation is broad and cannot be reduced to a single family of crystals.

Finally we can consider the evolution of the total number of nuclei with the flow rate. As expected, as we can see in figure 9.22 the faster the CO 2 is, the higher the number of nuclei will be. This can be seen in the value of Γ. Indeed, the faster the drying is, the stronger the nucleation has to be in order for the crystal growth to become efficient and compensate the evaporation.

Upscaling and poromechanics

The presented simulation allows to calculate the behavior and the crystallization pressure for a unique pore. In order to compare the results to the previous models, we now have to upscale this behavior and apply the laws of poromechanics as presented in figure 9.23. To do so we consider that the REV is composed of parallel pores of same size. Accordingly, the behavior of every pore will be the same. As a consequence, we can directly use the poromechanical calculation and the equivalent tensile stress as usual: In this model, we have neglected the capillary effects and thus the liquid pressure. The previous equation reduces then to:

̟ = b 2 + K + (4/3)µ N (S l + S s ) p l -p 0 l + S s ω p ( 9
̟ = b 2 + K + (4/3)µ N S s ω p (9.49)
The resulting equivalent tensile stress is presented in figure 9.24. We see two different regimes of ̟. The end of phase 2 corresponds to a higher value than the phase 3 because of the drop of supersaturation. The value during phase 2 is relatively stable because the decrease in crystallization pressure is compensated by an increase in the crystal saturation. On the contrary, the strain increases at first because of the increase of the crystal saturation, then decreases with the decrease of the crystallization pressure.

Conclusion on the modellings

In section 3.4, we have described three kinds of residual water after the end of the multiphase transport: wetting water, capillary trapped water and blob water. Blob water is defined as the brine remaining in the small pores while carbon dioxide is flowing in the big pores. Capillary trapped water is defined as the brine being trapped by capillary forces in the corners of the big pores. It is obvious that these two kind of brine will have very different behavior. Corners filled with capillary trapped water will dry very quickly because of the strong evaporation rate imposed by the flowing carbon dioxide. On the contrary, blob water drying is much slower because of the smaller interfacial area between CO 2 and brine. The modellings presented here describe the behavior of these two kind of remaining water: the macroscopic modellings (sections 9.1.1 and 9.1.2) describe the behavior of blob water, while the microscopic modelling of section 9.2 describes the behavior of the corners of the big pores. As expected the duration of drying is very different, from 1s for the corners up to several minutes or even hours for the blob water. The global behavior of the REV will be a juxtaposition of these two behaviors. At first, at small time scales, occurs the drying and the stress creation at the corner level, then at medium to long time scales occurs the crystallization and stress creation at the blob water level. Both phenomena depend strongly on the carbon dioxide flow rate and on the characteristics of the medium (capillary curve, corner size. . . ).

It is finally important to remind that the scenarios presented do not have yet an experimental verification. They are based on experimental and theoretical works from the literature, but no direct measurement and verification of the models have been proceeded. Experimental study of the crystallization of salt in rock cores subjected to supercritical carbon dioxide injection was intended as a part of this PhD work. However, because of technical difficulties, these experiments could not be carried-out properly. The description of the experimental set-up and the different issues met is the topic of the next part.

Part IV

Experimentation of reactive percolation

CHAPTER 10. GENERAL PURPOSE OF THE EXPERIMENTS of a triaxial cell in which the rock plug is inserted. A complete set of pressure sensors and flowmeters allows to obtain the intrinsic and relative permeabilities along with the evolution of the core saturation as a function of time (see chapter 13). The axial strain is measured with a LVDT in order to compare it to the macroscopic strain calculated in the modelling (see part III). Combined with these direct measurements during the injection, other studies were planned in order to obtain a whole overview of the behavior: comparison of the µ-CT scans prior and after the injection in order to obtain the clogging and the possible fracture pattern; chemical study by dosing the concentration of ions in the brine and comparison with the injected brine; SEM observation of the crystallization pattern and the crystal shapes. . .

During the initial development of the prototype, the emphasis was put on the decoupling of the different processes involved. Indeed, as we have seen in chapter 6.3, injection of supercritical carbon dioxide in a brine saturated rock leads to a strongly coupled THMC behavior. To solve this problem and consider the processes separately, different solutions were implemented: modification of the rock core and its chemical reactivity (sintered glass, sandstones and limestones; chapter 12); modification of the gas (nitrogen, gaseous or supercritical carbon dioxide, saturation of CO 2 with water before injection to decrease the influence of evaporation, section 11.1.2.2); modification of the composition of the brine (pure water, NaCl, CaSO 4 . . . ).

Finally, in order to consider a real case of CCS, different scenarios of experiment were planned: continuous injection of CO 2 only; of CO 2 and brine; or discontinuous scenarios alternating injection of both CO 2 and brine in order to study the effect of a interruption of the injection and the following resaturation by the aquifer brine.

Unfortunately, because of the novelty of the prototype, the initial development delays and the resolution of the different issues faced during the experiments forbade to follow this experimental plan. Among the most important issues faced is the clogging of the tubing with salt crystals (section 14.2). This problem, still unsolved at the moment of the redaction of this thesis made impossible to proceed with drying-out experiments and to study the different scenarios presented here. Nevertheless, the prototype was tested and is found to perform well on intrinsic and relative permeability measurements combined with capillary curve measurement, which gives interesting results as presented in chapter 13.

Description and sketch of the prototype

The final system built is shown in figure 11.1. In the center is the triaxial cell. In this triaxial cell is inserted a rock core within a sleeve which will transmit the confining pressure. Below the triaxial cell is the bottle to mix the fluids. This mixing bottle is designed as a storage for the water and as a mixing tool (see 11.1.2.2). Carbon dioxide is injected by pump ST2, while pressure in the bottle is controlled by pump ST1. Confining pressure is set by pump ST3.

The set-up is designed to be controlled with a computer software. The control panel is presented in figure 11.2. On the top left is the pump ST2 (11.1.2.4) and on the bottom left is the reservoir of brine connected to the bottle (11.1.2.2) in the bottom center. Above the bottle is represented the triaxial cell (11.1.2.1) with the two manual injection valves V1 and V3 (and the corresponding purges V2 and V4). The light gray rectangle represents the climatic chamber (11.1.2.3). Finally on the bottom right are the pumps ST1 and ST3 respectively controlling the piston of the bottle and the confining pressure. On the top right is the micro-pump controlling the back pressure valve (11.1.2.5) (represented on the close right of the triaxial cell). Left to this micro-pump is the fluid measurement system with the gasometer and the weighting scale (11.1.2.6). The color code is described in table 11.1. In addition to the manual valves like V1, there are several electrical valves named EV and controlled by a simple click on the control panel.

Finally, the software has several built-in procedures allowing to automate some actions like injection, pressure regulation . . .

Material

Metallic parts

The set-up is supposed to work properly in a very aggressive environment. Indeed, as expressed above, the system is supposed to work with a binary mixture of supercritical carbon dioxide and brine. Strong corrosion is to be expected with this mixture and the material constituents of the set-up had to be chosen accordingly. For the triaxial cell and the bottle which suffer the strongest aggression of the mixture, we chose hastelloy. The other parts, in particular the tubing, are build in inox. As we can see from pictures 11.3 corrosion still happens on the inox. In particular, in the brine reservoir, strong corrosion occurred. The bottom was then replaced by Teflon.

Sleeve

The aim of the sleeve is to transmit the confining pressure of the confining fluid in the triaxial cell. The material of the sleeve was difficult to choose. Indeed, it has to sustain high pressure and temperature, and CO 2 presence. The different polymer available are Viton, Nitrile and Kalrez. As these are polymers, carbon dioxide acts as a solvent and can penetrate within the sleeve. While there is sufficient confining pressure this penetration has no consequence, but as soon as the confining pressure is released, the polymer swells strongly (see figures 11.4). Viton is the one which reacts the worst with carbon dioxide, Nitrile is a little better and Kalrez is the best. Regarding temperature, Viton can bear temperatures up to 200℃ while Nitrile can withstand only 120 ℃. Kalrez is the best polymer to make sleeves, however, because of its price, we have decided to stay with Viton and Nitrile sleeves, changing them regularly after experiments.

Description of the different parts of the set-up 11.1.2.1 Triaxial cell

The design of the cell is presented in figure 11.5. The confining pressure is transmitted to the core through the impermeable sleeve and the fluid filling the cell. The part in dark gray in the left of the picture is mobile and is used to transmit the deviatoric pressure controlled by the pressure of the fluid in the small chamber on the left. There are two fluid inlets (one for water and one for gas) but only one can be seen on the figure, the other being behind. The cell characteristics are presented in table 11.2.

Fluid bottle

The mixing bottle is an important part of the set-up. It has two purposes: it serves as a storage for the water to be injected, and as a mixing device. Design of the bottle is presented in figure 11.6. It consists of a body and a piston which is controlled by pump ST1. The particularity and aim of this bottle is the three inlet/outlet which allow a prepartionning of the two fluids. Indeed, gas can be injected through the middle inlet to make it bubble in the water. As a result, from the bottom outlet we can produce water saturated with carbon dioxide, and from the top outlet we can produce carbon dioxide saturated with water. However, in most experiments reported in this dissertation, only pure water was used. Bottle characteristics are reported in table 11.3.

Oven

The oven is a classic climatic chamber (MEMMERT) allowing to set the temperature of the whole set-up (triaxial, bottle and back pressure valve). It is however possible as shown in figure 11.7 to work outside the oven by adding additional tubing and by using the cart. This allows an easier access to the different parts of the set-up. This characteristic however was not very useful as the additionnal tubing generates a change in dead volumes (see 11.2). Finally, the circuit of CO 2 includes a coil in order to heat the carbon dioxide which is stored outside the oven in the pump and to allow its phase change (from liquid to supercritical) before the injection. The temperature can be set up to 250 ℃ (remembering that several parts of the set-up have lower maximum working temperatures) and the precision on temperature is 0.5 ℃ from ambient to 100℃ and 1℃ above.

Pumps

There are three pumps in the system which control fluid injection and axial/deviatoric stress: ST1 controls the piston of the bottle and thus the injection of brine in the cell. ST2 controls the carbon dioxide injection and ST3 controls the confinement of the cell. These STIGMA1000 pumps are specific of Sanchez Technologies. They can work either in pressure or flow rate control but are more efficient for the latter. STIGMA1000 pumps are shown in figure 11.8 and their technical characteristics are reported in table 11.4.

Back Pressure Valve

The back pressure valve (BPV) controls the pressure at the outlet of the system. This valve is controlled by a micro-pump and acts as a classic differential valve: if the pressure of the incoming fluid is higher than the BPV instructed value, the valve opens and the fluids flows through. On the contrary if the fluid pressure at the exit of the rock core is lower than this instructed value the valve stays closed. The BPV is sketched in figure 11.9 and its characteristics are reported in table 11.5.

Fluid separation and measurement

After the BPV, the two phases are separated and their quantity measured with time. The separator presented in figure 11.10 uses simply the difference of density, brine staying at the bottom while the gas is going up. At the bottom of the separator is a valve controlled by a sensor which detects the passage of the meniscus. Before the experimentation the separator is filled with water until the meniscus reaches the sensor. As a result, each drop of liquid coming in the separator moves the meniscus up and triggers the sensor which opens the valve. A corresponding quantity of brine is then extracted from the system and added to the beaker on the weighting scale to be measured. The valve closes when the meniscus goes below the detection range of the sensor.

The gas quantity is measured by a gasometer. Its principle is simple: the incoming gas forces a wheel to rotate and the number of turns the wheel does informs on the volume of gas. In order to have the proper resistance on the wheel rotation, the gasometer is half filled with a precise volume of demineralized water. Its precision is 5.10 -5 dm 3 and it works at atmospheric pressure and ambient temperature. Note that between the separator and the gasometer is added a condensation trap in order to remove most of the water from the gas. Finally the separator is covered with a heating system able to heat it to temperatures up to 60 ℃.

Dead volumes measurement

As can be seen in the figures and the description of the system, between the active parts (pumps, rock core, separator...) are present numerous tubings adding volume to the system. This volume which we will call dead volume has to be measured in order to get precisely the volume of injected fluids and calculate the saturation of the core during the injection (see section 13.2). In figure 11.11 are represented the considered dead volumes.

We divide the dead volumes in three categories:

• dead volume in: volume of tubing followed by the injected fluids under pressure (before BPV)

• dead volume out: volume of tubing followed by the injected fluids at atmospheric pressure

• dead volume sensors: volume of tubing leading to the pressure sensors P6 and P10

At first we measure the dead volume under pressure i.e. the sum of the dead volume in and the dead volume sensors. There are several methods to measure these dead volumes. The method considered here consists in varying the volume with the pump ST2 and measuring the resulting pressure (which is considered homogeneous in the part considered). As the matter quantity is constant, the dead volume can be simply deduced from the intercept of the curve using the perfect gas law. Indeed, we have:

RT p = V measured + V core n + V ST 2 n (11.1)
Figure 11.12: Steel tube for dead volume measurement.

considering n the gas quantity, V ST 2 the volume of gas in the pump, V core the empty volume of the core used in the triaxial cell and V measured the volume we are looking for.

Method

In order to get the volume corresponding to the dead volume in and the dead volume sensors, we did two different measures: we first measured the volume of tubing between the pump ST2 and the valve V3, and then we measured the total volume between the pump ST2 and the BPV. We then have the desired dead volume by substracting both values. The measurement is proceeded as follows: the pump ST2 is set at different volumes and equilibrium pressure is recorded. The equilibration time can range from a few minutes to several dozens. The smaller the measured volume is, the longer it takes to the pressure to equilibrate.

In this measurement, it is preferable not to use a real rock core but rather an empty cylinder of steel in order to know precisely the volume V core . The cylinder used is shown in figure 11.12

Results

Results of the experiment are presented in figure 11.13. The slope of the curves allows to get the matter quantity in the volume considered 11.2. DEAD VOLUMES MEASUREMENT ST2/BPV ST2/V3 V3/BPV 44.2 mL 37.9 mL 6.29 mL Table 11.6: Dead volumes in the system. and the Y-intercept gives the dead volume in the system according to equation (11.1). It is interesting to remark that when volumes are small (for example the beginning of 11.13(a)), the equilibration time is very long and the result is very sensitive to the conditions. This explains why the first value in figure 11.13(a) is a little off the fit.

Resulting dead volumes are given in table 11.6. The intrinsic dead volume of the triaxial cell according to Sanchez Technologies is 0.3mL. Consequently, the dead volumes are mainly composed of the tubing linking the different parts of the set-up.

The dead volume measured here are not only the tubing which will be used by the injected fluids but are also the tubing leading to the different pressure sensors (P10 and P6) (dead volume sensors). In most of our experiments, theses volumes will not be considered as the fluid filling them will not change (no flow in these tubings). The method presented here cannot then gives us the useful dead volume i.e. the dead volume actually followed by the moving fluids. We then had to resort to a less precise method to obtain the dead volume we are looking for: we measured the length of tubing, the inner diameter being 1.6mm. The considered dead volume in is then approximatively 2mL, while the dead volume out is 2.5mL. This solution is unfortunately less precise because of the numerous bends of the tubings. The measurement of the length of the tubings leading to the sensors however gives a volume for the dead volume sensor of approximatively 4.5 mL confirming the validity of the measure (the addition of the measured dead volume in and the dead volume sensor gives approximatively the 6.29 mL measured with the pressure-step method).

Rock cores

The system is designed to work with cores of porous materials. In our experiments we have considered two sorts of porous materials: natural rocks and sintered glass. Natural rocks are easily accessible but present a more complex structure and chemical composition than sintered glass. However, the sintering technique leads to very high permeabilities and the pressure drop was then not clearly discernible. As a result, experiments were rather carried out with natural rocks. The rocks used were reservoir rocks of two types: limestones and sandstones. Limestones are rocks mainly composed of calcite (CaCO 3 ) while sandstones are mostly made of quartz SiO 2 . For feasibility all the rock cores were drilled in quarried slabs.

Pierre de Lens

The Pierre de Lens is an oolithic limestone with 99 % of calcite [START_REF] Saad | Influence du changement climatique et des conditions extremes sur les massifs fracturés, Role des fluides dans leur processus d'altération[END_REF]. It is quarried in the Gard region in France. The rock were drilled in slabs by F. Martineau (IFSTTAR). Porosity was measured by saturating the sample with water and weighting it. We found a mean value of 13.5 % for the different cores used (values range from 10 to 17 % in the literature). In figure 12.1(a) is shown a rock core of Pierre de Lens

Grès des Vosges

The Grès des Vosges is a sandstone quarried in the Rothbach quarry in the Vosges Mountains (France). Its composition is almost exclusively quartz. The measured porosity for the core used for the experiments is 18.05 %. The results of the analysis of the rock made by the BRGM are reported in appendix F (in french). In figure 12.1(b) is shown a core of Grès des Vosges

Rocks from the Dogger aquifer of the Paris Basin

In order to study real reservoir rocks, we sampled some rock core in the TOTAL core library in Boussens (a warehouse where samples of rocks drilled in France by TOTAL are stored). We decided to sample several rocks from the Dogger aquifer of the Paris Basin. The different rock drilled and their description can be found in [START_REF] Badinier | Proposition de méthodologie de caractérisation géologique et pétrophysique d' échantillons de roches réservoir[END_REF]. However, because of the lack of time and the experiment difficulties, we did not use these samples in the experiments. Nevertheless, the modelling has been proceeded with the characteristics of two of these rocks: the Mondeville sandstone (Boissy sandstone) and the Bois-Brûlé

Sintered glass beads

In order to study a model material, we decided to produce sintered glass. Sintered glass is a model material composed of monodisperse and aggregated glass beads created by heating the SiO 2 beads in an oven at a temperature below the melting temperature of the glass but high enough to be in a plastic regime for the glass (usually between 600 ℃ and 660 ℃). The surface of the beads partially melts and the beads are then partially fused together. During the cooling, the glass solidificates and the beads are then tightly and cohesively held together forming an ideal and perfectly known porous medium.

Sintering temperature and duration

Petrophysical properties of the produced material (porosity and permeability) depend on three parameters: size of the beads, sintering duration and sintering temperature. Indeed, the higher the temperature, the more the plastification of the beads is pronounced leading to a more compact and dense material and thus to a low permeability and porosity. In the same way, sintering duration increases the density of the material allowing more time for the beads to fuse. Charts [START_REF] Gauthier | Courbes de porosité sur frittés[END_REF] allow to determine the properties of the final material as a function of these three parameters. However, they are not comprehensive and depend particularly on the oven used and the position of the beads in this oven (for example through the homogeneity of the temperature in the oven and its heating capacity). . . Empirically, a temperature program has been determined in order to lower the density variations, the shrinkage of the material during the heating and to reduce the thermal strains. It is presented in figure 12.2.

Moulding of the material

Glass beads are set in a mould which has to resist the high sintering temperatures. In particular, the mould should not present any notable thermal deformation. Two kinds of moulds are possible which respect these conditions: refractory plaster and quartz tubes.

Refractory plaster mould

At first we tried a plaster mould created with a negative of the wanted shape of sintered material. The plaster used was MOLDA 3 Normal (St Gobain) and the ratio water/plaster used was 0.5.

Plaster is poured in a plexiglass mould with the negative covered in oil to avoid any adhesion within. This plexiglass mould is also made of two distinct parts which are sealed together with parafilm and aluminum scotch. At first we wanted to obtain a reusable mould by inserting a plastic sheet around the negative in the plexiglass mould in order to obtain two separate halves of the plaster mould (see figure 12.3). However, the plastic sheet was difficult to insert and the two halves after demoulding did not fit perfectly. We tried then to glue the two parts with a refractory mastic (Rubson Haute-Temperature), but the differential strain of the mastic and the plaster lead to fracturing of the mould and to the leakage of the beads. This forbade the use of micrometer sized beads and was thus abandoned. The plaster mould was then created in only one part which would be destroyed in order to retrieve the sintered glass core. This technique is efficient and allows to create the wanted material. However, it is time consuming as one mould has to be created for each core. Moreover, the moulding of the plaster around the negative does not create a smooth surface. As a result, the core is not well shaped and its surface is particularly rough. It then has to be trimmed before use. Finally, another issue with this technique is that the plaster after curing still contains a lot of free water. As a result, the mould has to be put in an oven before use in order to remove all this water. If not, vaporisation of this water will fracturate the mould and lead to the leakage of beads. Figure 12.4 shows a core obtained with a plaster mould. We clearly see the rough surface, and a small leakage of the beads on the top face of the core. Therefore, we decided not to continue with plaster moulds but rather use quartz tube which are much easier to handle.

Quartz tubes

Quartz is a transparent material which resists well the temperatures used for sintered glass. Quartz tubes are really easy to handle: the tubes are put on a steel sheet and beads are poured in it. At the end of the sintering program, the cores are simply removed from the tubes which can then be used for another sintering. This technique is then much more faster and easier to use than the plaster moulds. However, it is less flexible as the diameter of the tubes and their height is fixed. It is also much more expensive. Figure 12.5 shows the result of the tube moulding.

Shrinkage and density profile

During sintering the density of the beads assemblage increases. However, this density increase is not homogeneous because of gravity: the lower layer supporting all the weight of the upper layers, will have a higher density than these latter. In the case of small cores with small height, this density effect can be neglected. However, if we want a very controlled and homogeneous density profile, a solution is to produce a huge core and then to drill small ones horizontally in a layer of homogeneous density. Another solution is to put on top of the poured beads a weight which will render negligible the effect of the beads layer. However, this create radial stresses on the tube which almost always fracture it. In our case, the height of the sample being small, we neglect this density effect.

The beads compaction increase is not radially homogeneous in the tube and is more pronounced on the center creating a convex shape on the top of the core. Trimming is then necessary to obtain a real cylindrical shape. As expected, the higher temperature and the longer the sintering, the more pronounced this effect is.

Finally, another shrinking effect justifies the choice of the quartz tube: for plaster mould, it often happens that the diameter of the core is reduced during sintering, this effect being more pronounced in the middle of the core. )) by forcing mercury to penetrate into the sample with increasing pressure steps. As a pressure corresponds to a radius of a meniscus and thus to a pore radius (the contact angle between mercury and the rock matrix being known), the intrusion volume added during the pressure increase gives the porous volume corresponding to the said pore radius. By sweeping over a large range of pressure, we can then obtain the volume corresponding to each pore size present in the sample.

Material and methods

The MIP is processed on a AutoPore IV from MicroMeritics Instrument Corporation. The mercury used comes from Ophram Laboratories. There are several routine measurement already built in the system, we have chosen the "carbonate" method. It does 2 series of measurement, one at low pressure and the other at high pressure on two different parts of the set-up. The mercury parameters considered are 130°for the contact angle and 485 J/m 2 for the surface tension. Equilibration time (time to wait between the pressure step and the intrusion volume measurement) is 10s for low pressure and 60s for high pressure. A standard blank measurement is used to remove the various offsets of the machine. Measurement of the intrusion volume is made with a capacitive method with a constant for the penetrometer of 10.898 µL/pF. The measuring pressure range from 3.45.10 -3 MPa to 227.5MPa, which corresponds to pore radii from 2.7nm to 180 µm. To be inserted in the porometer, the rock to be measured is broken in small chunks. Several chunks are

Results

MIP has been processed for the different rocks from the core library (figure 12.7), for a sample of Grès des Vosges and, a sample of Pierre de Lens (figure 12.8) and two samples of sintered glass (figure 12.9). Distribution of pore sizes for natural rocks ranges from 0.1 µ to 20-50 µm. As can be expected from the results of porosity and permeability measurement (see chapter 13), the Grès des Vosges presents much bigger pores than the Pierre de Lens, which explains both the intrinsic permeability difference and the residual liquid/gas saturations. Finally, in figure 12.9, we can clearly see the influence of the glass beads diameter on the pore size distribution (both core have a sintering temperature of 640℃ and a sintering time of 60min). Unfortunately, the technical limitations of the instruments do not allow to obtain the whole pore size distribution for the core with the biggest beads.

Intrinsic permeability 13.1.1 Method

Measurement of intrinsic permeability is based on Darcy's law as presented in section 3.2.2. The idea is to inject a fluid at different flow rates and to measure the difference of pressure between the inlet and the outlet. The slope of the straight line obtained gives the permeability, knowing the viscosity and the geometry of the core sample of the considered rock. Intrinsic permeability has been measured for several rocks with gaseous nitrogen and pure water. Depending on the permeability of the medium considered, the flow rates have been chosen to obtain a significant pressure difference. The back pressure valve is set at 2 bars (the BPV pressure is not precise at 1 bar) for Pierre de Lens and 100 bar for the Grès des Vosges. The equilibrating time for the pressure difference at each flow rate depends on the permeability of the medium: the smaller the permeability is, the longer it takes to reach the equilibrium value. Consequently, in order to obtain the right equilibrium difference (i.e. at infinite time), an exponential fit has been processed as presented in figure 13.1 for each flow rate. The equilibrium value of the pressure difference is then the value at infinite time of the obtained fit.

Results

Results of measurement of intrinsic permeability show that the system can obtain accurate values of pressure difference under a constant flow rate. Indeed, the curves obtained are perfectly linear. Here are presented the different results obtained with limestone Pierre de Lens and sandstone Grès des Vosges.

In figure 13.2 are plotted the evolutions of the pressure difference with the flow rate for pure water 13.2(a) and gaseous nitrogen 13.2(b) and their linear fit. Note the difference in pressure difference for the liquid and for the gas.

Considering now Darcy's law (3.10) for saturated case with no gravity effect, the slope of the figure is given by the relation:

Q ∆P = k 0 A ηL (13.1)
Q being the flow rate, ∆P , the pressure difference, k 0 the intrinsic permeability, A the section of the sample, L its length, and η is the viscosity of the fluid depending on the temperature of the study. For the limestone, measurement for water has been performed at 26.5 ℃ leading to a viscosity of 0.860 mPa.s and measurement for nitrogen have been performed at 15 ℃ leading to a viscosity of 1.8.10 -5 Pa.s. For the sandstone, measurements have been processed at 60 ℃ and 100 bar (in order to be in the conditions of the measure of relative permeability water/CO 2 , see below). The rock sample of limestone was a cylinder or radius 1.5" and of length 3" while the sample of sandstone had a diameter of 1.5" and a length of 8.5cm. Unfortunately, measurements of permeability for the sandstone with gas could not be made as the pressure difference was to low because of the high permeability. Results are presented in table 13.1.

We see that both results for Pierre de Lens are close (4 % variation). This is an interesting validation of the accuracy of the set-up. In [START_REF] Saad | Influence du changement climatique et des conditions extremes sur les massifs fracturés, Role des fluides dans leur processus d'altération[END_REF] is measured the permeability of a Pierre de Lens rock core. They obtain a value of 409 µD. The order of magnitude is correct, the difference being easily explained by the difference of provenance of the two cores.

Methods

Measures of relative permeability are much more difficult to carry-out than intrinsic permeability because of the presence of the two fluids in the porosity. However, the principle is the same, injection of the two fluids at constant flow rate and measure of the pressure difference between the inlet and the outlet. The values of the two flow rates are chosen in order to sweep the whole range of water saturations. However, this technique is limited at the low saturation and values close to the residual saturation cannot easily be obtained. Experiments have been proceeded with the core initially saturated either with water of with the gas phase. The latter allowing to reach more easily very low water saturation.

Three relative permeability measurements were realized, one with Pierre de Lens limestone with the couple N 2 /water, another with the sandstone Grès des Vosges with gaseous carbon dioxide and water, and finally the last measurement was made with the same sample of Grès des Vosges with supercritical carbon dioxide and water in order to determine the influence of the supercritical state on the relative permeability (measurement of relative permeability with CO 2 and water is impossible for the Pierre de Lens because of the dissolution of calcite during the injection). The first two experiments have been proceeded at low pressure (BPV set to 5 bars and 2 bars respectively) and ambient temperature, while the second experiment needed high temperature and pressure to remain in the supercritical domain (BPV set at 100 bars and temperature set at 60 ℃). Whereas the first two experiments do not present experimental issues, the third is much more complex. Indeed, at the end of the experiment, the whole system has to be put at 1 bar in order to retrieve the sample. This decrease of pressure transforms the supercritical carbon dioxide into gaseous carbon dioxide. This phase change and the volume change coming along can be the cause of variation of the fluid content in the core. Moreover, the density of water is different at 100 bar and 60 ℃ and the saturation has to be corrected accordingly.

Finally, using the definition of relative permeability recalled in section 3.2.2 and the intrinsic permeability measured in the previous paragraph, we can obtain the following relation for the relative permeability:

k r J (S l ) = Q J A η J k 0 L ∆P (13.2)
where J is the considered phase. In table 13.2 are presented the different flow rates used in the experiment and the corresponding water saturation. The saturation of the core was measured by dismantling and weighting it.

Q G (cc/min) Q l (cc/

Grès des Vosges

In table 13.3 are presented the flow rates used in the measurement of relative permeability of the sandstone with gaseous carbon dioxide and water.

In order to reach the very low water saturation for the last injection for Grès des Vosges we injected before the experiment a small and precise quantity of water. Then carbon dioxide only was injected and the pressure difference recorded. As the relative permeability of water is really low in this case (low water saturation), the loss of information on the water relative permeability by not injecting water is not a problem. This allows to avoid the problem of spontaneous water imbibition of the sample when water is injected along with carbon dioxide. In the case of supercritical conditions, measurement of the relative permeability is much more difficult because of the high pressure and temperature conditions. In table 13.4 are presented the flow rates and the weighted water saturation used in the experiment.

The water saturation presented in table 13.4 is the weighted water saturation and is different from the real water saturation at 100 bar and 60 ℃ because of the change of density of the two phases during the depressurization and additionally because of the evaporation associated with the depressurisation. In order to estimate the real water saturation during the injection, we consider the volume of water injected during the experiment and compare it to the volume extracted during the experiment.

Let us consider the conservation of water in the cell and the dead volumes at the beginning of the injection (0) and at the end of the injection (1). We then have the following equation:

m 0 core + m 0 dead -m out + m inj = m 1 core + m 1 dead (13.3)
where m core is the mass of water in the rock core, m dead is the mass of water in the dead volumes, m out is the mass of water obtained at the exit of the system and m inj is the mass of water injected during the experiment. Expressing the different masses, we obtain:

S 1 l φπr 2 Lρ * w -S 0 l φπr 2 Lρ * w + ∆m dead = V inj ρ * w -m out (13.4)
ρ * w is the density of water under 60 ℃ and 100 bar, φ is the porosity of the rock core, r its radius and L its length, ∆m dead = m 1 dead -m 0 dead is the variation of the mass of water in the dead volume during the experiment.

As described in the section 11.2 the dead volumes are divided in three parts: the tubing leading to the pressure sensors, the dead volumes between the triaxial cell and the BPV, and finally the dead volumes after the BPV. Considering the injection process, CHAPTER 13 we assume that the dead volumes leading to the pressure sensor are not modified by the injection process: if they were filled with water (respectively gas) they will remain filled with water (respectively gas). Due to the expansion of the gas phase after the BPV, we observe that this dead volume (dead volume out) is filled only with carbon dioxide at the end of the experiment. However, between the core sample and the BPV, the two fluids penetrate alternatively the tubing. Neglecting the variation of pressure in the tubings, the water saturation in the tubing between the rock core and the BPV can be approximated by the comparison of the flow rates of liquid and gas

S l = Q l Q G +Q l .
As a result, we have for the dead volumes:

∆m dead = Q l Q G + Q l V bef ore ρ * w -S initial (V af ter ρ w + V bef ore ρ * w ) (13.5)
with V af ter the dead volume out after the BPV, V bef ore the volume between the rock core and the BPV, S intial is the initial water saturation in the dead volumes and ρ w the density of water at ambiant pressure and temperature. S initial = 1 if the rock was initially saturated with water and 0 if not.

Finally, using equations (13.4) and (13.5), we obtain the water saturation in the rock core during the experiment:

S 1 l = S 0 l + 1 φπr 2 L S initial V af ter ρ w + V bef ore ρ * w S initial - Q l Q G + Q l + V inj ρ * w -m out ( 13
.6) In table 13.5 are presented the corrected water saturation.

It is important to remark that the corrected water saturation are always higher than the weighted water saturation, which is perfectly normal, mainly because of evaporation. Paradoxically this method is not accurate to measure the water saturation for the low BPV pressures. Indeed, the conception of the mixing bottle does not allow to remove all the gas from it, leaving a gas cap at the top of the reservoir. As a result, the volume change associated with the pressure increase in the bottle is in the most part due to the gas compression and not only to water displacement. At high pressure (100 bar), the gas in the bottle is already compressed and the pressure increase because of the injection is only of 1% for a 1 bar pressure gradient, compared to an increase of 50 % for the 2 bar BPV. The injected volume method is then possible at high pressure and impossible at low BPV pressures.

Results

Pierre de Lens

In figure 13.3 is presented the resulting relative permeability of the Pierre de Lens to nitrogen and pure water. Relative error on the permeability, considering that the error on the pressure difference is 0.01bar (minimal possible variation), is of an order of magnitude between 0.1 % and 1%. A fit of the curve with the Corey and van Genuchten's equations (see 3.2.2) is also plotted in the figure. The best fit is obtained with m = 0.5 and S lr = 0.41. The last points at low saturation for the gas relative permeability do not belong to the Corey equation fit because they are below the residual water saturation.

Grès des Vosges

In figure 13.4 are presented the relative permeability of the Grès des Vosges with supercritical carbon dioxide and pure water 13.4(a) and with gaseous carbon dioxide and pure water 13.4(b), and the van Genuchten/Corey fits. For supercritical carbon dioxide, the fitting parameters are m = 0.65, S lr = 0.05 and S gr = 0.17, while for gaseous carbon dioxide, the fitting parameters are m = 0.62, S lr = 0 and S gr = 0.23. However, the fit for the gas relative permeability is not entirely accurate at high water saturations. Moreover, the lack of points at low water saturations (because of the difficulty to reach these points) creates a strong uncertainty on the determination of the residual water saturation. Indeed, in [START_REF] André | Impact of porous medium dessication during anhydrous CO2 injection in deep saline aquifer: up scaling from experimental results at laboratory scale to near-well region[END_REF] is presented a curve of relative permeability for a sample of Grès des Vosges obtained from the capillary pressure curve. They found a residual water saturation S lr = 0.19. However, their determination of the coefficient m = 0.675 is close to our value (relative difference of 8 %) which confirms our measurement.

From figure 13.4 and the fitting curve it is clear that the two curves are similar, especially the two water relative permeabilities. As a result, it is possible to conclude the values of the relative permeabilities measured with gaseous carbon dioxide can be used for supercritical CO 2 with good accuracy. Finally, we can see that the curve for the sandstones differs strongly from the limestone as the residual water saturation are strongly different: the limestone retains more water that the sandstone, which can be explained by the different pore size distribution, the pores of the Pierre de Lens are much smaller than those of the Grès des Vosges (see section 12.3). Capillarity effect are then less pronounced for the latter leading to a smaller water residual saturation. Everything else being the same, it seems that the Grès des Vosges presents better characteristics that the Pierre de Lens for storage purposes because of the lower water retention S lr and the higher water mobility: it is then possible to remove more water and store more CO 2 in rocks like Grès des Vosges than in Pierre de Lens.

Capillary pressure measurement

The system can also be for the measurement of the capillary pressure curve (see section 3.3.2) of the Grès des Vosges sandstone with carbon dioxide and pure water. The set-up was not initially designed for this kind of measurement, however, a good approximation of the capillary curve (initial drainage) could be obtained.

Method

In order to get the capillary pressure curve the measurement was proceeded as follows: at first the rock was saturated with pure water and the BPV set to 5 bars. The valve V1 was then closed and the valve V3 opened to connect the rock plug to the gas system. The carbon dioxide pressure is then increased by small increments (0.1 bar) from 5 bars and the extracted water quantity at each increment is recorded. Capillary pressure is obtained by the measurement of the pressure difference between the fixed carbon dioxide pressure and the BPV pressure. We can then sweep the whole range of capillary pressure and thus obtain the capillary pressure curve. However, this methods is based on several questionable hypotheses: we consider that the water pressure is homogeneous and equal to the BPV pressure and we consider that the water saturation is homogeneous in the sample. Both of these hypotheses are partially wrong as argued in [START_REF] Dana | Experimental study of two-phase flow in three sandstones. I. measuring relative permeabilities during two-phase steady-state experiments[END_REF][START_REF] Dana | Experimental study of two-phase flow in three sandstones. II. capillary curve measurement and relative permeability pore space capillary models[END_REF]. Nevertheless, as we will see in the next paragraph, the resulting curve is relatively accurate with the van Genuchten fit obtained from the measurement of the relative permeability of the considered sandstone (see figure 13.4).

The liquid saturation is obtained similarly to the calculation for the relative permeability: we substract from the extracted water volume the dead volume in order to get core water saturation, the tubing being initially full of water. However, for the second measure (∆P = 0.06 bar, see table 13.6), there is no percolation: the gas penetrates the rock core, displacing some water, but does not percolate through the entire core. As a result, the dead volume is not substracted from the extracted water volume and the water saturation is directly calculated with the extracted volume as there is no change in the fluid occupation of the dead volume.

Results

Results of the measurement are presented in the table 13.6 and are plotted in figure 13.5. In dots are the experimental values and in dotted line is the curve obtained with the van Genuchten equation (3.20) of the same core (m = 0.62, S lr = 0.05), with a value for the entry pressure p 0 = 0.06 bar. This value is obtained from the situation of the second measured point i.e. penetration but no percolation.

As we can see from the results, the measured capillary pressure corresponds accurately to the curve obtained with the coefficients extracted from the relative permeability curve. The slight dispersion of the points is caused by both the hypotheses of the measurement (homogeneous water saturation and pressure in the core) and the imprecision on the value of the dead volume. However, the results show that a capillary pressure estimation of a core can be accurately and easily obtained with the experimental set-up.

Issues and solutions

As the system is a new prototype, there were several technical issues with the system which needed to be solved. The main problem during the measurement of intrinsic and relative permeability was the separator and fluid quantity measurement (gasometer and weighting scale).

Carbon dioxide pressure and leakage

The gasometer used in the set-up is supposed to work at a pressure not exceeding 65 mbar above the atmospheric pressure. However, depending on the injection rate of supercritical carbon dioxide, a big quantity of gas can arrive to the gasometer as there is a phase change from supercritical to gas after the BPV (the volume change from supercritical carbon dioxide at 60 ℃ and 100 bar to gaseous carbon dioxide at 25 ℃ and 1 bar is 1 for 161). As a result, if the carbon dioxide flow rate is too high (above 1 to 2 cc/min), the gasometer is under too much pressure and there is a leakage of the water filling the gasometer (see section 11.1.2.6). This results in a loss of gas which can also leak with the water, but also on a decrease of the resistance on the wheel of the gasometer leading to erroneous gas measurement. This phenomenon happens also during the depressurization 13.4. ISSUES AND SOLUTIONS for the relative permeability of supercritical CO 2 and water. A too fast depressurisation leads to the same problem.

For very high amount of gas (very high flow rate or BPV oscillations (see below)), the gasometer is not the only one to be impacted. Indeed, the plastic tubing at the bottom of the separator (the one with the valve for water) can also leak if the carbon dioxide pressure in the separator is too high, leading to a loss of both water and gas. Finally, in rare cases of very strong gas quantity, the plastic plug closing the condensation trap is ejected.

As a result, the flow rate of gas exiting the BPV has to be controlled either by setting a low CO 2 flow rate, or by depressurizing slowly the system at the end of the experiment.

BPV oscillations

As said above, there is a phase change after the BPV as the carbon dioxide is not under supercritical conditions anymore. This leads to a volume change of approximatively 1 for 161 and also to a change of compressibility. The BPV works as a differential pressure valve. This means that it uses the incoming pressure of the fluid to build-up its own pressure. When the fluids begins to be injected, the pressure before the BPV increases. This increases the pressure of the BPV itself which then works to decrease this pressure to reach its instructed value. The valve then opens and fluid flows. But as there is a phase change and because of the compressibility of the new phase, the BPV finds itself pushing not on a dense phase as it had, but on a gaseous phase. This phase change modifying the environment of the BPV creates a drop in the BPV pressure, this latter working then to increase it again to the instructed value. The valve then closes and the whole cycle restarts. For high permeability systems (Grès des Vosges) these oscillations of the BPV are of small magnitude (a few bars at maximum) and damped, and the system is stable. However, for low permeability material (Pierre de Lens) these oscillations are much larger (several bars or even ten or more bars) and the damping is less strong. For high carbon dioxide flow rates (1 cc/min or above), this leads to the phenomenon described in the previous section (gasometer and separator leakage, and often condensation trap plug ejection). The only solution to this effect is to limit the supercritical carbon dioxide flow rate for low permeability materials.

Meniscus and sensor in the separator

We have faced numerous issues with the triggering of the sensor controlling the valve at the bottom of the separator. The first issue is that if there is too much water suddenly arriving in the separator, the meniscus is going up too fast and the sensor is not triggering. As a result the sensor cannot trigger anymore because the meniscus is already above it. As the quantity of water which is added to the separator before the experiment is not precisely known, this forbids to measure the water quantity. The worst case scenario is that if the water flow rate is too strong and the problem not solved (during an overnight experiment for example), water fills not only the separator, but also the condensation trap and even goes within the gasometer. To solve this issue, a small polymer foam was added in the separator. It floats at the meniscus and thus increases the spatial range of detection of the sensor (figure 13.6).

The other problem is the non desired trigger of the valve. It happened during hot days (during summer the room temperature reached easily 30 ℃). As a result, there was evaporation of the water filling the separator and some condensation occurred on the glass of the separator. These condensation droplets could sometimes trigger the sensor if they condensed right on it. As a result, the valve is opened but as it is not triggered by a moving meniscus but rather by an immobile droplet, the valve could not close. The consequence is a loss of water and gas measurement. Indeed, as the valve at the bottom of the separator does not close, a part of the gas follows this path and is not measured

Precipitation in the gasometer

The pure water used in the gasometer and in the experiments is not completely deionized. It is made by use of an ion exchange resin. As a result there are still ions in it (the resistance is 18 MΩ). It has no real consequence on the experimentation because of the very low amount of ions, but it has a consequence on the gasometer. Indeed, there was a lot of precipitation of calcite in the gasometer because of the continuous supply of carbon dioxide on the water filling it. This however had no real consequence apart a slight modification of the mass and the moment of the gasometer wheel.

Purpose

The main target of the presented set-up is to measure the evolution of intrinsic and relative permeability of a rock core sample saturated with a saline solution and submitted to a flow of carbon dioxide. Additionally, the measure of the axial deformation of the sample with a LVDT gives information on the stress acting on the rock sample, which can be related through poromechanics to the crystallization of salt and to capillary pressure. The set-up is designed to allow different types of drying scenario: it is possible to inject prepartitioned fluids (thanks to the mixing bottle (section 11.1.2.2)), to inject simultaneously brine and carbon dioxide to simulate the capillary currents, or to alternate brine and carbon dioxide injection in order to simulate resaturation of the aquifer after a stop in the injection for example. However, because of several issues developed here and because of the lack of time, these scenarii were not experimented.

Issues and solutions

Although the experimental set-up was designed to dry-out rock cores saturated by saline solution by injection of supercritical carbon dioxide, we faced several issues in carrying-out these experiments. Besides the issues described in section 13.4, the main problem during the drying is the precipitation of salt out of the rock core sample. Inox is a water-wet material and thus after the desaturation of the rock core by the flow of carbon dioxide and during the evaporation, droplets of brine can stay within the tubing after the triaxial cell and the rock core. Evaporation and precipitation then occurs in these droplets. It has often no consequence, but at some accumulation points, the precipitation leads to a clogging of the path followed by the carbon dioxide. Clogging of the tubing is easily seen on the control panel during the experiment, as P6 pressure increases while BPV pressure stays at the instructed value and no flow is coming out of the system. As the injection continues, P6 increases until it is high enough to break the salt plug. The consequences of this are dramatic. A very strong flux of supercritical carbon dioxide arrives at the BPV, as the pressure gradient between P6 and the BPV can reach 10 or more bars. This strong flux of supercritical CO 2 added to the phase change at the exit of the BPV make it oscillate with oscillation amplitude of 20 or even 30 bars. It causes leakage in the separator, the gasometer and ejection of the plug of the condensation trap. Finally, because of the strong oscillations, the seal of the BPV is under too much stress and can be deteriorated.

Another annoying precipitation occurs at the mouth of the tubing leading to P10. Clogging of this tubing does not lead to disastrous consequences but disconnects the pressure sensor from the main system giving erroneous pressure measurement.

In order to avoid this unwanted precipitation we needed to determine the main accumulation points. Two critical parts were identified: the T junction between the main path and P6 sensor, and the BPV seal. Solution for the BPV seal was to increase its diameter from 0.1mm to 1mm. The solution for the T junction was to replace it by a cross. By this way, salt and brine accumulation does not occur in the main pathway but just outside of the T-junction as represented in red in figure 14.1.

Finally, the protocol was also modified. Initially, saturation, desaturation and drying-out were made in a single step at working temperature and pressure. We then decided to uncouple the saturation/desaturation step and the drying step. The first one is made at atmospheric pressure and ambient temperature in order to avoid any parasitic precipitation. Then the rock core is dismounted and the whole system is flushed with pure water in order to remove all salt from the tubings. In particular, brine filling the tubing to the sensor P6 is replaced by pure water. Finally, the rock core sample is reinserted and the system put at working temperature and pressure (60 ℃ and 100 bar).

Addendum

Between the first submission of this dissertation and the defense, an experiment of dryingout was successfully carried-out using the improvements described in the previous section. The core used was the same core of Grès des Vosges which was used in the relative permeability experiments.

The experiment was carried-out with gaseous carbon dioxide instead of supercritical carbon dioxide. Indeed, the molar fraction of water in the carbon dioxide is higher in the gaseous state than in the supercritical state [START_REF] Spycher | CO2-H2O mixtures in the geological sequestration of CO2. I. assessment and calculation of mutual solubilities from 12 to 100C and up to 600 bar[END_REF], which allows a faster evaporation and thus allows to reduce the experiment time. The change between gaseous and supercritical carbon dioxide does not affect the behavior or the fluids in the core as shown in the chapter 13: relative permeabilities of supercritical and gaseous carbon dioxide are not different.

As a result the back pressure valve was set at 55 bar, the temperature 60℃ and the carbon dioxide flow rate at 5 cc/min. The rock core was initially saturated with a halite solution of 150g/L. The experiment lasted 8 days of carbon dioxide injection. The initial brine saturation after the desaturation with carbon dioxide was S 0 l =0.425. The evolution of the saturation in the core is directly obtained with the carbon dioxide flow rate (similarly to equation (9.44)):

S l (t) = S 0 l - V evap (t) V pores (14.1a) S l (t) = S 0 l -q CO2 ∆t 1 V pores M H 2 O M CO 2 y H 2 O ρ CO 2 ρ H 2 O (14.1b)
with V evap (t) the quantity of water evaporated after an injection of duration t, V pores the volume of pores of the considered rock core, q CO2 the carbon dioxide flow rate, ∆t the elapsed time since the beginning of the injection, M CO 2 and M H 2 O the molar mass of respectively carbon dioxide and water, y

H 2 O = n H 2 O
n H 2 O +n CO2 the molar fraction of water vapor in the carbon dioxide (n being the molar quantity of both species, and y H2O ≪ 1), ρ CO 2 and ρ H 2 O the density respectively of carbon dioxide and brine. Thanks to equation 14.1b, we can rescale the evolution of the pressure difference between the inlet and the outlet (P10 and P6 in figure 11.2) as a function of the water saturation. We then obtain figure 14.2.

We can see on figure 14.2 that the pressure difference is decreasing at first. It corresponds to the diminution of the water saturation and thus on the increase of the relative permeability of the carbon dioxide. However, at a saturation of S p l =0.18, there is a slight increase of the pressure difference. If we consider the supersaturation of the solution we obtain (from equation (9. The initial supersaturation corresponding to a halite concentration of 150 g/L is ap- As a result, if we consider that the precipitation occurs at a supersaturation equal to 1 (considering that the heterogeneous nucleation lowers sufficiently the nucleation threshold), we can calculate from equation (14.2) that precipitation of salt will occur at a water saturation of S p l ≈ 0.16.

The increase in the pressure difference corresponds then approximatively to the beginning of the precipitation of salts. We can then conclude that we observe the clogging consecutive to this precipitation. Indeed, if we consider the final pressure difference (∆P = 0.075 bar) we then obtain an intrinsic permeability of 90mD. The clogging has then reduced the intrinsic permeability of 40%.

Unfortunately, the study of the deformation did not give any notable result. It seems that the mechanical impact of this precipitation is negligible, this observation being confirmed by the study of the core by X-ray microtomography which did not show any fracture linked to the precipitation of salt. Nevertheless, this experiment should be backed up with others, changing the concentration, the flow rate, the temperature and the pressure in order to study the clogging effect and the condition of apparition of possible fracture linked to the precipitation of salt.

General conclusion

Salt crystallization during CCS is a major issue. The two main consequences of the precipitation, clogging and crystallization pressure cannot be overlooked. The first one decreases permeability sometimes by several orders of magnitude while the second can have potentially dramatic effects such as opening leakage pathways and reactivating faults. The positive effect of salt crystallization on injectivity is to fracture the rock and open new percolation paths, thus rising the permeability in the aquifer. In this work we have presented a comprehensive study of crystallization pressure in the case of the drying of a rock by a flow of carbon dioxide. There are mainly two types of precipitation corresponding to the residual water pattern after CO 2 percolation. The first one is the precipitation within the so-called blob water which is able to create very high local crystallization pressure and occurs at long time-scales. The second is the precipitation in the brine trapped by capillary forces in the corner of the pores followed by the principal flow of carbon dioxide and which occurs at short time-scales. These corners are subjected to the fastest precipitation and thus create crystallization pressure in a very short period, the intensity of the stress depending on the size of the corner in which brine is trapped. All results show that crystallization pressure in the case of CCS is a phenomenon which cannot be neglected. Several attempts of modelling the mechanical behavior of aquifers subjected to carbon dioxide injection, mainly by the numerical coupling of a reactive transport code (TOUGH/TOUGHREACT) and a mechanical code like FLAC3D [START_REF] Rutqvist | Thermohydromechanics of partially saturated geological media: governing equations and formulation for four finite element models[END_REF][START_REF] Rutqvist | A modeling approach for analysis of coupled multiphase fluid flow, heat transfer and deformation in fractured porous media[END_REF][START_REF] Rutqvist | A comparative simulation study of coupled THM processes and their effect on fractured rock permeability around nuclear waste repositories[END_REF] simply overlook this effect missing a whole part of the possible mechanical response of the aquifer to carbon dioxide injection. However, as presented in part III, a complete description of crystallization pressure during CCS seems particularly complicated. The coupling with experiments is then mandatory. Reactive percolation experiments under geotechnical conditions will provide inestimable results on the behavior of a rock core sample and if they are coupled with on field results and microscopical observation like the study of the crystallization pattern with SEM, they should provide a solid experimental framework to model more precisely the behavior of a reservoir rock subjected to carbon dioxide injection. Unfortunately, reactive percolation experiments which were supposed to be an important part of this work could not be proceeded in a satisfactory manner because of the novelty of the experimental set-up and the long time spent to design it. However, interesting measurements of the relative permeabilities of a sandstone with the pairs supercritical CO 2 /water and gaseous CO 2 /water have been carried out, along with measurements of capillary pressure. The results show little difference of the influence of the state of the carbon dioxide on the multiphase behavior of the considered rock, bringing interesting insights for the underground unsaturated numerical simulations. Moreover, a first attempt of a drying-out experiment has been carried out shortly before the PhD defense. This experiment shows clearly the clogging effect of the salt crystallization. This bodes well for the future experiments on this subject as almost all experimental issues have been solved. The set-up is now ready to carry-out all sort of experiments and I hope that in a close future, experimental results will be obtained and compared to the different modellings presented in this work.

A.2. FUGACITIES Finally, the complete set of parameters defining the model are:

β (0) M X , β (1) M X , β (2)
M X , C φ M X for each cation-anion pair MX, θ ij for each cation-cation and anion-anion pair, ψ ijk for each cation-cation-anion and anion-anion-cation triplet, λ ni for ion-neutral pairs. The interaction between three ions of same sign have been neglected because they are very unlikely to happen.

A.2 Fugacities

The non-ideal characteristics of the CO 2 rich phase is described by the fugacity coefficients of water vapor and carbon dioxide in the mixture. Correspondingly to the Pitzer's treatment of the activity coefficient, fugacities can be obtained by a virial development of the compressibility factor Z = P V RT . If the gases are considered as ideal, the compressibility factor reduces to the mole quantity of gas considered. Considering the non-ideality of the mixture, the virial development of Z gives: .12) or in terms of pressure: .13) with B ′ = B/RT and C ′ = (C -B 2 )/(RT ) 2 . Considering a mixture we have:

Z = 1 + B/V + C/V 2 + D/V 3 . . . (A
Z = 1 + B ′ p + B ′ p 2 (A
B = i j y i y j B ij C = i j k y i y j y k C ijk (A.14)
with i, j, k, the components of the mixture and B ij and C ijk , the interaction parameters between respectively 2 and 3 components. Coefficients of the virial development are function only of the temperature and are tabulated (see [START_REF] Spycher | Fugacity coefficients for H2 CO2 CH4 H2O, and of H2O-CO2-CH4 mixtures: a virial equation treatment for moderate pressures and BIBLIOGRAPHY temperatures applicable to hydrothermal boiling[END_REF]).

Finally, fugacity coefficient can be calculated as:

ln Φ = p 0 Z -1 p dp (A.15)
Using the following dependency on temperature for the parameters of the compressibility factor B ≈ a/T + b + cT and C ≈ d/T 3 + e/T 2 + f /T + g, we obtain:

The term j σ j ∂A j ∂h i dh i corresponds to the variation of surface energy of all the faces A j impacted by the modification of the length h i .

Writing that dn (i) s = dV i υs , and noting that dV i = dh i × A i for an infinitesimal transformation, we can then rewrite the previous equation as:

df i = -(p s -p i )A i dh i + j σ j ∂A j ∂h i dh i + µ s A i dh i υ s (B.5)
Geometrically, we have the following relations for a pyramid:

A i = 1 2 j h j ∂A j ∂h i (B.6)
Using this expression of A i in the free energy we get: Considering now the global system at equilibrium, the total energy variation is zero which leads to:

df i = -(p s -p i ) 1
dF = 0 = i df i + dF l = i df i + df (i) l (B.9)
with df (i) l the liquid free energy variation associated with the growth of face i.

Combining now equations (B.1), (B.2), (B.8), (B.9) with dn s = i dn (i) s , we obtain:

0 = i dh i j ∂A j ∂h i - 1

4 . 3

 43 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.4 Schematics of the different targets for CCS (from Kaltediffusion.ch) . . . . 1.5 Relative efficiency of the different solutions . . . . . . . . . . . . . . . . . . 3.1 Illustration of the Lagrangian porosity and the deformation of the porous volume in a two-fluid case (from [20]). . . . . . . . . . . . . . . . . . . . . 3.2 Definition of the surface energy as the creation of new surfaces. . . . . . . 3.3 Sketch of the different wetting situations. . . . . . . . . . . . . . . . . . . . 3.4 Curvature between two fluids in a pore. . . . . . . . . . . . . . . . . . . . . 3.5 Evolution of capillary pressure with water saturation S w in drainage and imbibition conditions (from [40]). . . . . . . . . . . . . . . . . . . . . . . . 3.6 Determination of fluid behavior according to the comparison between capillary numbers and viscosity ratio. The grey zone represents the zones where the fluids have an undetermined (intermediate) behavior (from [35]. . . . . 3.7 Sketch of the percolation of CO 2 in a porous medium and the different kinds of residual water. The solid matrix is shown in black, water is blue and CO 2 is orange. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.1 Sillén diagram for carbon dioxide aqueous species. . . . . . . . . . . . . . . 4.2 Evolution of the Gibbs free energy of formation of a nuclei. . . . . . . . . . xv Nucleation rate of epsomite, mirabilite and halite as a function of sursaturation assuming heterogeneous nucleation (θ = 57°, T = 23℃, from [8]). . . 4.4 Modelisation of a conceptual pore (a) with the tubes-in-series model (b), from [80]. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6.1 Schematic of the aquifer during carbon dioxide injection (from [96]). . . . . 6.2 Distribution of the CO 2 saturation according to numerical simulation (from [97]). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8.1 Illustration of the Wulff shape for a cubic halite crystal. . . . . . . . . . . . 8.2 Vicinal surface of a crystal with an angle θ. . . . . . . . . . . . . . . . . . . 8.3 Overlap between the energy profile of two flat surfaces and the creation of the disjoining pressure, adapted from [19]. . . . . . . . . . . . . . . . . . . 8.4 Crystal growing in a cylindrical pore and the different curvatures (adapted from [107]). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

8. 7

 7 Halite crystal at homogeneous pressure p s confined in a pore with a square section. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8.8 Tilted halite crystal at homogeneous pressure confined in a pore with a square section. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8.9 Evolution of the equilibrium thickness δ with the supersaturation for halite at 40 ℃. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9.1 Modelling of the porous medium as a succession of spherical pores of decreasing size linked by cylindrical channels of decreasing radius. . . . . . . 9.2 Decomposition of the transmitted stress into crystallization pressure and liquid pressure for halite at 40 ℃ for the MDV capillary curve. . . . . . . . 9.3 Evolution of the transmitted stress with water saturation for different temperature for hydrohalite/halite 9.3(a) and9.3(c) and the couple anhydrite/gypsum 9.3(b) and 9.3(d) for the BBL and MDV rocks. . . . . . . . . 9.4 Evolution of the equivalent tensile stress with water saturation in a pure water condition. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9.5 Evolution of the difference in equivalent tensile stress with water saturation for different temperature for halite/hydrohalite 9.5(a) and 9.5(c) and gypsum/anhydrite 9.5(b) and 9.5(d) for BBL and MDV rocks. . . . . . . . 9.6 Profiles of liquid pressure 9.6(a), equivalent tensile stress 9.6(b), liquid 9.6(c) and crystal 9.6(d) saturations for BBL at 5mm/s at 196s, MDV at 5mm/s at 12s and MDV at 5µm/s at 14996s. . . . . . . . . . . . . . . . . . 9.7 Comparison of the three equivalent tensile stress obtained with the two macroscopic calculations. . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9.8 Presentation of the capillary trapped water in the corner subjected to evaporation by CO 2 flowing through the big pore. . . . . . . . . . . . . . . . . . 9.9 Speed of growth for a small crystal as a fonction of the radius for constant supersaturation. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9.10 Description of the evolution of the corner during evaporation and precipitation. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9.11 Algorithm used for the simulation of the behavior of the corner. . . . . . . 9.12 Pore model used for the calculation. . . . . . . . . . . . . . . . . . . . . . . 9.13 Creation of the crystallization pressure for a spherical crystal close to the pore wall. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9.14 Evolution of the water saturation and the supersaturation in the pore. . . . 9.15 Crystallization pressure corresponding to the supersaturation of figure 9.14. 9.16 Influence of the CO 2 flow rate (top) and the pore size (bottom) on the maximum crystallization pressure. . . . . . . . . . . . . . . . . . . . . . . . 9.17 Evolution of the nucleation rate and the critical radius of nucleation with time. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9.18 Variation of size for the nuclei at the moment of their creation. . . . . . . . 9.19 Evolution of the crystal sizes for all families with time. . . . . . . . . . . .

CHAPTER 1 .

 1 CONTEXT OF THE STUDY1.1 Carbon dioxide and greenhouse effectDuring the last century, anthropic emissions of greenhouse gases and particularly carbon dioxide have been drastically increasing as represented in figure1.1.

Figure 1 . 1 :

 11 Figure 1.1: Evolution of the anthropic emissions of carbon dioxide since 1800 (from [13])

Figure 1 . 2 :Figure 1 . 3 :

 1213 Figure 1.2: Schematic of the greenhouse effect (from [16])

1. 2 .

 2 Figure 1.4: Schematics of the different targets for CCS (from Kaltediffusion.ch)

Figure 1 . 5 :

 15 Figure 1.5: Relative efficiency of the different solutions

3 . 1 3 . 2

 3132 Porosities and phase saturations . . . . . . . . . . . . . . . . . 3.1.1 Definitions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.1.2 Diagenesis and the formation of natural porous media . . . . . Mass balance and fluid velocity . . . . . . . . . . . . . . . . . .

Figure 3 . 1 :

 31 Figure 3.1: Illustration of the Lagrangian porosity and the deformation of the porous volume in a two-fluid case (from [20]).

Figure 3 . 3 :

 33 Figure 3.3: Sketch of the different wetting situations.

Figure 3 . 4 :

 34 Figure 3.4: Curvature between two fluids in a pore.

Figure 3 . 6 :

 36 Figure 3.6: Determination of fluid behavior according to the comparison between capillary numbers and viscosity ratio. The grey zone represents the zones where the fluids have an undetermined (intermediate) behavior (from [35].

Figure 3 . 7 : 4 Chemical behavior of the aquifer Contents 4 . 1

 37441 Figure 3.7: Sketch of the percolation of CO 2 in a porous medium and the different kinds of residual water. The solid matrix is shown in black, water is blue and CO 2 is orange.

4. 1 .. 4 ) 4 . 1 . 1 . 2

 144112 CARBON DIOXIDE/BRINE PARTITIONING And we have the following relation: G(p, T ) = i n i µ i (Euler identity) (4Chemical Equilibria and Gibbs-Duhem relation:

CHAPTER 4 .

 4 CHEMICAL BEHAVIOR OF THE AQUIFER ∆H = ∆G + T ∆S (

CHAPTER 4 .Figure 4 . 1 :

 441 Figure 4.1: Sillén diagram for carbon dioxide aqueous species.

Figure 4 . 2 :

 42 Figure 4.2: Evolution of the Gibbs free energy of formation of a nuclei.

Figure 4 . 3 :

 43 Figure 4.3: Nucleation rate of epsomite, mirabilite and halite as a function of sursaturation assuming heterogeneous nucleation (θ = 57°, T = 23℃, from [8]).

  :

4. 2 .

 2 CHEMICAL REACTIONS INDUCED BY CARBON DIOXIDE DISSOLUTION (MINERAL TRAPPING)

Figure 4 . 4 : 5 Thermal behavior of the aquifer Contents 5 . 1 5 . 2 5 . 3 5 . 4

 44551525354 Figure 4.4: Modelisation of a conceptual pore (a) with the tubes-in-series model (b), from [80].

Figure 6 . 1 : 1

 611 Figure 6.1: Schematic of the aquifer during carbon dioxide injection (from [96]).

=

  ε and Σ = tr = Σ /3 are the volumetric dilation and the mean stress, deviatoric strain and stress, K and µ are respectively the bulk and shear moduli, while b J and N I,J are the generalized Biot coefficients and the generalized Biot coupling moduli.

Contents 8 . 1

 81 Ostwald-Freundlich equation and Wulff theorem . . . . . . . . 8.1.1 Chemical equilibrium of a crystal in solution . . . . . . . . . . 8.1.2 Ostwald-Freundlich equation for small crystallite . . . . . . . . 8.1.3 Wulff theorem and equilibrium shape . . . . . . . . . . . . . . 8.1.4 Non-flat surfaces . . . . . . . . . . . . . . . . . . . . . . . . . . 8.2 Interaction energy between two solid surfaces . . . . . . . . .

Figure 8 . 2 :

 82 Figure 8.2: Vicinal surface of a crystal with an angle θ.

Figure 8 . 3 :

 83 Figure 8.3: Overlap between the energy profile of two flat surfaces and the creation of the disjoining pressure, adapted from [19].

8. 3 In-pore growth of a single crystal and crystalliza- tion pressure 8 . 3 . 1

 3831 First case: crystallization pressure in a cylindrical poreLet us now consider a crystal in a pore in equilibrium with a supersaturation S. Let us consider the situation of a liquid-like crystal (i.e. a crystal able to grow with any possible curvature), in a cylindrical pore. This is the classical description by[START_REF] Scherer | Crystallization in pores[END_REF][START_REF] Scherer | Stress from crystallization of salt[END_REF][START_REF] Flatt | Salt damage in porous material: How high supersaturations are generated[END_REF][START_REF] Steiger | Crystal growth in porous materials -I: The crystallization pressure of large crystals[END_REF].

Figure 8 . 4 :

 84 Figure 8.4: Crystal growing in a cylindrical pore and the different curvatures (adapted from [107]).

Figure 8 . 5 :

 85 Figure 8.5: Crystal growing in a spherical pore of radius r p and invading the small entry channels of radius r s . The curvature of the crystal in the channel depends on the current supersaturation (adapted from[START_REF] Scherer | Crystallization in pores[END_REF]).

8. 3 .

 3 Figure 8.7: Halite crystal at homogeneous pressure p s confined in a pore with a square section.

Figure 8 . 9 :

 89 Figure 8.9: Evolution of the equilibrium thickness δ with the supersaturation for halite at 40 ℃.

Figure 9 . 2 :

 92 Figure 9.2: Decomposition of the transmitted stress into crystallization pressure and liquid pressure for halite at 40 ℃ for the MDV capillary curve.

Figure 9 . 3 :

 93 Figure 9.3: Evolution of the transmitted stress with water saturation for different temperature for hydrohalite/halite 9.3(a) and9.3(c) and the couple anhydrite/gypsum 9.3(b) and 9.3(d) for the BBL and MDV rocks.

Figure 9 . 4 :

 94 Figure 9.4: Evolution of the equivalent tensile stress with water saturation in a pure water condition.

Figure 9 . 5 :

 95 Figure 9.5: Evolution of the difference in equivalent tensile stress with water saturation for different temperature for halite/hydrohalite 9.5(a) and 9.5(c) and gypsum/anhydrite 9.5(b) and 9.5(d) for BBL and MDV rocks.

Figure 9 . 6 :

 96 Figure 9.6: Profiles of liquid pressure 9.6(a), equivalent tensile stress 9.6(b), liquid 9.6(c) and crystal 9.6(d) saturations for BBL at 5mm/s at 196s, MDV at 5mm/s at 12s and MDV at 5µm/s at 14996s.

Figure 9 . 7 :

 97 Figure 9.7: Comparison of the three equivalent tensile stress obtained with the two macroscopic calculations.

Figure 9 . 8 :

 98 Figure 9.8: Presentation of the capillary trapped water in the corner subjected to evaporation by CO 2 flowing through the big pore.

9. 2 .

 2 MICROSCOPIC BEHAVIOR: NUCLEATION AND STRESS CREATION AT THE PORE LEVEL can then rewrite equations (4.54) and (4.54b) as:

  r c ln S(t) 3 (9.41) with k ′ = 4kπ, f (S l ) a function of the brine saturation corresponding to the surface avail-9.2. MICROSCOPIC BEHAVIOR: NUCLEATION AND STRESS CREATION AT THE PORE LEVEL able for nucleation, and r c = 2συs

9. 2 .Figure 9 . 11 :

 2911 Figure 9.11: Algorithm used for the simulation of the behavior of the corner.

Figure 9 . 12 :

 912 Figure 9.12: Pore model used for the calculation.

. 45 )Figure 9 . 13 :

 45913 Figure 9.13: Creation of the crystallization pressure for a spherical crystal close to the pore wall.

Figure 9 . 15 :

 915 Figure 9.15: Crystallization pressure corresponding to the supersaturation of figure 9.14.

Figure 9 .

 9 Figure 9.16: Influence of the CO 2 flow rate (top) and the pore size (bottom) on the maximum crystallization pressure.

Figure 9 . 17 :

 917 Figure 9.17: Evolution of the nucleation rate and the critical radius of nucleation with time.

Figure 9 . 18 :

 918 Figure 9.18: Variation of size for the nuclei at the moment of their creation.

Figure 9 .

 9 19(b) presents a zoom-in at the beginning of the nucleation, while 9.19(a) presents the whole families.

9. 2 .

 2 MICROSCOPIC BEHAVIOR: NUCLEATION AND STRESS CREATION AT THE PORE LEVEL

Figure 9 . 19 :

 919 Figure 9.19: Evolution of the crystal sizes for all families with time.

Figure 9 . 20 :

 920 Figure 9.20: Evolution of the normalized number of crystal and crystal radius with time.

Figure 9 . 21 :

 921 Figure 9.21: Evolution of the distribution of occupied volume with crystal radius.

  Figure 9.22: Evolution of the number of nuclei with the flow rate.

Figure 9 . 23 :

 923 Figure 9.23: Upscaling to an REV of the corner.

Figure 11 . 1 : 1 :Figure 11 . 3 :

 1111113 Figure 11.1: Picture of the experimental set-up.

Figure 11 . 4 : 2 :

 1142 Figure 11.4: Reaction of sleeves with carbon dioxide. 11.4(a) is the initial sleeve. 11.4(b) and 11.4(c) show sleeves (Viton) after carbon dioxide injection.

Figure 11 . 5 :Figure 11 . 6 : 3 :

 1151163 Figure 11.5: Cross-section of the triaxial cell.

Figure 11 . 7 :

 117 Figure 11.7: Set-up for external experimentation. Note the additional tubing.

Figure 11 . 8 :

 118 Figure 11.8: STIGMA1000 pump (here ST1 and ST3).

Figure 11 . 9 :

 119 Figure 11.9: Sketch of the Back Pressure Valve.

11. 1 .Figure 11 . 10 :

 11110 Figure 11.10: Phase separator.

Figure 11 . 11 :

 1111 Figure 11.11: Representation of the considered dead volumes.

Figure 12 . 1 :

 121 Figure 12.1: Rock core of Pierre de Lens 12.1(a) and Grès des Vosges 12.1(b)

12. 2 .

 2 Figure 12.2: Temperature program for sintering.

Figure 12 . 3 :

 123 Figure 12.3: Plaster mould used for glass sintering.

Figure 12 . 4 :

 124 Figure 12.4: Sintered glass core obtained with a plaster mould.

Figure 12 . 5 :

 125 Figure 12.5: Result of the glass sintering with quartz tubes: 12.5(a) represents the final cores after trimming, and 12.5(b) represents the cores in their tubes.

Figure 12 . 6 :

 126 Figure 12.6: Schematics of the mercury intrusion porosimetry (MIP).

12. 3

 3 Mercury Intrusion Porosimetry of the different cores 12.3.1 Principle of the measurement Mercury Intrusion porosimetry (MIP) is an experimental technique which allows to determine the distribution of pore sizes in a sample of porous medium. The measurement is based on Laplace law (equation (3.19

12. 3 .

 3 Figure 12.7: MIP results for the different rock core from the core library. BB = Bois-Brûlé, MDV = Mondeville, CF1= Clos-Fontaine

Figure 12 . 8 :

 128 Figure 12.8: MIP results for Pierre de Lens and Grès des Vosges.

Figure 12 . 9 :

 129 Figure 12.9: MIP results for two samples of sintered glass.

Figure 13 . 1 : 1 :

 1311 Figure 13.1: Fit of the evolution of pressure difference after a step increase of injection flow rate.

Figure 13 . 2 :

 132 Figure 13.2: Evolution of the pressure difference with the injection flow rate for a rock core of Pierre de Lens for pure water 13.2(a) and nitrogen 13.2(b).

Figure 13 . 3 :

 133 Figure 13.3: Relative permeability of Pierre de Lens to nitrogen and pure water.

Figure 13 . 4 :

 134 Figure 13.4: Relative permeability of Grès des Vosges to supercritical 13.4(a) and gaseous 13.4(b) carbon dioxide and pure water (60 ℃ and 100 bar for supercritical carbon dioxide, 2 bar 20 ℃ for gaseous carbon dioxide).

Figure 13 . 5 :

 135 Figure 13.5: Capillary pressure curve for the Grès des Vosges sandstone.

Figure 13 . 6 :

 136 Figure 13.6: Polymer foam to improve the meniscus detection

13. 4 .

 4 Figure 13.7: Precipitation of calcite on the wheel of the gasometer because of the ions in the filling water.

Figure 14 . 1 :

 141 Figure 14.1: Modification of the T junction for salt precipitation.

  4)):S(S l ) = S(S 0 l )

Figure 14 . 2 :

 142 Figure 14.2: Evolution of the pressure difference during the evaporation with brine saturation. The red line is the Savitzki-Golay smoothing on 20 000 points of window

4.2 Chemical reactions induced by carbon dioxide dissolution (mineral trapping) . . . . . . . . . . . . . . . . . . . . . . . . . .

	4.2.1 Carbon dioxide dissolution and pH evolution . . . . . . . . . .
	4.2.2 Mineral reactions and carbon dioxide trapping . . . . . . . . .
	4.2.3 Kinetics of mineral reactions: nucleation and crystal growth . .
	4.2.4 Porosity and permeability variations induced from dissolu-
	tion/precipitation . . . . . . . . . . . . . . . . . . . . . . . . .

  4.1. CARBON DIOXIDE/BRINE PARTITIONING K H (T )Φ CO 2 p CO 2 = γ CO 2 x CO 2 (4.25) K H is the Henry's constant at temperature T , p CO2 the partial pressure of CO 2 (in the case of a gas at pressure p tot , if y CO 2 is the mole fraction of gaseous carbon dioxide, p CO 2 = y CO 2 p tot (Dalton's law)); Φ CO 2 is the fugacity coefficient, γ CO 2 the activity coefficient of aqueous CO 2 , and x is the molar fraction in the aqueous solution. The quantity of dissolved CO 2 can also be measured by mass or mole concentration, or molality.

  σ

	CHAPTER 4. CHEMICAL BEHAVIOR OF THE AQUIFER
	Ion		• a
	Cl -	3	• A
	Na +	4	• A
	Ca 2+	6	• A
	SO 4	2-	4	• A
	Table 4.1: Values of the mean electrostatic radius for different solute species.
	with σ[x] = 3 x 3 1 + x -1 1+x -2 ln(1 + x) . The value of and is represented in table 4.1 [58].	• a depends on the considered ion
					• aB γ	√	I -BI 2	(4.35)

  Mg 5 Al 2 Si 3 O 10 (OH) 8 Al 2 Si 2 O 5 (OH) 4

	KAlSi 3 O 8 K-f eldspar	+ 5 2 M g-chlorite	+ 25 2 CO 2 (aq)
		= KAlSi 3 O 10 (OH) 2 muscovite	+ 3 2 kaolinite	magnesite 2 MgCO 3 + 25	+
					9 2 SiO 2 qtz/chal + 6 H 2 O (4.41)

4.2. CHEMICAL REACTIONS INDUCED BY CARBON DIOXIDE DISSOLUTION (MINERAL TRAPPING) are represented by Mg-Chlorite, Fe-Mg carbonates are represented by magnesite):

Table 4 .

 4

	3

2: Shape factor for different regular geometries.

8.3 In-pore growth of a single crystal and crystallization pres- sure

  8.2.1 Equilibrium of a crystal surface close to the pore wall . . . . . 8.2.2 Expression of the interaction energy . . . . . . . . . . . . . . .

. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8.3.1 First case: crystallization pressure in a cylindrical pore . . . . 8.3.2 Second case: crystallization pressure in a spherical pore with small entry channels . . . . . . . . . . . . . . . . . . . . . . . . 8.3.3 Crystallization pressure and Wulff shape . . . . . . . . . . . . . 8.3.4 Thickness δ and interaction forces . . . . . . . . . . . . . . . . 8.3.5 Effectively transmitted stress . . . . . . . . . . . . . . . . . . . 8.4 Final remark . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

8.

1 Ostwald-Freundlich equation and Wulff theorem 8.1.1 Chemical equilibrium of a crystal in solution

  Let us consider the dissolution of an hydrated saltA ν A B ν B . . . • ν w H 2 O, A, Bare the ions composing the crystal, H 2 O is also denoted as w and ν are the stoichiometric coefficients.

2 Microscopic behavior: nucleation and stress creation at the pore level

	9.2.4 Nucleation and crystal growth . . . . . . . . . . . . . . . . . .
	9.2.5 Upscaling and poromechanics . . . . . . . . . . . . . . . . . . .

Contents 9.1 Macroscopic behavior: poromechanics at the REV scale . . . 9.1.1 Modelling at constant concentration . . . . . . . . . . . . . . . 9.1.2 Simulation of a REV submitted to carbon dioxide evaporation 9.1.3 Comparison of the two modellings . . . . . . . . . . . . . . . . 9.. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9.2.1 Evolution of the salt quantity in the corner . . . . . . . . . . . 9.2.2 Algorithm . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9.2.3 Results of the simulation . . . . . . . . . . . . . . . . . . . . . 9.3 Conclusion on the modellings . . . . . . . . . . . . . . . . . . . OF CCS

Table 11 . 5 :

 115 Technical specifications of the back pressure valve (and the micropump).

	Technical specifications
	Volume	1L
	Flow rate	0.001 cc/min to 360 cc/min
	Precision	0.1 % of the consigne
	Maximum working pressure	300 bar
	Precision on the pressure	0.05 % Full scale
	Material	Inox 316
	Working temperature	5 to 40 ℃
	Table 11.4: Technical specifications of the STIGMA 1000 pump.

Table 13 .

 13 2: Flow rates used in the measurement of relative permeability of Pierre de Lens and the corresponding water saturation. The gas used here is nitrogen.

		min)	S l
	10	0.01	0.1518
	10	0.01	0.4688
	5	0.01	0.6577
	5	0.01	0.6726
	0.5	0.01	0.7976
	0.1	0.5	0.9107
	0.1	0.5	0.9464
	0.05	0.5	0.9464
	Q G (cc/min) Q l (cc/min)	S l
	1	1	0.8636
	1	0.1	0.7061
	5	0.1	0.4415
	5	0.01	0.4350
	10	0	0.256

Table 13 .

 13

3: Flow rates used in the measurement of relative permeability of Grès des Vosges and the corresponding water saturation. The gas used here is carbon dioxide

13.2.1.1 Pierre de Lens

Table 13 .

 13 13.2. RELATIVE PERMEABILITYQ G (cc/min) Q l (cc/min) Weighted S l Initial S l 4: Flow rates used in the measurement of relative permeability with supercritical carbon dioxide of Grès des Vosges and the corresponding weighted water saturation.

	0.05	0.5	0.6200	1
	0.02	0.5	0.6382	1
	0.07	0.5	0.5907	1
	0.1	0.02	0.5790	1
	0.1	0.02	0.4022	0
	5	0.02	0.1077	0

  . INTRINSIC AND RELATIVE PERMEABILITY MEASUREMENT Weighted S l Corrected S l

	0.6200	0.7111
	0.6382	0.8236
	0.5907	0.73505
	0.5790	0.6037
	0.4022	0.4480
	0.1078	0.2661
	Table 13.5: Corrected water saturation.

Table 13 .

 13 and coefficients extracted from the relative permeability curve 13.3. CAPILLARY PRESSURE MEASUREMENT Pressure difference (bar) Extracted water (g) Calculated S l Percolation 6: Values of capillary pressure and corresponding liquid saturation for the Grès des Vosges.

	0	0	0	no
	0.06	3.6	0.79	no
	0.13	14.26	0.50	yes
	0.22	14.5	0.49	yes
	0.29	16.0	0.42	yes

Table D .

 D 2 h j (p s -p i ) + σ j + 1: Correlation coefficients for the capillary pressure curve (p cap in MPa).

	Temperature (℃)		a(T).10 9 [Pa -1 ]	b(T)
			Hydrohalite
	-10		2.25258		1.04933
			Halite
	-10		3.23255 Mondeville		1.46624
	10	a	2.1217	1.71665	1.55749
	30	b	1.65655	-20.72356	1.60692
	40	c	1.52922	45.28661	1.62281
	50	d	1.44492	-41.64102	1.6345
	80	e	1.36284	12.35925	1.65011
			Bois-Brûlé
		a	Anhydrite	1.40432
	-10	b	11.8	-15.93202	-4.1987
	10	c	9.48357	130.04237	-4.28815
	30	d	8.36036	-537.98057	-4.48369
	40	e	8.0048	987.70878	-4.60369
	50	f	7.77385	-824.27113	-4.73448
	80	g	7.31559	255.97455	-5.17445
			Gypsum
	-10		10.4		-4.6545
	10		8.16458		-4.59564
	30		7.10101		-4.60153
	40		6.7646		-4.61502
	µ s 2υ s Table C.1: Correlation coefficients for variation of the equilibrium constant with liquid h j -α µ (l) α ν α i dV i υ s (B.10) pressure.

leading to (with dV

i = A i dh i ):

Other formulations of relative permeabilities can be found in[START_REF] Dana | Experimental study of two-phase flow in three sandstones. I. measuring relative permeabilities during two-phase steady-state experiments[END_REF]. However, these expressions are less used than the van Genuchten/Mualem relation

[START_REF] Appelo | Multicomponent diffusion modeling in clay systems with application to the diffusion of tritium, iodode, and sodium in opalinus clay[END_REF] suggests that the activity (see 4.1.4.1) should be used instead of the molality for ionic strength greater than 1 or 2

For experimental data of carbon dioxide/brine interfacial energy, see for example[START_REF] Chiquet | CO2/water interfacial tensions under pressure and temperature conditions of CO2 geological storage[END_REF] and[START_REF] Espinoza | Water-CO2-mineral systems: Interfacial tension, contact angle, and diffusions-implication to CO2 geological storage[END_REF].

See[START_REF] Bachu | Effects of in-situ conditions on relative permeability characteristics of CO2-brine systems[END_REF] for experimental measures of interfacial tensions, capillary pressure curves and relative permeabilities in rocks from the Alberta Basin. Other relative permeabilities and capillary pressure curves measurement can be found in[START_REF] Dana | Experimental study of two-phase flow in three sandstones. I. measuring relative permeabilities during two-phase steady-state experiments[END_REF][START_REF] Bennion | Relative permeability characteristics for supercritical CO2 displacing water in a variety of potential sequestration zones in the Western Canada sedimentary basin[END_REF][START_REF] Dana | Experimental study of two-phase flow in three sandstones. II. capillary curve measurement and relative permeability pore space capillary models[END_REF].

See[START_REF] Schneider | Wettablity determination of core samples through visual rock and fluid imaging during fluid injection[END_REF] for an experimental observation of the influence of wettability on hydrodynamic regimes

Viscous fingering is also called DLA: diffusion-limited aggregation. DLA is a process defining the formation of aggregates (fingers) by diffusion process of particles. DLA refers to a larger field of phenomena than just hydrodynamic regimes describing for example the electrodeposition of minerals or dendritic growth of polymers. Anti-DLA behavior is characterised by a stable front of the injected fluid, while capillary fingering presents similar patterns to DLA but dominated by capillary forces i.e. wettability.

For more elements on thermodynamics see[START_REF] Coussy | Mechanics and Physics of Porous Solids[END_REF] or[START_REF] Gaskell | Introduction to the thermodynamic of materials[END_REF]

For the effect of negative pressure on activity, see[START_REF] Mercury | Thermodynamic properties of solutions in metastable systems under negative or positive pressures[END_REF]

If we write u AA and u BB the energy of interaction between molecules of species A and B in the pure constituents, then u AB corresponding to the interaction between molecule A and B in the mixture has to be of the same magnitude. If we define w = 2u ABu AAu BB , then for an ideal mixture w is equal to 0, which is equivalent to ∆H = 0.

Supercritical carbon dioxide can also react with concrete (concrete carbonation), see[START_REF] Garcia-Gonzalez | Modification of composition and microstructure of portland cement pastes as a results of natural and supercritical carbonation procedures[END_REF].

The theory mostly comes from the condensation of vapor to liquid which explains the use of a contact angle. Of course in the case of crystals, there is no real contact angle. However, a value of θ can be obtained by Young-Dupré law (3.17) knowing the different surface energies.

See[START_REF] Royne | Rim formation on crystal faces growing in confinement[END_REF] for the effect of the diffusion of salt in the film on the shape of confined crystals.

(a) (b)Figure 9.14: Evolution of the water saturation and the supersaturation in the pore.

(a) (b)Figure 9.24: Results of the upscaling: equivalent macroscopic stress and strain.

Remerciements

List of Tables

Chapter 10

General purpose of the experiments

In chapter 9 we have presented three modellings of the crystallization in pores allowing to calculate the stress and strain of a REV subjected to carbon dioxide drying. These models are based on several fundamental hypotheses such as the residual water pattern, the crystallization processes. . . In order to verify these hypotheses and to obtain experimental data on the drying of a REV saturated with brine, this work includes an important experimental part. The initial idea was to study the evolution of intrinsic and relative permeability along with deformation of the sample during injection of supercritical carbon dioxide within a rock core saturated with brine under geotechnical conditions (axial and deviatoric stresses and controlled temperature). To proceed with these experiments, we developed with Sanchez Technologies the specifications and the characteristics of a new experimental prototype of reactive percolation. This set-up was intended as an important equipment for the laboratory and will be used after the end of this PhD work. The review of the development and the setting-up of the prototype is reported in table 10.1.

The system which is more precisely described in the next chapter 11

Appendices

Appendix A

Pitzer model and fugacity coefficient calculation A.1 Pitzer's model for activity coefficients

In the following we consider an aqueous solution containing several ions and non charged species (like aqueous CO 2 ). The Pitzer's model allows to calculate the activity of the different species in solution up to very high concentrations (I > 5 or more).

Let us start with the activity of the water (the solvent):

where φ is the osmotic coefficient and M w = 1.801528.10 -2 kg.mol -1 is the molar mass of water. The osmotic coefficient characterises the deviation of a solvent from ideal behavior.

i represents all solutes (anions, cations and neutral).

The calculation of the solubility equilibrium requires the knowledge of the activity and the osmotic coefficient as a function of composition, pressure and temperature. In order to do so, we need to express the excess Gibbs energy of the solution defined in (4.13). One way is to develop this excess Gibbs energy as a extended Debye-Hückel limiting law as presented section 4.1.4.1, and add a virial expansion in order to represent the short range ionic interactions. Only second and third virial coefficients are considered and the ionic strength dependence of the third virial coefficient is neglected. Expressions of activity and osmotic coefficients are then obtained by appropriates derivatives of the Gibbs free energy equation. The virial coefficients cannot be determined theoretically and have thus to be measured experimentally. CALCULATION Let us now express the virial expansion of the excess Gibbs energy:

W w is the mass of water in kilograms, m i , m l . . . are the molalities of the ions and I is the ionic strength of the solution. The first term represents the Debye-Hückel limiting law, λ ij (I) represents the short range interactions in the presence of solvent between solute particles i and j. This binary interaction parameter of second virial coefficient depends on ionic strength, on the considered species, and on the temperature and pressure. µ ijk is the third virial coefficient and represents interactions between three particles. In theory the third virial coefficient depends on the ionic strength but it is ususal to neglect this dependency. The parameters λ and µ are not measurable directly. We must then change the parameters to obtain a set which is experimentally measurable. We then define:

with M and N ions of same sign, and M and X ions of opposed sign. z X is the charge of ion X. The value of B can be obtained experimentally with pure electrolyte solutions and θ exists only for mixtures of electrolytes.

The excess Gibbs energy can then be rewritten:

The functions Z and C M X are defined as:

The virial coefficients depends on the ionic strength of the solution:

with the following definitions for g and g ′ :

when either cation M or anion X is monovalent, α M X = 2.0 kg 1/2 .mol -1/2 . For 2-2 or higher valence pairs, α M X = 1.4. In most cases β (2) is equal to zero for monovalent pairs whereas it is usually non-zero for 2-2 electrolytes.

Similarly we have:

where E θ ij (I) and I E θ ′ ij (I) are functions only of ionic strength and electrolyte pair type. The constant θ ij is a parameter of the model.

The second virial coefficient λ ni represent the interaction between ions and neutral species and is assumed constant. The third virial coefficients C φ M X and ψ ijk are also assumed to be independent of ionic strength. C φ M X is a single electrolyte parameter while ψ ijk refers to mixed electrolytes (i, j, k do not necessarily refers to ions of same sign). By convention, the first two subscripts are chosen as the like-charged ions, and the last one is chosen as the oppositely-charged ion). CALCULATION

neglecting the higher order factors (the values of a, b, c, d, e, d are not the same as in the expressions of B and C).

The same reasoning can be applied for a mixture. The fugacity coefficients are then:

leading for a 2 component mixture:

with a similar equation for Φ 2 .

Another way to model the fugacities is to use the Redlich-Kwong equation (A. [START_REF] Coussy | Mechanics and Physics of Porous Solids[END_REF]) instead of the virial development of the compressibility factor.

where a = k 0 + k 1 T and b are parameters derived from critical constraints. For mixtures, we have:

which is similar to the description of B and C in the virial development.

The fugacity coefficient is then:

A simplification can be made considering that the proportion of water in the carbon dioxide is very small. This allows to write y H 2 O = 0 and y CO 2 = 1, leading to a = a CO 2 and b = b CO 2 .

A.2. FUGACITIES

As we have seen, the fugacity coefficients depend directly on the equation of state used to describe the behavior of the gas/supercritical mixture. For example [START_REF] Duan | An improved model calculating CO2 solubility in pure water and aqueous NaCl solutions from 273K to 533K and from 0 to 2000bar[END_REF] used another equation of state and thus obtained another relation for the fugacity coefficients. However, all these relations present a good accuracy in describing the behavior of the mixture.

Appendix B Wulff construction and Wulff theorem

The Wulff construction is the construction of the crystal shape considering their surface energies. The construction is decomposed as follows:

• choose a random point which will be the geometrical center of the future crystal

APPENDIX B. WULFF CONSTRUCTION AND WULFF THEOREM

This is a summary of the 2D Wulff construction. In 3D, the construction is the same but with planes instead of edges. In order to proceed with this construction we need to get the Wulff length of all faces. To do so we will consider a crystal in its Wulff shape in equilibrium with a solution of infinite volume, such as there is a different pressure on its faces. We also consider that the pressure inside the crystal is homogeneous. Let us now consider an infinitesimal transformation which makes all faces i displace of a distance dh i from their earlier position as showed in figure B.2. The Helmholtz energy variation following this transformation is:

considering a constant temperature and an infinite volume. dn (l) α is the mole quantity of ion α which participated to the growth of the crystal. We can then write:

with dn s the mole quantity variation of the crystal and ν α the stoichiometric coefficient of the ion α in the precipitation equation.

Let us now consider the variation of the Helmholtz energy of the crystal. We divide the crystal in pyramids having as base one face of the crystal and as apex the center of the crystal. By additivity we can write:

with df i the energy variation of the considered pyramid. We can then express this variation as follows:

The chemical equilibrium gives: α µ (l) α ν α = µ s (equation (8.3)) which leads to:

If the crystal is subjected to an homogeneous pressure, we have p i = p l for all faces. In this case, a solution is to consider that all the terms in the sum are nil. We then get:

which is the Wulff theorem.

Let us now consider the case where each face is submitted to a different pressure. We can then divide the sum over i in two parts and use equation (B.6) leading to:

We can then switch the two summation signs to get:

∂h i is the exact differential of A j which leads to:

The index of the second summation can be replaced by an index i. We can then extend the summation to the whole expression. The volume of a pyramid having for base the i face is:

We then have two ways to calculate this volume:

This leads to:

Once again one solution of this equation is to make all the terms of the summation nil:

The Wulff theorem can then be extended to a case where all faces are subjected to a different pressure. This is particularly important because it allows to use this relation to consider the crystallization pressure of a crystal.

Appendix C Correlation of the equilibrium constant

Values of equilibrium constants have been obtained by the THERMO-ZNS code [START_REF] Lassin | Geochemistry of unsaturated soil systems: Aqueous speciation and solubility of minerals and gases in capillary solutions[END_REF]. This code gives (among other things) the value of the thermodynamic constant for the chosen mineral at the considered temperature and pressure. In order to calculate easily the evolution of the equilibrium constants with temperature and pressure, we have used the punctual data of THERMO-ZNS and fit them on the considered domain. If we consider the variation of K s with the liquid pressure we obtain:

In the considered range of liquid pressure, the molar volume of reaction is considered as constant. Accordingly, integration of equation (C.1) leads to:

with p * a reference pressure.

Variation of the equilibrium constant is then linear with the liquid pressure. The fit of data from THERMO-ZNS gives a linear correlation between the decimal logarithm of We then have by identification:

In table C.1 are referenced the correlation coefficients for the minerals considered in this work (halite, hydrohalite, anhydrite and gypsum).

Appendix D Capillary pressure curve estimation from mercury intrusion porosimetry

Measures of capillary pressure below the residual saturation are quite difficult. There are even more difficult if we consider the system brine/supercritical CO 2 at high temperature and pressure. The van Genuchten-Mualem correlation is also not well defined below this residual saturation. On way to estimate the capillary pressure curve is to calculate it from the Mercury Intrusion Porosimetry results as presented in chapter 13 and following the method presented in [START_REF] Fabbri | Indirect measurement of the ice content curve of partially frozen cement based materials[END_REF].

The idea is to compare the in-pore interface energy of the two systems in order to link their capillary pressure for a given water saturation:

with θ the wetting angle of the mercury on the pore wall, σ M and p M cap are respectively the mercury surface energy and capillary pressure, while δ is a function of the thickness of the adsorbed brine layer e and is equal to:

Using laplace law (3.19) we obtain the pore radius of penetration of the mercury: In

RESULTATS DES ANALYSES

-DETERMINATION DES MASSES VOLUMIQUES ET DE LA POROSITE :

Les résultats sont présentés dans le tableau ci-dessous.

Définitions des paramètres présentés :

-Masse volumique réelle : correspond à la densité de la matière solide + la densité de tous les vides fermés dans l'échantillon.

-Masse volumique apparente : correspond à la densité de la matière solide + la densité de tous les vides présents dans l'échantillon.

-Porosité ouverte : correspond au volume des vides ouverts (accessibles par la surface de l'échantillon) présents dans l'échantillon, rapportés en pourcentage du volume total de l'échantillon.

Référence des échantillons

Masse volumique réelle

En g / cm

-ANALYSE GRANULOMETRIQUE EN VOIE SECHE :

Pour la réalisation de l'analyse granulométrique, l'échantillon de grès a été déconsolidé par broyage ménagé au mortier manuel afin de ne pas réduire la taille des grains qui composent la roche.

Le « sable » obtenu est ensuite analysé par tamisage en voie sèche.

Les résultats de cette analyse sont présentés dans le tableau et graphique ci-dessous :