Benjamin Malvoisin
email: benjamin.malvoisin@univ-grenoble-alpes.fr

Lukas P Baumgartner

Mineral Dissolution and Precipitation Under Stress: Model Formulation and Application to Metamorphic Reactions

model of reaction by dissolution-precipitation under stress providing new creep laws for pressure solution and replacement reaction • in the presence of fluid, dissolution-precipitation creep is the dominant deformation mechanism in the Earth's crust • grain shape preferred orientation can only develop near thermodynamic equilibrium

Introduction

Stylolites, pressure shadows or flattened pebbles are common observations in the first 15 km of the Earth's crust [START_REF] Gratier | The Role of Pressure Solution Creep in the Ductility of the Earth's Upper Crust[END_REF]. They are markers of anisotropic mineral dissolution and mass transfer. Dissolution occurs when fluid is undersaturated with respect to equilibrium. As many mineral solubilities increase with pressure [START_REF] Eugster | Mineral solubilities and speciation in supercritical metamorphic fluids[END_REF][START_REF] Dolejš | Thermodynamic Model for Mineral Solubility in Aqueous Fluids: Theory, Calibration and Application to Model Fluid-Flow Systems[END_REF], grain-to-grain surfaces under high normal stress can locally dissolve. The dissolved aqueous species then diffuse and precipitate at surfaces under low normal stress. This process of mass transfer called pressure solution induces compaction by creep at large scale [START_REF] Durney | Solution-transfer, an Important Geological Deformation Mechanism[END_REF]. In the presence of fluid and at temperature below 500 • C, pressure solution is thought to be the dominant deformation mechanism in Earth's crust [START_REF] Rutter | The Kinetics of Rock Deformation by Pressure Solution [and Discussion[END_REF][START_REF] Gratier | Transition between seismic and aseismic deformation in the upper crust[END_REF]. Pressure solution theory is based on the presence at grain contacts of a thin aqueous film to provide the dissolved aqueous species necessary for crystal growth and transmit the normal stress across the grain contact [START_REF] Rutter | The Kinetics of Rock Deformation by Pressure Solution [and Discussion[END_REF][START_REF] Rutter | Pressure solution in nature, theory and experiment[END_REF]. The thinner estimates of aqueous film thickness ranges from 1 to 100 nm [START_REF] Renard | Pressure solution in sandstones: influence of clays and dependence on temperature and stress[END_REF][START_REF] Dysthe | Fluid in mineral interfaces -Molecular simulations of structure and diffusion[END_REF][START_REF] Desarnaud | The Pressure induced by salt crystallization in confinement[END_REF].

The existing models of pressure solution consider three steps: dissolution at grainto-grain contacts (source), diffusive transport in the thin aqueous film and precipitation in the pore space (sink) [START_REF] Rutter | The Kinetics of Rock Deformation by Pressure Solution [and Discussion[END_REF][START_REF] De Meer | On Mechanisms and Kinetics of Creep by Intergranular Pressure Solution[END_REF][START_REF] Niemeijer | Compaction creep of quartz sand at 400-600 • C: experimental evidence for dissolution-controlled pressure solution[END_REF][START_REF] Gundersen | Coupling between pressure solution creep and diffusive mass transport in porous rocks[END_REF]. They predict strain rate as a function of stress, dissolution and precipitation rates, diffusivity and geometrical parameters. Three main diffusive pathway geometries are used [START_REF] Den Brok | Effect of microcracking on pressure-solution strain rate: The Gratz grain-boundary model[END_REF]: a thin fluid film model with equally spaced grain boundaries, an island model with fluid-filled channels surrounding grains in contact [START_REF] Raj | Creep in polycrystalline aggregates by matter transport through a liquid phase[END_REF], and a micro-cracking model with complex fluid pathways continuously transformed by stress-corrosion at grain boundaries [START_REF] Gratz | Solution-transfer compaction of quartzites: Progress toward a rate law[END_REF]. These geometries play a first-order role on the estimated strain rates [START_REF] Xia | Structural and rheological evolution of the Laramide subduction channel in southern California[END_REF]. The pressure solution models proved their efficiency in reproducing observations in natural samples (sedimentary rocks) and in experiments [START_REF] Gratier | A pressure solution creep law for quartz from indentation experiments[END_REF][START_REF] Niemeijer | Compaction creep of quartz sand at 400-600 • C: experimental evidence for dissolution-controlled pressure solution[END_REF]. However, they are not directly applicable to metamorphic or hydrothermal systems since they consider a single phase. Observations in natural samples (A. [START_REF] Putnis | Fluid-Induced Processes: Metasomatism and Metamorphism[END_REF] and in experiments (C. [START_REF] Putnis | Direct observations of pseudomorphism: Compositional and textural evolution at a fluid-solid interface[END_REF][START_REF] Putnis | Mineral Replacement Reactions[END_REF][START_REF] Raufaste | The mechanism of porosity formation during solvent-mediated phase transformations[END_REF] indeed indicate that reactions occur through a replacement process in which primary min-eral dissolution and secondary mineral precipitation are coupled at the mineral interface.

A critical point to model dissolution-precipitation process in metamorphic rocks is thus to allow for precipitation at the stressed grain boundaries.

Structural features such as strain shadow, grain shape preferred orientation (GSPO) and crenulation cleavage provide evidence for the occurrence of dissolution-precipitation creep in metamorphic rocks. They are used to understand how deformation is accommodated from the scale of a single fault or shear zone [START_REF] Renard | Kinetics of cracksealing, intergranular pressure solution, and compaction around active faults[END_REF][START_REF] Gratier | Aseismic sliding of active faults by pressure solution creep: Evidence from the San Andreas Fault Observatory at Depth[END_REF][START_REF] Mcaleer | Reaction softening by dissolutionprecipitation creep in a retrograde greenschist facies ductile shear zone, New Hampshire, USA[END_REF] to the scale of a subduction zone [START_REF] Shimizu | Kinetics of pressure solution creep in quartz: theoretical considerations[END_REF][START_REF] Wassmann | Rheology of the plate interface -Dissolution precipitation creep in high pressure metamorphic rocks[END_REF][START_REF] Platt | Rheology and stress in subduction zones around the aseismic/seismic transition[END_REF]. The acquisition of GSPO during dissolutionprecipitation requires mass redistribution from highly to low stressed surfaces. It was recognized in a variety of rocks including marble [START_REF] Bestmann | Development of single-crystal σ-shaped quartz porphyroclasts by dissolutionprecipitation creep in a calcite marble shear zone[END_REF], amphibolite [START_REF] Imon | Development of shapeand lattice-preferred orientations of amphibole grains during initial cataclastic deformation and subsequent deformation by dissolutionprecipitation creep in amphibolites from the Ryoke metamorphic belt, SW Japan[END_REF][START_REF] Stokes | Deformation of amphibolites via dissolution-precipitation creep in the middle and lower crust[END_REF], serpentinite [START_REF] Wassmann | Dissolution precipitation creep versus crystalline plasticity in high-pressure metamorphic serpentinites[END_REF][START_REF] Bial | The microstructural record of porphyroclasts and matrix of partly serpentinized peridotite mylonites from brittle and crystal-plastic deformation to dissolutionprecipitation creep[END_REF] and micaschist (Wassmann & Stöckhert, 2013a). However, other processes can lead to the acquisition of GSPO [START_REF] Kretz | Interpretation of the shape of mineral grains in metamorphic rocks[END_REF][START_REF] Means | Experiments on preferred orientation of platy minerals[END_REF][START_REF] Willis | A kinematic model of preferred orientation[END_REF][START_REF] Brady | Intergranular diffusion in metamorphic rocks[END_REF][START_REF] Awad | Forsteritic olivine: effect of crystallographic direction on dissolution kinetics[END_REF]. The possible occurrence of all these mechanisms during metamorphism complicates the interpretation of GSPO in rocks and calls for a better understanding of the conditions required to form a GSPO during dissolution-precipitation.

Dissolution and precipitation under stress could play a fundamental role on the thermodynamics of reactions in metamorphic environments. [START_REF] Wheeler | Dramatic effects of stress on metamorphic reactions[END_REF][START_REF] Wheeler | The effects of stress on reactions in the Earth: Sometimes rather mean, usually normal, always important[END_REF] proposes that the reaction pathway followed during reaction by pressure solution does not allow thermodynamic equilibrium in stressed environments to be defined. However, the stability fields of the various mineralogical assemblages can still be determined and the reaction lines can be shifted by several hundreds of MPa under stress. These results are debated in the metamorphic community [START_REF] Fletcher | Dramatic effects of stress on metamorphic reactions[END_REF][START_REF] Hobbs | Does non-hydrostatic stress influence the equilibrium of metamorphic reactions?[END_REF][START_REF] Moulas | Relation between mean stress, thermodynamic, and lithostatic pressure[END_REF] and developing a quantitative model able to reproduce the reaction pathway proposed in [START_REF] Wheeler | Dramatic effects of stress on metamorphic reactions[END_REF] would help determining the effect of stress on thermodynamic equilibrium.

Here we develop a model for reaction under stress by considering both reaction kinetics and dissolved aqueous species transport. The model allows creep laws for pressure solution and dissolution-precipitation to be determined. These laws give insights into the importance of dissolution-precipitation creep for deformation in the Earth's crust.

We also use the model to determine the conditions favouring the acquisition of GSPO during metamorphism and to constrain the pressure controlling the metamorphic reactions observed in natural systems. The model predictions are compared to observations in metasediments from the Nufenen Pass (Swiss Alps).

2 Modelling replacement under stress

Theoretical framework

Fluid promotes metamorphic reactions such as the quartz-coesite transition, the albite = jadeite + quartz geobarometer or the eclogitisation of the lower crust. We develop in the following a model to estimate the first-order effect of stress on such reactions. Parameters used in the model are listed in Table 1.

Let us consider a stressed porous rock matrix reacting at constant temperature under the simplified geometry of Figure 1. The elementary unit of the model is a cubic volume initially composed of the reactant mineralogical assemblage engulfed in a porous network. We study the deformation of this elementary unit during mass redistribution and volume change associated with a replacement reaction under stress. The geometry used here allows for calculations at small porosity and simple upscaling with a face-centered cuboid packing. It differs from the face-centered cubic packing of equally sized spherical grains generally used in pressure solution studies [START_REF] Yasuhara | A mechanistic model for compaction of granular aggregates moderated by pressure solution[END_REF]. This latter cannot be used below a critical porosity of 3.49% corresponding to the porosity of the densely

L 3 L 2 L 1 σ 1 σ 2 σ 3 σ 1 > σ 2 > σ 3 A B L tot1 L tot2 L tot3
x 1

x 2

x 3

Figure 1. Geometry of the reaction under stress model. A: considered elementary unit. The reacting grain (light grey) is replaced on its stressed surfaces by reaction products (dark grey).

The reactants and the products are separated on each face by an aqueous film in which dissolved aqueous species are transported by diffusion (the aqueous films are displayed in blue). Exchanges of dissolved aqueous species between the grain surfaces occur through the porous network occupying the edges and the corners of the cuboid. B: grain packing (left) and associated pore network (right) for an assemblage of 8 cuboids.

packed spheres when the faces of adjacent spheres touch. The fluid in this rock occupies both the pores and the thin films at grain boundaries perpendicular to three stress components (σ α such that α = 1, 2 or 3) respecting the following inequality:

σ 1 > σ 2 > σ 3 .
The reference frame is defined by the stress components with coordinate axis vectors, x α , parallel to each stress component having the same subscript. The fluid pressure is equal to a constant in the pores (P f,0 ). In the aqueous film, thermodynamic calculations are performed by using the fluid pressure [START_REF] Dahlen | Metamorphism of nonhydrostatically stressed rocks[END_REF] which is equal to the normal stress σ α [START_REF] Israelachvili | Intermolecular and Surface Forces[END_REF]. This implies that we use a thin film grain boundary model. For simplicity, the fluid pressure in each film is denoted in the following with the stress it is equal to (σ α ).

We decompose the replacement reaction in the dissolution and then the precipitation of two different mineralogical assemblages: the initially rock-forming mineral(s)

(denoted with subscript r in the following) and the reaction product(s) (denoted with subscript p). At P f,0 , the initially rock-forming minerals are more soluble than the replacing minerals and they will thus dissolve whereas the replacing mineral will precipitate. The aqueous film is considered to be located in between the reacting and produced mineralogical assemblages (Figure 1). The chemical compositions of the reacting and produced mineral assemblages are identical in the model. (De-)volatilization reactions are thus not taken into account here. The involved dissolved aqueous species are identical for the reacting and precipitating phases. We also assume equal transport properties for all the mineral assemblage components in the fluid. This is in line at first-order with the calculations of [START_REF] Oelkers | Calculation of the thermodynamic and transport properties of aqueous species at high pressures and temperatures: Aqueous tracer diffusion coefficients of ions to 1000C and 5 kb[END_REF] showing diffusivity of the same order of magnitude for ions in bulk fluid. However, this does not allow to distinguish single minerals in each mineralogical assemblage. Such a level of complexity is considered in [START_REF] Wheeler | Importance of pressure solution and coble creep in the deformation of polymineralic rocks[END_REF] and [START_REF] Ford | Modelling interface diffusion creep in two-phase materials[END_REF] with more complex geometries involving several grains but it cannot be modelled quantitatively due to the high uncertainty on diffusion coefficients. In summary, the aim here is to provide first-order quantitative estimates for replacement under stress by focusing on mass redistribution in one single cuboid of our model. This requires assuming constant chemical composition in the cuboid, the same chemical composition for reactants and products, and the same diffusion coefficients for all the dissolved aqueous species.

For each mineralogical assemblage, the rate of dissolution/precipitation (in mol.s

where k {r,p} is a kinetic constant (with unit of mol.m -2 .s -1 ), S {r,p},α is the geometrical reactive surface area of the considered mineralogical assemblage in the aqueous film normal to σ α , R is the gas constant and T is the temperature. A {r,p},α is the chemical affinity of reaction of dissolution or precipitation under stress σ α . It is defined at the fluid/solid interface as the difference of the chemical potentials of the minerals under their solid state (µ s,{r,p},α ) and of the minerals dissolved in the fluid as aqueous species (µ f,α ): A {r,p},α = µ s,{r,p},α -µ f,α . A {r,p},α is first expressed as a function of the activity of the minerals dissolved as aqueous species in the fluid (a α ) and of their activity in the fluid at the equilibrium (a e,{r,p},α ). Then activities are replaced by the concentration of the dissolved minerals (c α ) by assuming that the activity coefficients do not depend on concentration. c e,{r,p},α is the concentration of the dissolving or precipitating phases at the equilibrium under stress σ α . The composition of the dissolved minerals formed by dissolution or consumed during precipitation are identical. There is thus no need to distinguish the activity (or concentration) of the dissolved minerals between reactants and products, and a α and c α are not mineral-specific.

Dissolved aqueous species migrate in the thin aqueous film as a result of chemical potential gradients between the aqueous film and the pores. The transport of dissolved aqueous species in the planar thin film normal to x α is modelled by diffusion. The release and consumption of aqueous species by dissolution or precipitation, respectively, are modelled with two source terms for dissolution and precipitation occurring at the minerals/fluid interfaces. These terms are obtained by dividing the reaction rate of equation 1 by the volume of the thin aqueous film. We assume here that the reactive surface area of the dissolving and precipitating phases in each fluid film are identical (S r,α = S p,α ) and that the aqueous film thickness (w) is constant so the thin film volume is S p,α w.

The equation describing the diffusion of dissolved aqueous species and their release/consumption by reaction is finally:

∂c α ∂t = D ∂ 2 c α ∂x 2 β + ∂ 2 c α ∂x 2 γ + k r w 1 - c α c e,r,α + k p w 1 - c α c e,p,α (2)
where D is the diffusion coefficient of dissolved aqueous species in the aqueous film, β = γ = α = 1, 2 or 3, and the subscript α refers to the normal stress σ α at which the concentrations are calculated.

To solve the differential equation 2, we use as boundary conditions the concentration of dissolved mineral in the pore (c 0 ). Diffusion in the pore proceeds rapidly compared to diffusion at grain boundaries [START_REF] Brady | Intergranular diffusion in metamorphic rocks[END_REF]. Therefore, we assume no spatial variations of c 0 . In the following, c 0 is either fixed by assuming a connection with an infinite source of dissolved aqueous species or calculated by assuming closed system conditions:

∂c 0 ∂t = - 4Dw V 0 α=1,2,3;β =α;γ =α,β L β ∂c α ∂x α bβ + L γ ∂c α ∂x α bγ + k r S 0,r V 0 1 - c 0 c e,r,0 + k p S 0,p V 0 1 - c 0 c e,p,0 (3)
where L is the length of the thin film, and the subscript 0 refers to the fluid pressure P f,0 at which the concentrations are calculated. In the geometry used here, each pore is shared between 8 elementary units (Figure 1) and each cuboid summit contains a eight part of a pore. V 0 is therefore the volume of one pore, and S 0,r and S 0,p are the reactive surface area in this pore of the reactants and the products, respectively. The two last terms on the right hand side correspond to precipitation/dissolution in the pore.

The first term is the Fickian flux of matter from the thin films to the pores at boundaries bα, bβ and bγ parallel to x α , x β and x γ , respectively.

Analytical solution in one dimension

It is instructive to study an analytical solution of equations ( 2) and ( 3). The geometry is simplified to recover a solution at steady state in one dimension. The aqueous films are connected to pores in x = 0 and x = L with x the distance along the aqueous film, so that c α (x = 0) = c α (x = L) = c 0 . The composition in the films only varies in one direction (infinitely long sheet-like films in the other direction). The solution of equation 2 at steady state in one dimension is:

c α (x) = c b,α + (c 0 -c b,α ) exp Da α x L + exp Da α 1 -x L 1 + exp (Da α ) (4)
where c b,α = kr+kp kr ce,r,α + kp ce,p,α is the concentration at steady state in the plane normal to x α when diffusion is fast compared to reaction (i.e. the amount of dissolving minerals is exactly balanced by the amount of precipitating minerals). Da α is the Damköhler number defined here as:

Da α = kr ce,r,α + kp ce,p,α L 2 Dw (5)
Da α compares the timescales of the processes involved in the system, that is reaction and diffusion. It depends on stress as the reaction rates are stress dependent. At high Da α , diffusion is slower than reaction and the system is thus controlled by diffusion. Conversely, the system is controlled by reaction at low Da α .

Equation 4 reveals that dissolution and precipitation reactions release and consume dissolved aqueous species in order to reach a concentration of c b,α in the aqueous film at which input by dissolution and output by precipitation are exactly balanced. However, this concentration is never reached due to diffusion from or towards the pores in x = 0 and x = L. The direction of diffusive transport depends on the relative values of c 0 and c b,α .

The flux of concentration at boundary x = 0 can be derived from equation 4: It is apparent from equations ( 4) and ( 7) that the concentrations are controlled by the Damköhler number, that is by the slowest process between diffusion and surface-controlled kinetics in the pores and in the thin aqueous films. Table 2 gives the borderline cases for c α and c 0 . The variations of Da 0 compared to Da α are mainly due to differences between the reactive surface area of the pores and the grain boundary, respectively.

∂c α ∂x x=0 = (c 0 -c b,α ) Da α L 1 -exp (Da α ) 1 + exp (Da α ) (6)

Expression of strain rate

Strain rate ( ˙ α ) allows the various deformation mechanisms in the Earth to be compared. As it is usually done in pressure solution studies, we derive in the following strain rates for dissolution-precipitation creep when the system is controlled by diffusion ( ˙ αD ) and reaction ( ˙ αR ). We also introduce strain rate for a system controlled by volume change during reaction ( ˙ αV ), which is specific to replacement reactions. As in most pressure solution studies, strain rate is calculated here as positive when the rock compacts and it is calculated with the solution for mass transfer in one dimension [START_REF] Gratier | The Role of Pressure Solution Creep in the Ductility of the Earth's Upper Crust[END_REF].

It can be directly recovered from the reaction rates in the thin films:

˙ α = 2k r V r L 1 - c α c e,r,α + 2k p V p L 1 - c α c e,p,α (8)
where V r and V p are the molar volumes of the replaced and replacing solid phases, respectively. These volumes are involved here since a replacement reaction can induce deformation through volume change. The presence of L in the denominator is related to the distribution of volume change in the considered volume.

At high Da α , equation 4 leads to c α = c b,α and equation 8 can thus be re-written as:

˙ αV = -2 ∆ r V L k m (9)
where k m = c e,{r,p},α depends on stress [START_REF] Paterson | Nonhydrostatic thermodynamics and its geologic applications[END_REF] as: c e,{r,p},α = c e,{r,p},0 exp (σ α -P f,0 ) V {r,p} RT

where c e,{r,p},0 is the concentration of the dissolved or precipitating minerals at the equilibrium at P f,0 . We obtain with equations 9 and 10 for k r k p :

˙ αV = -2 k r ∆ r V L 1 - c e,p,0 c e,r,0 exp (σ α -P f,0 ) ∆ r V RT (11)
and for k p k r :

˙ αV = 2 k p ∆ r V L 1 - c e,r,0 c e,p,0 exp - (σ α -P f,0 ) ∆ r V RT (12)
We consider in the following that equilibrium is achieved at P f,0 (c e,r,α = c e,p,α when σ α = P f,0 ). By definition, the dissolving assemblage has a higher solubility than the precipitating one (c e,r,α > c e,p,α ). These two latter statements as well as equation

10 imply that V r < V p (∆ r V > 0) if σ α < P f,0 . Conversely, ∆ r V is negative when σ α > P f,0 . As a result, the product (σ α -P f,0 ) ∆ r V is
˙ α = 8D ∂c α ∂x x=0 wV r L 2 (13)
Equation 13 can then be combined with equation 6 to obtain:

˙ α = 8Da α DwV r L 3 (c 0 -c b,α ) 1 -exp (Da α ) 1 + exp (Da α ) (14)
At Da α > 5, 1-exp (Daα) 1+exp (Daα) ∼ -1 and equation 14 can be simplified to recover the strain rate for a diffusion-controlled system:

˙ αD = - 8Da α DwV r L 3 (c 0 -c b,α ) (15)
At low Da α (control by reaction), the strain rate can be retrieved from equation 8 by substituting c α with the expression of c 0 for the different borderline cases (see Table 2).

Pressure solution

Pressure solution models consider only one phase which dissolves at grain contacts.

This corresponds to a model proposed here if only considering the reaction of dissolution. There is therefore no replacement reaction and the terms of the previous section with p subscript can all be dropped. Moreover, the geometry is slightly modified since the aqueous film needs to be shared between two elementary units. The strain rate of equation 14 needs thus to be divided by 2.

The strain rate expressed in equation 8 becomes:

˙ α = 2k r V r L 1 - c α c e,r,α (16)
This latter equation can be re-written with the limits given in Table 2 at low Da α to express the strain rate for a system controlled by reaction:

˙ αR = 2k r V r L 1 - c 0 c e,r,α (17)
The consideration of a one phase system also leads to c b,α = c e,r,α which allows equation 14 to be re-written as:

˙ α = - 4k r V r L 1 - c 0 c e,r,α 1 Da α 1 -exp (Da α ) 1 + exp (Da α ) (18)
For a system controlled by diffusion (high Da α ), equation 18 becomes:

˙ αD = 4Da α c e,r,α DwV r L 3 1 - c 0 c e,r,α (19)
Equations 17 and 19 provide strain rates in the reaction-controlled and diffusioncontrolled cases, respectively. The only unknown is c 0 in these equations. We determine this unknown under two assumptions (small or large porosity) in the following.

Firstly, at small porosity (Da 0 Da α ; closed system), the concentration in the pores is c 0 = α ce,r,αDaα α Daα (Table 2). For non-hydrostatic stresses, grain surfaces normal to σ 1 will dissolve whereas phase precipitation will occur on the grain surfaces normal to σ 3 . The rate of this mass transfer at the grain scale can be quantified with equations 17 and 19.

Secondly, we can assume that pores are connected to a large sink/source of aqueous species buffering the system composition at an equilibrium value of c 0 = c e,r,0 . This is the assumption used in most pressure solution studies when open system conditions are considered. For the reaction-controlled case, combining equation 17 with equation 10 under the assumption c 0 = c e,r,0 leads to:

˙ αR = 2k r V r L 1 -exp -(σ α -P f,0 ) V r RT ≈ 2k r V 2 r (σ α -P f,0 ) RT L (20)
When the system is controlled by diffusion (high Da α ) under the assumption of a large aqueous species reservoir, the strain rate is recovered from equation 19 with c 0 = c e,r,0 :

˙ αD = - 4Da α c e,r,0 DwV r L 3 1 -exp (σ α -P f,0 ) V r RT ≈ 4Da α c e,r,0 DwV 2 r (σ α -P f,0 ) RT L 3 (21)

Deformation mechanism maps

Equation 18 can be used to generate a deformation mechanism map, allowing pressure solution to be compared with other deformation mechanisms, that is dislocation and diffusion creeps. The parameters needed to plot such a map are the diffusion coefficient, D, the reaction rate, k r , the aqueous film thickness, w, the grain size, L, and the concentration at the equilibrium in the aqueous film, c e,r,α . Experiments of pressure solution [START_REF] Van Noort | Effects of orientation on the diffusive properties of fluid-filled grain boundaries during pressure solution[END_REF][START_REF] Van Noort | Structure and properties of loaded silica contacts during pressure solution: impedance spectroscopy measurements under hydrothermal conditions[END_REF], molecular dynamics simulations [START_REF] Dysthe | Fluid in mineral interfaces -Molecular simulations of structure and diffusion[END_REF] and theoretical studies [START_REF] Revil | Pervasive pressure solution transfer in a quartz sand[END_REF] give estimates of D of ∼ 10 -10 m 2 .s -1 , that is one order of magnitude lower than diffusivity in bulk fluid [START_REF] Oelkers | Calculation of the thermodynamic and transport properties of aqueous species at high pressures and temperatures: Aqueous tracer diffusion coefficients of ions to 1000C and 5 kb[END_REF]. In practice, we use in the following a diffusion coefficient depending on temperature such as:

D = 1.2.10 -9 exp - 10.10 3 RT (22)
The dependency in temperature in this latter equation is derived from the diffusivity estimates of [START_REF] Oelkers | Calculation of the thermodynamic and transport properties of aqueous species at high pressures and temperatures: Aqueous tracer diffusion coefficients of ions to 1000C and 5 kb[END_REF]. We use an aqueous film thickness of w = 2.10 -9 m based on the study of [START_REF] Dysthe | Fluid in mineral interfaces -Molecular simulations of structure and diffusion[END_REF]. The grain size is fixed to L = 100 µm. c e,r,α is estimated with equation 10 including the aqueous species concentration at the equilibrium at P f 0 (c e,r,0 ) calculated with Perple X [START_REF] Connolly | Computation of phase equilibria by linear programming: A tool for geodynamic modeling and its application to subduction zone decarbonation[END_REF] (black dots). The black line corresponds to the best fit (kr = 2.9.10 -5 mol.m -2 .s -1 ). The strain rate computed with the reaction rate determined in a bulk fluid (kr = k rbulk = 0.43 mol.m -2 .s -1 ) [START_REF] Tester | Correlating quartz dissolution kinetics in pure water from 25 to 625 • C[END_REF][START_REF] Worley | Quartz dissolution kinetics from 100-200 • C as a function of pH and ionic strength[END_REF] is also displayed with a red line. kr deduced from the fit of the experimental data is four of magnitude lower than k rbulk . B: Strain rate as a function of the ratio between actual and initial porosities (φ/φ0). The data from [START_REF] Niemeijer | Compaction creep of quartz sand at 400-600 • C: experimental evidence for dissolution-controlled pressure solution[END_REF] are displayed with dots and one colour for each experiment: blue, red, black, yellow, cyan, magenta and green for their experiments CPf3, CPf4, CPf5, CPf6, CPf7, CPf8 and CPf9, respectively. These experiments are performed at different temperature and effective pressure (see [START_REF] Niemeijer | Compaction creep of quartz sand at 400-600 • C: experimental evidence for dissolution-controlled pressure solution[END_REF] for details). The lines correspond to the best fit of the data for which kr = k rbulk /15. Dashed lines correspond to the strain rates calculated with kr = k rbulk . The colors of the lines are identical to the color used for the experimental results (dots) when the stress and the temperature used for the computation are identical to those used in the experiments.

The distribution of the dominant deformation mechanisms is similar for the two reaction rates used for calculating the deformation mechanism maps (Figure 3). Diffusion creep is not a dominant deformation mechanism in the investigated stress-temperature domain. Dislocation creeps dominates at stress above ∼ 10 8 P a and temperature above The pressure solution model proposed above assumes no precipitation in the aqueous film. This does not allow dissolution-precipitation creep in small porosity rocks to be modelled as it is commonly observed in hydrothermal and metamorphic systems. We use in the following the equations of section 2.2 to extend our understanding of polymineralic assemblages by calculating first-order estimates of strain rate for dissolution-precipitation creep. Such estimates allow to compute deformation mechanism maps including dissolutionprecipitation creep (Figure 4). We consider as an example the metamorphic reaction of muscovite + quartz <=> andalusite + sanidine + H 2 O. We assume isotropic stress and small porosity (Da 0 → 0) so that equation 9 can be used to calculate strain rate at both high and low Da α . This strain rate corresponds to the strain rate associated with volume change during reaction ( ˙ αV ). The muscovite dehydration involves water exchange between the solid and the fluid phases which is not considered in our model. However, it is one of the few metamorphic reactions for which reaction rates estimates are available as a function of pressure and temperature [START_REF] Schramke | The reaction muscovite + quartz = andalusite + K-feldspar + water; Part 1, Growth kinetics and mechanism[END_REF]. The use of this reaction as a general kinetic model for metamorphic reactions is justified by the fact that the kinetic data for many metamorphic reactions (including or not volatile species) can be rather accurately fitted with a single rate law [START_REF] Wood | Rates of hydrothermal reactions[END_REF]. Calculations are performed by considering that andalusite + sanidine + H 2 O react to produce muscovite + quartz at temperatures lower than the equilibrium whereas the opposite reaction occurs at higher temperature. Equilibrium is thus never achieved in Figure 4 except on the reaction line and reaction rates progressively increase away from this line in agreement with equation 1. This diagram is therefore only valid if disequilibrium conditions are achieved through changes in composition (e.g. fluid input) or in pressure and/or temperature. We fitted the data of [START_REF] Schramke | The reaction muscovite + quartz = andalusite + K-feldspar + water; Part 1, Growth kinetics and mechanism[END_REF] with a kinetic law as follows:

∼
k = k rbulk 1 -exp -A t,α,β RT n
where k rbulk = λ exp -Ea RT and A t,α,β is the distance from the equilibrium defined as the difference between the chemical potential of the dissolving minerals (µ s,r,α ) and the chemical potential of precipitating minerals (µ s,p,β ): A t,α,β = µ s,r,α -µ s,p,β . Two different fits are produced for hydration and dehydration. k r and k p are assumed equal. There are no experimental data for dissolution-precipitation creep which can be used to determine the relationship between k r and k rbulk as we did for pressure solution. We therefore used the same approach as for pressure solution and plotted deformation mechanism map with k r = k rbulk and k r = k rbulk /(10 4 ) (Figure 4).

The concentrations at the equilibrium (c e,r,α and c e,p,α ) are calculated with Perple X [START_REF] Connolly | Computation of phase equilibria by linear programming: A tool for geodynamic modeling and its application to subduction zone decarbonation[END_REF] We first validate the model by reproducing the analytical results presented above.

At low Da α (∼ 0.12), the model predicts a constant concentration in the whole system (aqueous films and pores; c α = c 0 ) which can be approximately calculated with the equation for concentration in the pores derived in one dimension (Table 2; Figure S1 A). The strain rates derived analytically in one dimension at low Da α are also reproduced with a good precision (Equation 8; Figure S1 C). At high Da α (∼ 380), the calculated concentration and strain rates are also approximately equal to the analytical solutions (Table 2 and Equation 9; Figure S1 B and D).

Grain shape preferred orientation

To better understand the acquisition of GSPO during replacement under stress, we numerically investigate the influence of Da α on GSPO by running several simulations with different reaction rates. GSPO is quantified with the aspect ratio calculated after complete replacement and defined as the length of the grain side normal to σ 3 over the length of the grain side normal to σ 1 ( Ltot3 Ltot1 ). At Da α > 100, the aspect ratio is close to 1 and there is no change in shape during reaction (Figure 6). A strong change in GSPO occurs at Da α < 10 with an aspect ratio reaching 3.3 at Da α = 0.282 (Figure 6).

Pressure controlling the replacement process

The reaction rate depends on the distance from the equilibrium (A t,α,β ; equation 1). As the chemical potential depends on pressure, the reaction rates in two differently oriented grain surfaces under different normal stresses are different. Reaction in one of The reaction line where equilibrium is reached is shown in red and separates two fields in which we consider either hydration (left) or dehydration (right). The contours correspond to log10 ( ˙ ).

-1 3 3 1 -3 1 - -13 3 1 - 2 1 -2 1 - -1 2 1 1 -1 1 - -1 1 -1 1 - -10 - 1
The domains where one of the three studied deformation mechanism (dissolution-precipitation, dislocation and diffusion creeps) dominates are displayed. The yellow and orange areas correspond to the regions where prograde metamorphism and hydrothermal conditions are expected to occur, respectively.

-15-manuscript submitted to Geochemistry, Geophysics, Geosystems The aspect ratio is defined as the length of the grain side normal to σ3 over the length of the grain side normal to σ1 (

L tot,3
L tot,1 ). Each dot corresponds to one simulation. The differential stress (σ1 -σ3) is of 100 MPa (

L tot,3
L tot,1 would be equal to 1 if the grain would grow under isotropic stress). and σ3 are displayed with black, red and blue dots, respectively. We modify the distance from the equilibrium by increasing the difference between σ3 and the equilibrium pressure (Peq). The simulations are performed at high Da with kr = kp = 1.10 -6 mol.m -2 .s -1 and fixed differential stresses σ2 -σ3 and σ1 -σ3 of 50 and 100 MPa, respectively.

the thin aqueous film can thus dominate the replacement process, leading to a control of the reaction by the stress normal to this film.

To quantitatively investigate this process, we run simulations to determine the contribution of each face of the cuboid to the overall reaction progress as a function of the distance from the equilibrium (Figure 7). The distance from the equilibrium increases with pressure in the simulations (A t,1,1 > A t,2,2 > A t,3,3 ). This is equivalent to con-sider that the reaction involves a negative change in volume (∆ r V = V p -V r < 0; left side of the reaction line for a reaction with a positive Clapeyron slope or right side for a reaction with a negative Clapeyron slope). We run several simulations with an increasing difference between σ 3 and the equilibrium pressure (P eq ) while the differential stresses, σ 2 -σ 3 and σ 1 -σ 3 , are fixed to 50 and 100 MPa, respectively. At low Da α , diffusion proceeds more rapidly than reaction and mass is redistributed between the aqueous films with precipitation which can occur in a different aqueous film than dissolution. The simulations are thus run at Da α > 1 to be able to isolate the contribution of each face to reaction. Most of grain conversion occurs on the face normal to σ 1 near the equilibrium (Figure 7). The pressure controlling reaction is thus σ 1 . For reactions with ∆ r V > 0, the pressure controlling reaction is σ 3 . Far from the equilibrium, reaction is equally distributed among the faces (Figure 7) and the thermodynamic pressure is thus the mean stress.

4 Grain shape preferred orientation in natural samples from the Nufenen pass area

Here we present an application of the model developed above to evaluate GSPO data from a regional metamorphic terrain. The Mesozoic sedimentary cover of the Gotthard massif is well exposed in the Nufenen pass area (Swiss Central Alps). We collected along the road of the Griessee intensively foliated graphite-rich schists. We focused on two different types: one has garnet porphyroblasts (0.6 to 8 mm in diameter), the other one has plagioclase porphyroblasts (1 to 10 mm diameter). The peak metamorphic conditions were estimated by [START_REF] Klaper | Alpine metamorphism of pelitic schists in the Nufenen Pass area, Lepontine Alps[END_REF] and [START_REF] Kamber | Regional metamorphism and uplift along the southern margin of the Gotthard massif; results from the Nufenenpass area[END_REF] at 450-500 • C and 0.5 GPa. In addition to plagioclase and garnet, the rock contains quartz, carbonates, muscovite, biotite, clinozoisite, ilmenite and carbonaceous matter. Garnet and plagioclase are often slightly oval, and visual investigation of the direction of the long axis of the porphyroblasts show that they are not aligned in the pronounced schistosity visible in the field (Figure 8). The schistosity is defined by mica. The stratigraphy, deformation and metamorphism were studied in detail for this area by Liskay-Nagy (1965); [START_REF] Leu | Geologie der Sedimentzüge zwischen Griessee und Passa del Corno (Nufenengebiet, Wallis)[END_REF]; [START_REF] Klaper | Alpine metamorphism of pelitic schists in the Nufenen Pass area, Lepontine Alps[END_REF]; [START_REF] Kamber | Regional metamorphism and uplift along the southern margin of the Gotthard massif; results from the Nufenenpass area[END_REF]. A detailed review of the different nomenclature and tectonic interpretations made by these authors is beyond interest of this study. Nevertheless, we agree with the observations of these authors, that in the profile studied here the peak metamorphic minerals plagioclase and garnet grew after the second folding phase (D2, following the nomenclature of Liskay-Nagy (1965)), which was responsible for the major schistosity in the area. The final deformation phase (D3) did in places result in a rotation of the blasts, visible in a rotation of the schistosity. This was not observed in the samples presented here.

The 3D orientation of the porphyroblasts is determined in seven plagioclase-bearing and six garnet-bearing samples by X-ray microtomography (µcT) with a Skyscan 1172 device at the University of Lausanne. Different settings are used for the two types of schists since garnet and plagioclase have higher and comparable X-ray attenuation than the rock matrix, respectively. Garnet-bearing samples were analysed at 130 kV and 60 µA with a resolution comprised between 20 and 25 µm and a brass filter. Plagioclase-bearing samples are analysed at 70 kV and 114 µA with a resolution comprised between 11 and 16 µm and an aluminium filter. 3D images are generated from transmission images with the NRecon software. They are then segmented with blob3D [START_REF] Ketcham | Computational methods for quantitative analysis of three-dimensional features in geological specimens[END_REF]) by considering three phases: the porphyroblast, ilmenite and the matrix. We fit an ellipsoid to each segmented grain to determine the grain aspect ratio (length of the longest axis over length of the shortest axis), the orientations of the main elongation axis and of the main shortening plane. Spherical orientation data are projected in 2D and statistically analysed using Bingham axial distribution with the Stereonet software [START_REF] Cardozo | Spherical projections with OSXStereonet[END_REF]. Porphyroblasts at the edge of the core and/or coalescing with other grains were excluded from the analysis. Ilmenite grains are elongated in the foliation; their main shortening plane thus gives foliation orientation. The results of statistical analysis are given in Figure 9 and Table 3.

We determine the orientation and the aspect ratio of 13 to 42 plagioclases in 7 different cores and 40 to 96 garnets in 6 different cores (Table 3). Ilmenite crystals outline the orientation of the foliation. This orientation is similar to the schistosity measured in the field for garnet-bearing samples whereas plagioclase-bearing samples collected near a fold hinge display different schistosity and foliation orientations (Table 3 and Figure 9). Garnets have an ellipsoidal shape whereas plagioclases are only elongated in one direction (Figure 10). The mean aspect ratios are 2.3 ± 0.5 and 1.8 ± 0.5 for plagioclase and garnet, respectively (Figure S2). Their elongation axis is oblique to the sedimentary layering (S0) (Figure 9). Plagioclases appear as elongated in the foliation even though the dispersion of the data does not provide a precise estimate of their orientation. Garnets are flattened along a plane oblique to the schistosity/foliation at an angle between 15 and 39 • (Table 3 and Figure 9), confirming the field observations. Garnets contain numerous inclusions of quartz, epidote and ilmenite. These inclusions do not appear to be preferentially elongated (Figure S2). In each core except plagioclase-bearing sample LB 2018 4, the orientations of the elongation or of the flattening of the porphyroblasts are similar, indicating a GSPO at the centimetre scale (Table 3). Plagioclase orientation data are more scattered than garnet orientation data at this scale but also at the outcrop scale (Figure 9). Garnet-bearing cores have indeed similar GSPO with a main shortening axis ranging from 158/88 to 169/73 (Azimuth/Dip). The orientation of the main shortening axis of plagioclase varies between 5/62 and 325/80, indicating a variation of the GSPO at the tens of centimetre to the metre scale. These results are in agreement with the observations made in the field, showing that the foliation and the main orientation of plagioclase progressively rotate over a distance of several tens of centimetres (Figure 8 C). There is no relationship between the main orientation of plagioclases and their aspect ratio.

Discussion

The model proposed here implements two processes which were not considered in the previous models of pressure solution: the possibility for mineral precipitation in the thin aqueous films, that is at grain boundary, and a reaction rate depending on the distance from the equilibrium. These two features lead to new laws for pressure solution and dissolution-precipitation creeps. Moreover, they allow to predict the acquisition of GSPO and to discuss the pressure controlling metamorphic reactions. We discuss in the following the changes brought by these new laws for modelling reaction-driven deformation in the Earth's crust.

Comparison with previous creep laws for pressure solution

Equation 18 is a general creep law for pressure solution. It is obtained by solving for mass transfer in a thin aqueous film in one dimension, like in previous studies on pressure solution [START_REF] Rutter | The Kinetics of Rock Deformation by Pressure Solution [and Discussion[END_REF][START_REF] Niemeijer | Compaction creep of quartz sand at 400-600 • C: experimental evidence for dissolution-controlled pressure solution[END_REF][START_REF] Zhang | Compaction creep of wet granular calcite by pressure solution at 28 • C to 150[END_REF]. The grain size (L) does not evolve during reaction with this approach. Equation 18 is thus only valid if incremental deformation is assumed. Estimating the error associated with this assumption requires to solve for grain shape evolution during dissolution-precipitation. This is done with the numerical model where the strain rate is shown to decrease with time (Fig- ure S1). This decrease is limited to less than one order of magnitude for a reaction progress of less than 0.9. The incremental deformation assumption thus provides valid first-order strain rate estimates for a large range of reaction progresses.

Equation 18 is valid for all Da α values and thus allows strain rate to be calculated in both the reaction-and the diffusion-controlled cases. At low Da α (reaction-controlled regime), equation 18 can be approximated with equation 20, equivalent, at low stress ((σ α -P f,0 ) V r RT ), to the strain rate determined in reaction-controlled conditions in [START_REF] Raj | Creep in polycrystalline aggregates by matter transport through a liquid phase[END_REF]; [START_REF] Niemeijer | Compaction creep of quartz sand at 400-600 • C: experimental evidence for dissolution-controlled pressure solution[END_REF]; X. [START_REF] Zhang | Compaction creep of wet granular calcite by pressure solution at 28 • C to 150[END_REF]. However, at large stress, equation 20 tends to-wards krVr L . This is consistent with the transition state theory predicting reaction rate independent from the force driving reaction (A {r,p},α ) when A {r,p},α 0 [START_REF] Lasaga | Kinetic Theory in the Earth Sciences[END_REF].

The expression generally used in the literature for strain rate in the reaction-controlled case [START_REF] Raj | Creep in polycrystalline aggregates by matter transport through a liquid phase[END_REF][START_REF] Niemeijer | Compaction creep of quartz sand at 400-600 • C: experimental evidence for dissolution-controlled pressure solution[END_REF][START_REF] Zhang | Compaction creep of wet granular calcite by pressure solution at 28 • C to 150[END_REF] does not follow the transition state theory since it depends on effective stress as exp (σα-P f 0 )Vr RT -1 and thus tends towards infinity at large stress. For the pressure solution of quartz under a stress of 100 MPa, the ratio between the strain rates calculated with the expression found in the literature and with equation 20 is 1.6. This ratio exponentially increases with stress since it is 10 at 500 MPa. It could therefore lead to strong change in theoretical estimates especially when stress concentration occurs as it is the case in compacting rock aggregates [START_REF] Niemeijer | Compaction creep of quartz sand at 400-600 • C: experimental evidence for dissolution-controlled pressure solution[END_REF].

At high Da α (diffusion-controlled regime), equation 18 can be approximated with equation 21. This latter equation can be compared to the one generally obtained in pressure solution studies [START_REF] Rutter | The Kinetics of Rock Deformation by Pressure Solution [and Discussion[END_REF] which can be expressed for the geometry and the formalism used here as (see Appendix C for a detailed derivation):

˙ αS = - 24c e,r,0 DwV r L 3 1 -exp (σ α -P f,0 ) V r RT ≈ 24c e,r,0 DwV 2 r (σ α -P f,0 ) RT L 3 ( 23
)
This equation is the same as the one generally used for strain rate in the diffusionlimited case with the exception of the geometric factor 24 associated with the cubic geometry used here [START_REF] Rutter | The Kinetics of Rock Deformation by Pressure Solution [and Discussion[END_REF][START_REF] Gratier | The Role of Pressure Solution Creep in the Ductility of the Earth's Upper Crust[END_REF]. Its derivation requires two main assumptions (Appendix C): a constant reaction rate and thermodynamic equilibrium in the thin aqueous film [START_REF] Rutter | The Kinetics of Rock Deformation by Pressure Solution [and Discussion[END_REF]. These assumptions may seem contradictory since the driving force for reaction is disequilibrium in the thin aqueous film, and the reaction rate should thus be zero and not a constant at the equilibrium. The assumption of a constant reaction rate leads to a parabolic equation for aqueous species concentration (Appendix C) whereas we derive in equation 4 an hyperbolic cosine function by using a reaction rate depending on the distance from the equilibrium. The assumption of thermodynamic equilibrium in the aqueous film is responsible for the main difference between ˙ αD and ˙ αS , that is the absence of the Da α factor in equation 23 compared to equation 21 (Appendix C). The strain rate is underestimated by one order of magnitude at Da α = 100 with equation 23 compared to equation 21 (Figure 11). The dependency of ˙ αD on Da α indicates that the strain rate depends on reaction rate (k r ) even if the system is controlled by diffusion. An increase in reaction rate indeed leads to a higher dissolved aqueous species concentration in the aqueous film. As the gradient in concentration between the aqueous film and the pores drives diffusion, the transport rate increases leading to a higher strain rate.

The new laws provided here have thus not only different mathematical forms than the previous creep laws but they also lead to orders of magnitude differences in the estimated strain rates. This has to be considered when discussing the importance of pressure solution for deformation in the Earth's crust. The black plain line corresponds to the general expressions for strain rate (equation 18). The black dash-dot and dashed lines are the borderline cases when considering non-constant rates for low Daα (equation 20) and high Daα (equation 21), respectively. The red plain line is the strain rate generally used for diffusion-controlled system. It is derived with constant reaction rate and a condition for thermodynamic equilibrium in the aqueous film (equation 23). The calculations are performed for a typical grain size of L = 100 µm, σα = 100 MPa, the molar volume of quartz Vr = 22.10 -6 m 3 .mol -1 , w = 2.10 -9 m, D = 10 -10 m 2 .s -1 , ce,r,0 = 1mol.m -3 and T = 573K.

Dominant deformation mechanism in the Earth's crust

Strain rates are calculated for different Daα by varying kr from 10 -13 to 10 -5 mol.m -2 .s -1 based on the range given in the works of [START_REF] Tester | Correlating quartz dissolution kinetics in pure water from 25 to 625 • C[END_REF] and [START_REF] Worley | Quartz dissolution kinetics from 100-200 • C as a function of pH and ionic strength[END_REF]. jer et al., 2002;[START_REF] Dysthe | Fluid in mineral interfaces -Molecular simulations of structure and diffusion[END_REF][START_REF] Gratier | A pressure solution creep law for quartz from indentation experiments[END_REF][START_REF] Van Noort | Effects of orientation on the diffusive properties of fluid-filled grain boundaries during pressure solution[END_REF][START_REF] Van Noort | Structure and properties of loaded silica contacts during pressure solution: impedance spectroscopy measurements under hydrothermal conditions[END_REF].

However, the reaction rate in the aqueous film, k r , is poorly constrained. Previous studies generally assume that the reaction rate determined on unconfined surfaces (bulk fluid) can be used. However, the presence of catalytic dissolved aqueous species [START_REF] Rimstidt | Rate equations for sodium catalyzed quartz dissolution[END_REF] or the decrease in water activity [START_REF] Schott | The Link Between Mineral Dissolution/Precipitation Kinetics and Solution Chemistry[END_REF][START_REF] Lamadrid | Effect of water activity on rates of serpentinization of olivine[END_REF][START_REF] Rubie | The catalysis of mineral reactions by water and restrictions on the presence of aqueous fluid during metamorphism[END_REF] can modify the reaction rate in the aqueous film by several orders of magnitude compared to rates measured in bulk fluid. We fit here strain rate data from [START_REF] Niemeijer | Compaction creep of quartz sand at 400-600 • C: experimental evidence for dissolution-controlled pressure solution[END_REF] and [START_REF] Gratier | A pressure solution creep law for quartz from indentation experiments[END_REF] with equation 18 to estimate k r . The determined values are one to four orders of magnitude lower than the reaction rates measured in the bulk fluid. We therefore plot two deformation mechanism maps for pressure solution (Fig- [START_REF] Gratier | The Role of Pressure Solution Creep in the Ductility of the Earth's Upper Crust[END_REF]. The computed maps indicate orders of magnitude lower strain rates for pressure solution and dissolution-precipitation with k r = k rbulk / 10 4 than with k r = k rbulk (Figures 3 C and4 C). The distribution of the dominant deformation mechanisms is similar for the two reaction rates. However, it is different from the deformation mechanism map proposed by [START_REF] Rutter | The Kinetics of Rock Deformation by Pressure Solution [and Discussion[END_REF] in which diffusion creep dominates at temperatures higher than ∼ 600 • C. The use of equation 18 instead of equation 23 plays a key role in this switch from a map with three possible deformation mechanism to a map where pressure solution dominates. Indeed, recovering with equation 18 a deformation mechanism map similar to [START_REF] Rutter | The Kinetics of Rock Deformation by Pressure Solution [and Discussion[END_REF] requires a reaction rate eleven orders of magnitude smaller than the one measured in a bulk fluid (k r = 10 -16 mol.m -2 .s -1 for quartz at 350 • C). Figures 3 and4 thus provide robust distributions of the deformation mechanisms even though the absolute strain rate values will probably need to be re-evaluated when disposing of better constraints on reaction rate. These maps indicate that pressure solution and dissolution-precipitation creeps are the main deformation mechanisms in the Earth's crust.

One has to keep in mind the limitations of the model when using the deformation mechanism maps. First, we use a simple geometry with all the minerals dissolving or precipitating in the same cuboid, and we consider that dissolved aqueous species are transported at similar rates. We can thus not reproduce the complex distribution of minerals in natural samples, and transport distances can be underestimated by several times the grain size. These are second-order issues regarding uncertainties in quantitative modelling which are mainly due to the variation in reaction rate estimates by several orders of magnitude. Moreover, other models are available to understand mineral distribution or to compare the deformation of monomineralic and polymineralic assemblages [START_REF] Wheeler | Importance of pressure solution and coble creep in the deformation of polymineralic rocks[END_REF][START_REF] Ford | Modelling interface diffusion creep in two-phase materials[END_REF]. Mineral distribution can become a first-order process to take into account only if it strongly modifies diffusivity or reaction rate. [START_REF] Greene | Role of electrochemical reactions in pressure solution[END_REF] measure an increase in quartz dissolution rate when it is in contact with mica rather than with another quartz crystal. The role of sheet minerals is not fully understood but is related to a change in electrical potential. These experimental results are consistent with observations in natural samples showing the activation of pressure solution by the presence of sheet silicates [START_REF] Heald | Stylolites in Sandstones[END_REF][START_REF] Engelder | Disjunctive cleavage formed at shallow depths in sedimentary rocks[END_REF]. Mineral distribution could thus be a key process to be introduced in further model development but this would require to better capture the physics underlying the catalytic effect of sheet sil-combined with predictions of our model suggest that large sections of the Earth's crust are dry or at least not wet enough to allow for pressure solution. It is critical to estimate the range of fluid content necessary to activate pressure solution. This first requires to better constrain the geometry of the fluid pathways involved in pressure solution. Evidence for dissolution-precipitation creep are mainly collected in metasediments containing hydrous minerals [START_REF] Stöckhert | Low effective viscosity during high pressure metamorphism due to dissolution precipitation creep: the record of HP-LT metamorphic carbonates and siliciclastic rocks from Crete[END_REF]Wassmann & Stöckhert, 2013a;[START_REF] Xia | Structural and rheological evolution of the Laramide subduction channel in southern California[END_REF]. Fluid availability can be transient for example in metamorphic environments where fluid expelled during dehydration is transported towards the surface [START_REF] Connolly | Devolatilization-generated fluid pressure and deformationpropagated fluid flow during prograde regional metamorphism[END_REF][START_REF] Omlin | Pore Fluid Extrac-tion by Reactive Solitary Waves in 3-D[END_REF]. Fluid transport mechanisms over several kilometres involve large fluid pathways whereas the stress sustained by the rock can only be transmitted to the fluid in nanometre-size aqueous films [START_REF] Israelachvili | Intermolecular and Surface Forces[END_REF]. High fluid availability may therefore not favour pressure solution but the quantification of this effect requires estimation of pore distribution up to the nanoscale.

Finally, dissolution-precipitation creep is a transient mechanism since it occurs while reaching equilibrium. It thus requires a driving force for reaction which can be a change in temperature, pressure or chemistry.

Acquisition of GSPO by dissolution-precipitation creep

The observation of a GSPO is commonly used as a marker of dissolution-precipitation creep in rocks [START_REF] Gratier | The Role of Pressure Solution Creep in the Ductility of the Earth's Upper Crust[END_REF]Wassmann & Stöckhert, 2013a). The formation of a GSPO requires that minerals form with an anisotropic shape and that the shapes of the grains are orientated similarly. Several mechanisms can lead to a GSPO: 1) the tortuosity is lower for a diffusion path parallel rather than perpendicular to the schistosity leading to an anisotropic growth if the mineral grows by diffusion. The long axis of the crystals is in the schistosity [START_REF] Brady | Intergranular diffusion in metamorphic rocks[END_REF]; 2) the grain growth kinetics depend on crystallographic orientation [START_REF] Awad | Forsteritic olivine: effect of crystallographic direction on dissolution kinetics[END_REF]; 3) the interfacial energies at the contact between the growing and the surrounding minerals can systematically vary due to a crystallographic preferred orientation of the surrounding minerals [START_REF] Kretz | Interpretation of the shape of mineral grains in metamorphic rocks[END_REF]. This leads to a growth preferred orientation inherited from the host mineral or matrix; 4) the grains can rotate and/or plastically elongate for a significant finite strain of the matrix [START_REF] Willis | A kinematic model of preferred orientation[END_REF][START_REF] Means | Experiments on preferred orientation of platy minerals[END_REF]; 5) as demonstrated here with numerical modelling or in experiments [START_REF] Gratier | Experimental evidence for rock layering development by pressure solution[END_REF], dissolution-precipitation under stress can also favour grain shape preferred orientation. A combination of these mechanisms will promote GSPO as shown for example with a positive feedback between pressure solution and crystallographic preferred orientation [START_REF] Bons | Crystallographic preferred orientation development by dissolution-precipitation creep[END_REF].

We document a pronounced GSPO for garnet and plagioclase in the samples from the Nufenen Pass area with aspect ratios ranging from 1.6 to 2.6 and a similar orientation at the centimetre scale.

Garnet grows isotropically, hence its crystallographic orientation cannot lead to the observed formation of ellipses. Since maximum temperatures never exceeded 500 • C, garnet elongation by plastic deformation is very unlikely [START_REF] Karato | Plastic deformation of garnets: systematics and implications for the rheology of the mantle transition zone[END_REF]. Moreover, the ellipses are obliquely oriented with respect to the schistosity and they all have a similar orientation (Figures 9 and10 andTable 3). The variation of the diffusion coefficient [START_REF] Brady | Intergranular diffusion in metamorphic rocks[END_REF] and/or of the interfacial energy [START_REF] Kretz | Interpretation of the shape of mineral grains in metamorphic rocks[END_REF] with the orientation relative to the schistosity are thus not responsible for the aspect ratios measured here. This suggests a formation after the main foliation event (D2), consistent with the conclusions of [START_REF] Klaper | Alpine metamorphism of pelitic schists in the Nufenen Pass area, Lepontine Alps[END_REF] based on textural observations. The most probable mechanism for the acquisition of GSPO is thus growth under stress. Garnet deformation by pressure solution has already been observed in metamorphic rocks [START_REF] Bell | Sequential porphyroblast growth and crenulation cleavage development during progressive deformation[END_REF][START_REF] Hara | Pressure solution of plagioclase and garnet during Nagahama folding in the Sambagawa belt of Central Shikoku[END_REF][START_REF] Bell | Porphyroblast nucleation, growth and dissolution in regional metamorphic rocks as a function of deformation partitioning during foliation development[END_REF][START_REF] Duebendorfer | Retrogressive dissolution of garnet: Effect on garnet-biotite geothermometry[END_REF][START_REF] Álvarez-Valero | Formation of elliptical garnet in a metapelitic enclave by melt-assisted dissolution and reprecipitation[END_REF][START_REF] Smit | Creep of garnet in eclogite: Mechanisms and implications[END_REF]Wassmann & Stöckhert, 2013a).

-25-Plagioclase commonly develops a GSPO. The growth and dissolution rates of this mineral depends on crystallographic orientation [START_REF] Arvidson | Single-crystal plagioclase feldspar dissolution rates measured by vertical scanning interferometry[END_REF]. However, anisotropic reaction rates are not sufficient to produce a GSPO since the grains also need to have the same orientation at the centimetre scale. GSPO in plagioclase can develop as a result of plastic deformation [START_REF] Stünitz | Dislocation generation, slip systems, and dynamic recrystallization in experimentally deformed plagioclase single crystals[END_REF], at temperature as low as 500 • C [START_REF] Shigematsu | Dislocation creep of fine-grained recrystallized plagioclase under low-temperature conditions[END_REF][START_REF] Rosenberg | Deformation and recrystallization of plagioclase along a temperature gradient: An example from the Bergell tonalite[END_REF]. At low metamorphic grades, plagioclase does not deform by plastic deformation but rather by dissolution-precipitation [START_REF] Shelley | Plagioclase and quartz preferred orientations in a lowgrade schist: the roles of primary growth and plastic deformation[END_REF]. The maximum temperature estimate in the Nufenen pass area is ∼ 500 • C [START_REF] Kamber | Regional metamorphism and uplift along the southern margin of the Gotthard massif; results from the Nufenenpass area[END_REF], that is close from the temperature at which plastic deformation can occur in plagioclase. However, there is no evidence for rotation in the studied samples, suggesting that elongation is not associated with plastic deformation. The plagioclases collected in the Nufenen pass area are elongated in the foliation, suggesting that GSPO is acquired during or after the main deformation event D2 [START_REF] Klaper | Alpine metamorphism of pelitic schists in the Nufenen Pass area, Lepontine Alps[END_REF]. Plagioclases are found in a fold at the contact between a massive carbonate strata and a finely-grained schist (Figure 8 C). This may explain the variation of the GSPO observed at the outcrop scale due to the refraction of schistosity during D2. As the plagioclases are oriented in the foliation, the mechanism responsible for GSPO acquisition is more difficult to determine than for garnet. It could be mass transfer during a dissolution-precipitation process but also anisotropic diffusion or growth in the foliation. There is no relationship between the orientation of the plagioclases and their aspect ratio (Table 3). If the mechanism responsible for GSPO is stress-driven mass transfer, this indicates that the differential stress is rather constant throughout the outcrop during reaction even though the orientation of the main stress components significantly varies at the centimetre scale.

The Dahmköhler number can be calculated with the fluid composition at the equilibrium, the reaction rate, the grain size, the diffusion coefficient and the aqueous film thickness (equation 5). There are some uncertainties associated with these parameters, in particular regarding the reaction rate which can vary over several orders of magnitude (see above). Taking into account these uncertainties, Da α is above 1 at T > 300 • C for the muscovite + quartz <=> andalusite + sanidine + H 2 O reaction, except near the reaction line [START_REF] Schramke | The reaction muscovite + quartz = andalusite + K-feldspar + water; Part 1, Growth kinetics and mechanism[END_REF] (Figure 4).

The acquisition of GSPO is not expected from numerical modelling at high Da α , that is in most of the investigated stress-temperature domain except near the reaction line (Figure 6). The dissolution-precipitation process is indeed controlled by diffusion in most of the stress-temperature domain. This implies that mass transfer between the faces is slow compared to reaction, leading to a precipitation of the secondary phases where the primary phases are dissolved. Thus, the lack of GSPO in rocks does not necessary indicate low stresses but rather suggests that reaction took place at high Da α , in far from the equilibrium conditions. Interestingly, rocks reacting in these conditions (e.g. hydrothermal alteration, contact metamorphism) do not generally display GSPO [START_REF] Roselle | Nucleationdominated crystallization of forsterite in the Ubehebe Peak contact aureole, California[END_REF][START_REF] Gerdes | Stable Isotopic Evidence for Limited Fluid Flow through Dolomitic Marble in the Adamello Contact Aureole, Cima Uzza, Italy[END_REF]. Hydrothermal alteration of peridotite at mid-ocean ridges forms rocks with a mesh texture without GSPO [START_REF] Malvoisin | Control of serpentinisation rate by reaction-induced cracking[END_REF]. Such a reaction occurs in environments where stresses can be high as expected near detachment faults along which mantle rocks are exhumed [START_REF] Olive | Rapid rotation of normal faults due to flexural stresses: An explanation for the global distribution of normal fault dips[END_REF]. Examples of reaction in far from the equilibrium conditions also exist during regional metamorphism. [START_REF] Moore | Stress orientationdependent reactions during metamorphism[END_REF] describe a primary mineralogical assemblage out of equilibrium reacting to form two different secondary mineral assemblages with one stable at higher pressure (A1) than the other (A2). A1 and A2 form under different orientations and have a similar thickness, leading to a limited GSPO. This indicates reaction at high Da α , which is consistent with the high differential stress of 0.2 GPa calculated in [START_REF] Moore | Stress orientationdependent reactions during metamorphism[END_REF] (Figures 4 A andB).

Low Da α can only be reached at low temperature (T < 150 • C) or near thermodynamic equilibrium where reaction rates significantly decrease (Figure 4). Figure 6 indicates that GSPO can only be formed at low Da α . Diffusion indeed proceeds faster than reaction, allowing mass to be transported from one face to the other before precipitation. The metasediments of the Nufenen pass area continuously dehydrate during the prograde path in Barrovian metamorphism. At the same time, temperature increases leading to the overstepping of metamorphic reactions. This overstepping is limited since the rock rapidly reacts by a dissolution-precipitation process in the presence of fluids. As a result, prograde metamorphic reactions are thought to take place close to the equilibrium line [START_REF] Wood | Rates of hydrothermal reactions[END_REF]. When fluid is no longer present, the reaction stops leaving a mineralogical assemblage which can be used to retrieve the pressure-temperature conditions of metamorphism [START_REF] Spear | Metamorphic phase equilibria and pressure-temperature-time paths[END_REF]. Plagioclase and garnet in the metasediments of the Nufenen pass area precipitate during or after the main foliation event (D2), that is close to the peak metamorphism conditions [START_REF] Klaper | Alpine metamorphism of pelitic schists in the Nufenen Pass area, Lepontine Alps[END_REF]. They have thus undergone prograde metamorphism in the presence of fluid during the Alpine collision. We thus interpret GSPO in these samples as the result of replacement under stress near thermodynamic equilibrium. This interpretation is consistent with the results of numerical modelling provided here, showing that the distance from the equilibrium is the key parameter controlling the acquisition of GSPO in rocks.

The acquisition of a GSPO during dissolution-precipitation requires the presence of fluid serving as a transport medium for the dissolved minerals. GSPO is observed here in volatile-rich rocks [START_REF] Bestmann | Development of single-crystal σ-shaped quartz porphyroclasts by dissolutionprecipitation creep in a calcite marble shear zone[END_REF][START_REF] Imon | Development of shapeand lattice-preferred orientations of amphibole grains during initial cataclastic deformation and subsequent deformation by dissolutionprecipitation creep in amphibolites from the Ryoke metamorphic belt, SW Japan[END_REF][START_REF] Wassmann | Dissolution precipitation creep versus crystalline plasticity in high-pressure metamorphic serpentinites[END_REF][START_REF] Stokes | Deformation of amphibolites via dissolution-precipitation creep in the middle and lower crust[END_REF][START_REF] Bial | The microstructural record of porphyroclasts and matrix of partly serpentinized peridotite mylonites from brittle and crystal-plastic deformation to dissolutionprecipitation creep[END_REF]Wassmann & Stöckhert, 2013a) and could be used as a marker for fluid presence during metamorphism. In addition, our numerical simulations clearly show that GSPO can only be formed near thermodynamic equilibrium. This suggests that GSPO could be used to probe the distance from the equilibrium during metamorphism. The distance from the equilibrium plays a key role on the kinetics of reaction [START_REF] Lasaga | Kinetic Theory in the Earth Sciences[END_REF] but also, as shown below, on the barometric record during metamorphism.

Pressure controlling metamorphic reactions

Determining the pressure of metamorphic reactions is crucial for interpreting observations in natural samples, thermodynamic calculations and numerical models [START_REF] Moulas | Relation between mean stress, thermodynamic, and lithostatic pressure[END_REF]. It serves as a basis in the metamorphic petrology community to reconstruct the geodynamic history of the rocks collected in the field [START_REF] Spear | Metamorphic phase equilibria and pressure-temperature-time paths[END_REF]. In the presence of fluid, experiments [START_REF] Llana-Fúnez | Metamorphic reaction rate controlled by fluid pressure not confining pressure: implications of dehydration experiments with gypsum[END_REF], theoretical studies [START_REF] Dahlen | Metamorphism of nonhydrostatically stressed rocks[END_REF] and field data [START_REF] Holdaway | Rock pressure vs. fluid pressure as a controlling influence on mineral stability: an example from New Mexico[END_REF] suggest that fluid pressure is the relevant pressure for determining thermodynamic equilibrium. Fluid pressure in the aqueous film varies at grain boundary due to interactions with solid surfaces [START_REF] Israelachvili | Intermolecular and Surface Forces[END_REF]. It is equal to the normal stress in thin fluid films leading to strong pressure gradients at the micrometre scale between differently oriented mineral surfaces [START_REF] Gratier | The Role of Pressure Solution Creep in the Ductility of the Earth's Upper Crust[END_REF]. [START_REF] Moore | Stress orientationdependent reactions during metamorphism[END_REF] document an orientation dependence for precipitation at grain boundary during metamorphism. [START_REF] Wheeler | Dramatic effects of stress on metamorphic reactions[END_REF] proposes to calculate the affinity between differently oriented grain boundaries (reaction pathway) to estimate the effect of stress on thermodynamic equilibrium. Mineral stability fields shift by several hundreds of megapascals with this approach. However, [START_REF] Fletcher | Dramatic effects of stress on metamorphic reactions[END_REF] indicates that one cannot freely select a reaction pathway without taking into account reaction kinetics and diffusion. The model proposed here considers these two processes and can thus be used to estimate the effect of stress on pressure record during metamorphism.

The effect of stress on reaction can first be evaluated by trying to conceptually reproduce the approach of Wheeler (2014) with our model. [START_REF] Wheeler | Dramatic effects of stress on metamorphic reactions[END_REF] considers an incongruent pressure solution pathway to estimate the displacement of reaction line with stress. This is equivalent to consider that dissolution and precipitation occur in the aqueous film normal to σ 1 and σ 3 , respectively (low Da α ). The reaction line is defined as the zero chemical affinity curve (A t,1,3 = 0 or c e,r,1 = c e,p,3 with our formalism). By definition, the dissolving phases have a higher solubility than the precipitating phases at the same stress, so that c e,r,1 > c e,p,1 and c e,r,3 > c e,p,3 . Moreover, solubilities always increase with pressure except near the critical point of water leading to c e,r,1 > c e,r,3 and c e,p,1 > c e,p,3 . Combining these four latter inegalities lead to c e,r,1 > c e,p,3 , which is not compatible with the condition for zero chemical affinity. As a result, the incongruent pressure solution pathway used in Wheeler (2014) cannot be reproduced with our pressure solution model and pressure solution is thus not expected to significantly modify thermodynamic equilibrium. The only way to reproduce Wheeler (2014)'s approach is to consider that dissolution and precipitation occur in the aqueous film normal to σ 3 and σ 1 , respectively (c e,p,1 = c e,r,3 ). This requires a mechanism generating a large supersaturation under σ 3 before transport and precipitation under σ 1 . Such a mechanism is known as force of crystallisation [START_REF] Steiger | Crystal growth in porous materialsI: The crystallization pressure of large crystals[END_REF].

The model proposed here also provides a way to quantitatively estimate the effect of stress on equilibrium. Mass transfer between differently oriented grain surfaces can only occur at low Da α (fast diffusion compared to reaction). The estimates of Da α in Figure 4 are high except near the reaction line. Moreover, the reaction rate does not depend on stress at high A t,α,β (equation 1). Differential stress is thus not expected to have an effect on thermodynamics in far from the equilibrium conditions (Figure 12). The numerical simulations performed here indicate that, near the equilibrium, reaction is affected by differential stress. σ 1 controls the reaction if A t,α,β increases with pressure (∆ r V < 0). σ 3 controls the reaction otherwise (∆ r V > 0). These results are in agreement with the experiments performed by [START_REF] Richter | Stresses and pressures at the quartz-to-coesite phase transformation in shear deformation experiments[END_REF], showing that the quartz to coesite transition (∆ r V < 0) depends on σ 1 . Thermodynamic equilibrium recorded in mineralogical assemblages is thus strongly dependent on the relationship between reaction rate and A t,α,β .

Conclusion

We propose here a model for pressure solution and dissolution-precipitation creeps.

This model provides revised deformation mechanism maps for pressure solution. The previous deformation mechanism maps were not incorporating the effect of concentration changes in the aqueous film on reaction rate. We show the critical role of this dependency for estimating strain rates. We also propose new deformation mechanism maps and creep laws for dissolution-precipitation creep. Dissolution-precipitation creep is the main deformation mechanism in the Earth's crust provided that fluid is available and disequilibrium is achieved. The refinement of these deformation mechanism maps requires a detailed investigation of the dependence of reaction rate on pressure and water activity in thin fluid films. We also use numerical modelling to gain insights on the acquisition of GSPO during metamorphism and the pressure controlling metamorphic reactions. We show here that the distance from the equilibrium plays a key role for these two topics.

The model is indeed controlled by diffusion in far from the equilibrium conditions, leading to limited mass transfer between the reacting mineral faces. As a result, GSPO does not develop and the effect of stress on thermodynamic pressure is small. In close from the equilibrium conditions, significant differences in reaction rate occur between the mineral faces, leading to a control of the reaction by σ 1 or σ 3 . Moreover, the mineral grains change of shape and a GSPO is formed, in agreement with observations in samples from the Nufenen pass area. These results can serve as a basis to interpret GSPO and barometric record in metamorphic rocks.

Appendix A Fit of the experiments of quartz compaction

We determined k r by fitting with equation 18 the dataset of [START_REF] Niemeijer | Compaction creep of quartz sand at 400-600 • C: experimental evidence for dissolution-controlled pressure solution[END_REF] for quartz grains compaction in pure water under stresses ranging from 50 to 150 MPa and at temperature comprised between 400 and 600 • C. Strain rate data are reproduced from their figure 4. k rbulk ranges from 3.3.10 -6 mol.m -2 .s -1 at 400 • C to 1.4.10 -4 mol.m -2 .s -1 at 600

• C [START_REF] Tester | Correlating quartz dissolution kinetics in pure water from 25 to 625 • C[END_REF]. The grain sizes are given in their Table 1. The experiments are loaded with sands and the grain contacts thus evolve during compaction. To take into account grain contact evolution, a stress concentration factor depending on porosity (B (φ)) is introduced in the equations in front of the (σ α -P f,0 ) term [START_REF] Dewers | Rate laws for water-assisted compaction and stress-induced water-rock interaction in sandstones[END_REF]. This factor theoretically tends towards infinity at the beginning of the compaction even though a small elastic distortion occurs leading to finite grain contact (X. [START_REF] Zhang | Compaction creep of wet granular calcite by pressure solution at 28 • C to 150[END_REF]. B is modelled here by using a body centered cubic packing of spherical grains for which:

B (φ) = d (φ) 2 π L 2 2 -d(φ) 2 2 (A1)
where d is the diameter of the grain-to-grain contact determined by numerically solving:

8 3 π L 2 3 - 4 3 π (L -d) 2 L + d 2 -(1 -φ) 2d √ 3 3 = 0 (A2)
B is approximately 40 for a ratio between actual porosity and initial porosity of 0.9. This corresponds to an effective pressure of 4 GPa at grain contact for an effective pressure of 100 MPa applied on the sample. The reaction rate constant theoretically varies with pressure as [START_REF] Schott | The Link Between Mineral Dissolution/Precipitation Kinetics and Solution Chemistry[END_REF] Initial porosity during compaction experiments none µ f,{r,p},α Chemical potential of the mineralogical components dissolved in the fluid J.mol -1 µ s,{r,p},α Chemical potential of the solid mineralogical assemblage J.mol -1 A {r,p},α Chemical affinity of the dissolution or precipitation reactions in the aqueous film normal to x α (µ s,{r,p},α -µ f,{r,p},α )

J.mol -1

A t,α,β Chemical affinity of the replacement reaction (µ s,r,α -µ s,p,β ) J.mol -1 R Gas constant J.K -1 .mol

  It is the opposite at boundary x = L. This flux can then be introduced into the equivalent of equation 3 in one dimension under the assumptions of steady state to obtain an expression of c 0 : that S 0 = S 0,r = S 0,p .

  and ∆ r V = V p -V r are, respectively, the overall reaction rate constant and the solid change in volume for the replacement reaction including both dissolution and precipitation. ˙ αV is thus the strain rate when the system is controlled by volume change during reaction. When the reaction rate is controlled by dissolution, the reaction rate constant of the dissolving phase is smaller than the one of the precipitating phase (k r k p ), and k m = k r 1 -ce,p,α ce,r,α . When the reaction rate is controlled by precipitation (k p k r ), k m = -k p 1 -ce,r,α ce,p,α . In the case of isotropic stresses, equation 9 is equivalent to a kinetically-delayed change in solid volume.

  and the merged high pressure databases of[START_REF] Sverjensky | Water in the deep Earth: The dielectric constant and the solubilities of quartz and corundum to 60 kb and 1200 • C[END_REF] and[START_REF] Holland | An improved and extended internally consistent thermodynamic dataset for phases of petrological interest, involving a new equation of state for solids[END_REF] for dissolved aqueous species and solid phases, respectively. k r is poorly constrained since kinetic parameters are generally determined in a bulk fluid with experiments on powders in batch reactors in which the reaction processes are probably strongly different than in an aqueous film. The quartz indentation experiments of[START_REF] Gratier | A pressure solution creep law for quartz from indentation experiments[END_REF] provide precise strain rate measurements at 350 • C in a NaOH-bearing fluid and at stress ranging from 25 to 350 MPa. Their experimental dataset was fitted with equation 18 to estimate k r (Figure2 A). We obtained a value of k r = 3.10 -5 mol.m -2 .s -1 which is four orders of magnitude lower that the rate of quartz dissolution in a bulk fluid estimated with the laws of[START_REF] Tester | Correlating quartz dissolution kinetics in pure water from 25 to 625 • C[END_REF] and[START_REF] Worley | Quartz dissolution kinetics from 100-200 • C as a function of pH and ionic strength[END_REF] in the same conditions (k rbulk ).The experiments of[START_REF] Niemeijer | Compaction creep of quartz sand at 400-600 • C: experimental evidence for dissolution-controlled pressure solution[END_REF] of quartz grain compaction can also be used to compare k r to the dissolution laws obtained with experiments on powders. They were performed with stresses ranging from 50 to 150 MPa and at temperature comprised between 400 and 600 • C. Estimating stress (and associated reaction progress) at grain contacts from the applied far-field stress in the experiments of[START_REF] Niemeijer | Compaction creep of quartz sand at 400-600 • C: experimental evidence for dissolution-controlled pressure solution[END_REF] requires introduction of a stress concentration factor (B)[START_REF] Dewers | Rate laws for water-assisted compaction and stress-induced water-rock interaction in sandstones[END_REF] (see details on the calculation in Appendix A). This factor decreases with porosity during compaction as more grain contacts are formed (Figure2 B). Therefore, estimates of reaction rate in the thin aqueous film are more indirect than with the dataset of[START_REF] Gratier | A pressure solution creep law for quartz from indentation experiments[END_REF]. They are at least one order of magnitude lower than in the bulk fluid (Fig-ure 2 B).The fit of experimental data thus indicate k r values below the one determined in bulk fluid. Since the uncertainty on k r is high (1 to 4 orders of magnitude below k rbulk ), we investigated the influence of this parameter on strain rate prediction by determining two deformation mechanism maps for quartz with either k r = k rbulk or k r = k rbulk / 10 4 (Figure3C and D). The Damköhler number calculated with equation 5 and the parameters given above are required to plot the deformation mechanism maps. It is displayed in Figure3A and B for distinguishing between diffusion-and reaction-controlled pressure solution creep. The deformation mechanism maps are obtained by considering, in addition to pressure solution, dislocation and diffusion creeps for which strain rates are calculated with the laws given in[START_REF] Tokle | Flow laws and fabric transitions in wet quartzite[END_REF] and[START_REF] Rutter | Experimental grain size-sensitive flow of hotpressed Brazilian quartz aggregates[END_REF], respectively (see Appendix B for details on these laws). .

Figure 2 .

 2 Figure 2. Fit of experimental data with equation 18 to determine kr. A: Strain rate as a function of applied normal stress in the quartz indentation experiments of[START_REF] Gratier | A pressure solution creep law for quartz from indentation experiments[END_REF]

  650 • C whereas pressure solution creep dominates elsewhere. The pressure solution creep domain extends as the reaction rate increases. The strain rate calculated for pressure solution are two orders of magnitude higher when calculated with k r = k rbulk than with k r = k rbulk / 10 4 , in agreement with equation 18. Pressure solution creep is controlled by reaction at low temperature whereas it is controlled by diffusion at high temperature. The temperature of transition from reaction-to diffusion-controlled pressure solution creep is ∼ 70 • C and ∼ 350 • C for k r = k rbulk and k r = k rbulk / 10 4 , respectively.

m 3 .Figure 3 .

 33 Figure 3. Damköhler number and deformation mechanism maps for quartz as a function of stress and temperature. The calculations are performed with an aqueous film thickness of w = 2.10 -9 m, a grain size of L = 100µm, a diffusion coefficient given in equation 22, and a solubility calculated with equation 10. Two different reaction rates are used: kr = k rbulk (A and C) and kr = k rbulk /(10 4 ) (B and D). A and B: calculated Damköhler number. The contours correspond to log10 (Daα). C and D: deformation mechanism maps. The contours correspond to log10 ( ˙ ). The domains where one of the four studied deformation mechanism (reaction-and diffusion-controlled pressure solution creep, dislocation creep and diffusion creep) dominates are displayed. The strain rates for dislocation and diffusion creeps are calculated with the equations given in Tokle et al. (2019) and[START_REF] Rutter | Experimental grain size-sensitive flow of hotpressed Brazilian quartz aggregates[END_REF], respectively. The strain rate for pressure solution creep is calculated with equation 18. Reaction-and diffusion-controlled domains correspond to Daα < 1 and Daα > 1, respectively.

Figure 4 D

 4 Figure4 D). The strain rate calculated for dissolution-precipitation creep are four orders of magnitude higher for k r = k rbulk than for k r = k rbulk /(10 4 ), in agreement with equation 9. The calculations do not predict diffusion creep to be a dominant deformation mechanism. The main difference with pressure solution (Figure3) is observed near thermodynamic equilibrium (i.e. near the reaction line) where dislocation creep can dominate due to the low reaction rate calculated in this region.

Figure 4 .

 4 Figure 4. Damköhler number and deformation mechanism map for a system in which the replacement reaction of muscovite + quartz <=> andalusite + sanidine + H2O occurs. The calculations are performed with an aqueous film thickness of w = 2.10 -9 m, a grain size of L = 100µm, a diffusion coefficient given in equation 22, and a solubility calculated with equation 10. Two different reaction rates are used: kr = k rbulk (A and C) and kr = k rbulk /(10 4 ) (B and D). A and B: calculated Damköhler number. The contours correspond to log10 (Daα). C and D: deformation mechanism map. The strain rates for dislocation and diffusion creeps are calculated with the equations given in Tokle et al. (2019) and[START_REF] Rutter | Experimental grain size-sensitive flow of hotpressed Brazilian quartz aggregates[END_REF], respectively. The strain rate for dissolution-precipitation creep is calculated with equation 9, which is relevant for high Daα and does not depend on differential stress (control by volume change during reaction).

Figure 5 .Figure 6 .

 56 Figure 5. Numerical modelling of replacement reaction under stress for kr = kp = 1.10 -10 mol.m -2 .s -1 . A: concentration in the aqueous film normal to σ1 at a reaction progress of ξ = 10%. B: concentration in the aqueous film normal to σ2 at ξ = 10%. C: concentration in the aqueous film normal to σ3 at ξ = 10%. D: grain shape at ξ = 10% deduced from reaction rate calculated with the numerical model for concentration in the thin fluid film. The initially rock-forming mineral is displayed in yellow and the reaction products are displayed in green. The orientation of the stress components are indicated. E: grain shape at ξ = 40%. F: grain shape at ξ = 90%.

Figure 7 .

 7 Figure 7. Proportion of the overall reaction progress on each group of mineral faces with the same orientation as a function of the distance from the equilibrium. The faces normal to σ1, σ2

Figure 9 .

 9 Figure 9. Plagioclase and garnet porphyroblasts orientation for the Nufenen pass area on an equal area lower hemisphere projection. The schistosity measured in the field is in green, the foliation measured with ilmenite orientation in black and the mean porphyroblast shortening plane in red. When the sedimentary layering (S0) is different from the schistosity, it is displayed in blue (samples LB 2018 1, LB 2018 2, LB 2018 3 and LB 2018 3 2). The white dots are the principal elongation axes of the porphyroblasts. The dot density is contoured. Sample LB 2018 10 was not collected in place and orientations are thus given relative to an arbitrary plane chosen as the planar bottom of the core.

Figure 10 .

 10 Figure 10. 3D images of segmented cores. A: plagioclase-bearing sample (LB 2018 3 2). B: garnet-bearing sample (LB 2018 6). Garnet, plagioclase and ilmenite are displayed in yellow, green and blue, respectively. The ilmenite crystals serve to determine the orientation of the foliation.

Figures 3 CFigure 11 .

 311 Figures 3 C and D and 4 C and D provide new deformation mechanism maps including the pressure solution and dissolution-precipitation creep laws determined here.As shown above, these creep laws involve less assumptions than the previous formalisms.Moreover, equation 18 encompasses all Da α values allowing strain rate to be computed without an a priori assumption on the process controlling deformation (diffusion or reaction). Plotting deformation mechanism maps requires quantitative estimates of the parameters used in the creep laws. Significant efforts have been devoted to the determination of the thickness of the aqueous film and the diffusion coefficient[START_REF] Revil | Pervasive pressure solution transfer in a quartz sand[END_REF] Niemei-

ure 3 C

 3 and D) and dissolution-precipitation (Figure 4 C and D) with reaction rates equal or four orders of magnitude smaller than the one determined in the bulk fluid. For quartz, the transition from reaction-to diffusion-controlled pressure solution is predicted to occur between ∼ 70 • C (k r = k rbulk ; Figure 3 C) and ∼ 350 • C (k r = k rbulk / 10 4 ; Figure 3 D). This is consistent with the previous temperature estimate of ∼ 150 • C for this transition

3 Figure 12 .

 312 Figure12. Sketch showing the thermodynamic pressure for metamorphic reactions under stress in a pressure-temperature domain. The plain line is the reaction line where equilibrium between reactant and products is achieved. Far from the equilibrium, the thermodynamic pressure is the mean stress (σmean). Near the equilibrium, the thermodynamic pressure depends on the total volume change during reaction. It is σ1 for ∆rV < 0 and σ3 for ∆rV > 0.

3 D 2 w 3 .

 323 -3 c e,{r,p}α Concentration at the equilibrium of dissolved aqueous species in the aqueous film normal to x α mol.m -3 c 0 Concentration of dissolved aqueous species in the pores mol.m -3 c e,{r,p},0 Concentration at the equilibrium of dissolved aqueous species in the pores mol.m -3 c b,α Concentration when balance between the amounts of dissolving and precipitating minerals in the aqueous film normal to x α . Defined below equation 4. mol.m -Diffusion coefficient in the aqueous fluid film m 2 .s -1 S {r,p},α Reactive surface area in one aqueous film normal to x α m 2 S 0,{r,p} Reactive surface area in one pore m Thickness of the aqueous film m x Length along the grain boundary in 1D m L Length of the grain boundary in 1D m d Diameter of the grain-to-grain contact defined in equation A2 m Da α Damköhler number in the aqueous film normal to x α . Defined in equation 5. none Da 0 Damköhler number in the pores. Defined below equation 7. none V {r,p} Molar volume m 3 .mol -1 ∆ r V Change in volume during reaction (V p -V r ) m

  α c b,α Daα α Daα c = c b,α c 0 = c b,0 c α = c b,α

n

  ilm and n pb are the number of analyzed ilmenite and porphyroblast grains, respectively. F/MSP is the angle between the foliation orientation measured with ilmenite and the main shortening axis. F/MEA is the angle between the foliation orientation measured with ilmenite and the main elongation axis. *: sample LB 2018 10 was not collected in place and orientations are thus given relative to an arbitrary plane chosen as the planar bottom of the core. An other equation for concentration is generally considered by assuming the achievement of thermodynamic equilibrium in the aqueous film by integrating equation 10 over the aqueous film length (equation 3 from Rutter (1976)): L 0 c α (x)dx = c e,r,α L (C3) Equations C2 and C3 deal with the same single variable, c α (x). They are generally combined to express one of the model parameters (usually the reaction rate, k r ) as a function of the others: of equation C2 at x = 0 into equation 13 and replacing k r with equation C4 provide the following equation for strain rate : ˙ αS = -24c e,r,0 DwV r L 3 1 -exp (σ α -P f,0 ) V r RT ≈ 24c e,r,0 DwV 2 r (σ α -P f,0 ) RT L 3 (C5)

  -1 ) is modelled under stress σ α with the kinetic law of Lasaga (2014):

	ξ{r,p},α = k {r,p} S {r,p},α 1 -exp	-A {r,p},α RT	= k {r,p} S {r,p},α 1 -	a α a e,{r,p},α
	= k {r,p} S {r,p},α 1 -	c α c e,{r,p},α

  from the equilibrium constants of muscovite + quartz or andalusite + sanidine + H 2 O dissolution reactions written with SiO 2,aq , AlO - 2 , K + and H + as dissolved aqueous species. The molar volumes are 162.68.10 -6 m 3 .mol -1 and 160.38.10 -6

  : Subscript referring to the dissolving phases p Subscript referring to the precipitating phases σ α Stress component α with α = 1, 2 or 3 and σ 1 > σ 2 > σ 3 Reaction rate in the aqueous film normal to x α mol.s -1 k {r,p} Reaction rate constant mol.m -2 .s -1 k m Mean reaction rate constant defined below equation 9 mol.m -2 .s -1 k rbulk Reaction rate constant in the bulk fluid mol.m -2 .s -1 a α Activity of dissolved aqueous species in the aqueous film normal to x α none a e,{r,p}α Activity at the equilibrium of dissolved aqueous species in the film normal to x α none c α Concentration of dissolved aqueous species in the aqueous film normal to x α mol.m

	Symbol Definition	Unit
	r
			P a
	P f,0	Pore pressure	P a
	x α	Axis vector parallel to σ α	m
	ξ{r,p},α

Table 1 .

 1 Strain rate in direction x α s -1 ˙ αV Strain rate in direction x α when control by volume change during reaction s -1 ˙ αD Strain rate in direction x α when control by diffusion s -1 ˙ αR Strain rate in direction x α when control by reaction s -1 ˙ αS Strain rate in direction x α derived in a similar way to Rutter (1976) s -1 ˙ Symbols used in the mathematical model

	-1

dif Pre-exponential factor for diffusion creep m 2 .P a -1 .s -1 Ea dif Energy of activation for diffusion creep J.mol -1

Table 2 .

 2 Borderline cases for cα and c0 as a function of the Da numbers in the pore (Da0) and in the aqueous film (Daα).

	LB 2018 1	Plagioclase	348/73 1126	140/74 32	156/41	5/62	62	1 2.0 ± 0.4
	LB 2018 2	Plagioclase	348/73 1417	150/90 40	181/66 169/86	20	12 2.2 ± 0.4
	LB 2018 3	Plagioclase	337/81 1182	112/85 40	266/67 222/70	72	5 2.4 ± 0.6
	LB 2018 3 2	Plagioclase	337/81	814	304/85 25	326/79 216/78	87	1 2.0 ± 0.5
	LB 2018 4	Plagioclase	138/89	839	279/87 13	270/59 247/81	32	7 2.3 ± 0.7
	LB 2018 5	Plagioclase	138/89 1539	108/90 17	282/86 227/86	61	0 2.4 ± 0.7
	LB 2018 11	Plagioclase	131/89	991	312/83 42	350/68 325/80	13	7 2.6 ± 0.4
	LB 2018 6	Garnet	352/78 1822	356/76 75	330/76 158/88	24	19 1.8 ± 0.4
	LB 2018 7	Garnet	352/81 1480	352/83 67	207/63 169/73	25	22 1.6 ± 0.3
	LB 2018 7 2	Garnet	352/81 1471	354/84 96	236/77 158/77	25	18 1.8 ± 0.4
	LB 2018 8	Garnet	350/79 1374	342/82 40	200/79 165/84	15	15 2.1 ± 0.5
	LB 2018 9	Garnet	350/79 1364	348/79 40	197/80 171/85	16	15 2.3 ± 0.9
	LB 2018 10	Garnet	NA*	654	126/30* 65 132/11*	5/15*	39	7 1.7 ± 0.2

Sample

Porphyroblast S1(field)

n ilm F(ilmenite) n pb MEA MSP F/MSP ( • ) F/MEA ( • ) Aspect ratio

Table 3 .

 3 Results of the statistical analysis for porphyroblasts orientation (values are given as Azimuth/Dip). S1(field), F(ilmenite), MEA and MSP refer to the trend/plunge of the Schistosity measured in the field, the Foliation measured with the ilmenite main shortening plane, the Mean main Elongation Axis of the porphyblast and the Mean main Shortening Plane of the porphyroblast, respectively. The mean orientations are determined with Bingham distribution.

icates.

Secondly, we assume that free fluid is available in the system during deformation.

Fluid is probably not commonly present in the Earth's crust explaining why quartz commonly deforms by dislocation creep in natural settings at temperatures above 300 • C (Stipp et al., 2002;Trepmann & Stöckhert, 2003;Kilian et al., 2011). Similarly, diffusion creep is expected to occur at high temperature if fluid availability is low. These observations -24-manuscript submitted to Geochemistry, Geophysics, Geosystems Acknowledgments B.M. acknowledges support from the Swiss National Science Foundation (Grant No. PZ00P2 168083).

The mathematical model can be reproduced following derivations specified in the text of the paper. We thank J.-P. Gratier and an anonymous referee for reviews, W. Behr for editorial handling, and J. Wheeler and J. Platt for discussions.

(C) Plagioclase-porphyroblast (Pl) bearing schist. Plagioclases are mainly form at the lithological contacts. The main elongation axis is displayed with white dashed lines for some porphyroblasts.

where ∆V * is the volume of activation (volume of the activated complex minus volume of the reactant). ∆V * is approximately -1.5 cm 3 .mol -1 at 200-300 • C and several tens of megapascals (R. [START_REF] Zhang | Dissolution kinetics of quartz in water at high temperatures across the critical state of water[END_REF]. Using this volume of activation, the reaction rate of quartz dissolution should increase with pressure by a factor ∼ 1.7, ∼ 13 and ∼ 3.10 4 for stresses of 200 MPa, 1 GPa and 4 GPa, respectively. Equation A3 allows to correct reaction progress for the effect of pressure.

Appendix B Dislocation and diffusion creeps for the deformation mechanism maps

Dislocation creep is modelled following [START_REF] Tokle | Flow laws and fabric transitions in wet quartzite[END_REF] as:

with A dis a material parameter, n is a stress exponent, f H2O is the water fugacity, h is the water fugacity exponent and Ea dis is the activation energy for dislocation.

We used the parameters from the extrapolated fits of Tokle et al. ( 2019) (their Table 1).

We take into account the change of stress exponent in the flow law (from n = 4 to n =

3) when stress increases. We use a linear transition from one flow law to the other between 100 and 400 MPa. Water fugacity is calculated with the Pitzer and Sterner (1994)'s equation of state.

Diffusion creep is modelled following [START_REF] Rutter | Experimental grain size-sensitive flow of hotpressed Brazilian quartz aggregates[END_REF]:

where A dis is a material parameter and Ea dif is the activation energy for diffusion creep.

Appendix C Creep law for pressure solution with a derivation similar to [START_REF] Rutter | The Kinetics of Rock Deformation by Pressure Solution [and Discussion[END_REF] We derive in the following a creep law in the diffusion-controlled case in a similar way to [START_REF] Rutter | The Kinetics of Rock Deformation by Pressure Solution [and Discussion[END_REF] but with the formalism and the geometry adopted here. The two main assumptions used in [START_REF] Rutter | The Kinetics of Rock Deformation by Pressure Solution [and Discussion[END_REF] are that the reaction rate is constant ( ξr,α = k r S r,α ) and that thermodynamic equilibrium is achieved in the thin aqueous film (c α = c e,r,α ).

The assumption of constant reaction rate simplifies equation 2 at steady state in one dimension as:

which can be solved to express c α as a parabolic function: