Equations for passivation

Two approaches can be considered to model a passivation process (Fig. 1). The first involves meshing the space for the explicit modeling of transport mechanisms. While useful for a deeper understanding of the transport process, such modeling becomes numerically complex when transport mechanisms in the fluid are taken into account. Indeed, the concentration profiles in the solid for each fluid cell position in space have to be held in memory by the code. The second approach only describes the transport flux arising in the fluid. Since it is numerically simpler, it enables computational effort to be focused on the effect of the transport flux on the fluid chemistry and the way in which that can drive the dissolution of the protective layer, which is the key long-term mechanism. Including the transport through the layer is essential for modeling glass dissolution experiments effectively as it has a major influence on both the evolution of pH as a function of time and glass surface area to solution-volume ratio.

Fig. 1 Two possible approaches to modeling passivation: defining a meshing space for an explicit transport model (1) or simply modeling the flux of elements arising in the fluid (2) Full size image

In the GRAAL model, the glass alteration rate is a function of the protective layer’s transport properties, that is, its thickness and apparent diffusion coefficient. The thickness of the protective layer is function of a creation term, r 1 , at the interface between the primary solid and the protective layer and a dissolution term, r 2 , at the external face of the protective layer (Fig. 2). The coupling between these terms was solved analytically by ref. 32. A first simplification is discussed in ref. 21, but it does not allow for any feedback effect of the advancing glass hydration front on the ion diffusion rate within the protective layer. Equations and numerical resolutions are presented in the next sections.

Fig. 2 Kinetics relative to the protective layer: Creation kinetic (r 1 ) and dissolution kinetic (r 2 ) of the protective layer (PL) Full size image

Experiments on borosilicate glasses show that the thicker the protective layer (PL) is, the slower is the alteration rate of the primary solid (PS)56 and that, consequently, the protective layer grows with the square root of time when not undergoing dissolution at its surface. This simple experimental observation is expressed as (Eq. 1). r 1 is the primary solid alteration rate controlled by diffusion (m s−1), x PL is the protective layer’s thickness (m), and A is a constant of proportionality (m2 s−1).

$$r_1 = \frac{A}{{x{_{{{\rm{PL}}}}}}}.$$ (1)

Equation 2 defines r 1 , thus:

$$r_1 = \frac{{{\mathrm{d}}x{_{{{\rm{PL}}}}}}}{{{\mathrm{d}}t}}.$$ (2)

Combining Eqs. 1 and 2, considering x = 0 at t = 0, and assuming no dissolution at the outer surface, Eq. 3 can be stated as:

$$x{_{{{\rm{PL}}}}} = \sqrt {2At}.$$ (3)

This equation is to be compared with the resolution of Fick’s second law (Eq. 4) with the following hypotheses:

(H 1 ) diffusion of a glass mobile element of concentration C (g/m 3 ),

(H 2 ) in a semi-infinite media defined on one side by the glass/fluid initial interface and on the other side by the infinite glass,

(H 3 ) considering that the driving force is the concentration gradient between the mobile element concentration in the pristine glass C 0 (g/m 3 of glass) and zero at the interface with the fluid,

(H 4 ) considering a constant diffusion coefficient in the semi-infinite media D (m2 s−1).

$$\frac{{\partial C}}{{\partial t}} = D\frac{{\partial ^2C}}{{\partial x^2}}.$$ (4)

With those initial conditions, Eq. 5 is the solution of Eq. 4.

$$C\left( {x,t} \right) = \,C_0\,erfc\left( {\frac{x}{{2\sqrt {Dt} }}} \right).$$ (5)

Following Fick’s first law, derivation of Eq. 5 with respect to x and multiplication by −D enables the flux of the mobile element to be calculated. Integration with time of the flux taken in x = 0, gives the concentration in the fluid C f (t) (g/m3 of solution) of the mobile element going through the external surface S (m2) and diluting in a perfectly steered solution of volume V (m3) (Eq. 6).

$$C_{\mathrm{f}}(t) = \frac{S}{V}\,2\,C_0\sqrt {\frac{{Dt}}{\pi }}.$$ (6)

The “equivalent thickness of alteration,” with respect to a glass constituent, is given by a simple mass balance (Eq. 7). It is the thickness of glass that has been altered to explain the concentration C f of a glass constituent in solution.

$$E{\mathrm{th}} = \frac{{C_{\mathrm{f}}V}}{{C_0S}}\,.$$ (7)

Equation 8 follows from Eqs. 6 and 7:

$$E{\mathrm{th}} = \frac{2}{{\sqrt \pi }}\sqrt {Dt}.$$ (8)

Equations 8 and 3 are equivalent. However, the proportionality constant, D in Eq. 8, is more meaningful than A in Eq. 3, in agreement with hypotheses H 1 to H 4 . Therefore, Eq. 8 is preferred for the GRAAL model. In agreement with Eq. 8, Eq. 1 becomes Eq. 9.

$$r_1 = \frac{{D\frac{\pi}{2}}}{{x_{{{{\rm{PL}}}}}}}.$$ (9)

According to Eq. 9, the time required to create a protective layer from glass when its thickness is zero, is infinite. However, there is no physical reason for the dissolution rate of mobile ions to be infinite. Therefore, a constant kinetic, r h , known as the hydration rate, is added in order to limit the creation kinetics of the protective layer (Eq. 10). This rate is meant to be higher than the initial dissolution rate in agreement with the higher dissolution rate of mobile ions in comparison to silicon ions.

$$r_1 = \frac{{r_{\mathrm{h}}}}{{1 + x{_{{{\rm{PL}}}}}\frac{{r_{\mathrm{h}}}}{{D\frac{\pi }{2}}}}}.$$ (10)

The hydration rate is negligible as soon as \(x_{\mathrm{{PL}}} > \frac{{D\frac{\pi }{2}}}{{r_{\mathrm{h}}}}\), and more so if \(x_{\mathrm{{PL}}} > \frac{{D\frac{\pi }{2}}}{{r_0}}\). Whatever the pH and temperature, numerical modeling using nuclear glass parameters shows that the hydration rate is largely negligible for a 1 nm thick layer. The conclusion is that the parameter is not measurable on nuclear glasses: whatever the value of r h > r 0 , modeling results after a few seconds of alteration do not change. However, the parameter is useful from a numerical point of view as it avoids the error that arises from dividing a number by zero.

Equation 11 expresses the pH and temperature dependence of the dissolution rate of the protective layer, and the exponential rate drop when the concentrations approach saturation, following the formalism proposed by Aagaard and Helgeson.20

$$r_2 = k_ + \left( {H^ + } \right)^ne^{\frac{{ - E_{\mathrm{a}}}}{{RT}}}\left( {1 - \frac{Q}{K}} \right),$$ (11)

where, k + is the dissolution rate of the protective layer in pure water (forward rate) (m d−1), n is the pH-dependence coefficient of the initial rate, E a is the apparent activation energy at the initial rate (J mol−1), Q is the activity product of the protective layer, K is the activity product at saturation of the protective layer, and R is the ideal gas constant (J mol−1 K−1).

Case of r 1 » r 2 : the thickness of the diffusion zone increases with the square root of tine and the flux of mobile species released from the glass decreases with the inverse square root of time. This occurs at initial time and once concentrations of elements belonging to the protective layer approach saturation in the fluid.

Case of r 1 ≥ r 2 : once r 2 is significantly different from 0, r 2 approaches r 1 resulting in a roughly constant thickness x PL , then r 1 has no time dependence except due to pH variations.

Implementation in a reactive transport code

The CHESS/HYTEC computational code (from versions 3.7 to 5.2) developed by MINES ParisTech was chosen for implementing the GRAAL model. One interesting feature is that the GRAAL code does not require any software improvement: the input file that describes the experimental conditions is to be written taking into account GRAAL equations. How the model equations translate into code equations is described in Table 1 and with the help of the thickness–concentration relationship and the mass balance relationship described in the following paragraphs.

Table 1 From GRAAL model to reactive transport code Full size table

The requirements for the code to handle GRAAL equations, are the availability of the Monod type equation (cf. Table 1), the ability to use the concentration of a solid as a variable of the equation, and an adaptive timescale (available in CHESS/HYTEC code since version 3.7): the time step has to be small at the beginning of the calculation when the protective layer is thin in order to prevent its dissolution within a single time step.

A spherical shrinking-core model is used to describe the reduction in surface area due to the reduction in the glass grain size, although it does not have any effect in the experimental conditions discussed in this paper.

C PL (mol m−3) is defined as the concentration of the protective layer within a calculation cell. The thickness of the protective layer x PL is proportional to C PL (Eq. 12).

$$x_{\mathrm{{PL}}} = \frac{{C_{{\mathrm{PL}}}}}{{\rho _{{\mathrm{PL}}}S{\mathrm{/}}V}},$$ (12)

where, ρ PL is the protective layer’s density in mol m−3 of the protective layer, S is the protective layer’s surface area in the modeling cell (m2 m−3); the layer being thin and protective, is by definition equal to the surface area of the primary solid. V is the volume of the water in the calculation cell (m3 m−3).

The surface area to solution-volume ratio can also be written as the product of the specific surface area SSA (m2/mol) and the concentration C PS (mol m−3) of the primary solid (Eq. 13).

$$x_{\mathrm{{PL}}} = \frac{{C_{\mathrm{{PL}}}}}{{\rho _{{\mathrm{PL}}}{\mathrm{SSA}}\,C_{{\mathrm{PS}}}}}.$$ (13)

A mass balance can be written for each element belonging to the primary solid and the protective layer (Fig. 3 and Eq. 14).

$$F_i = \frac{{i_{{\mathrm{PL}}}\rho _{{\mathrm{PL}}}\;x_{{\mathrm{PL}}}}}{{i_{{\mathrm{PS}}}\rho _{{\mathrm{PS}}}\;x_{\mathrm{{Px}}}}},$$ (14)

where, F i is the fraction of element i coming from the primary solid and recovered in the protective layer, ρ PS is the primary solid’s density in mol m−3, i PS is the molar fraction of element I in the primary solid, i PL is the molar fraction of element i in the protective layer, x PS is the thickness of primary solid that has been altered, and x PL is the thickness of protective layer (Fig. 3).

Fig. 3 Mass and volume balance in the alteration layer. PS primary solid, PL protective layer, NPL non-protective layer Full size image

Measuring the composition and the structure of the protective layer enables the calculation of i PL and ρ PL . F i can be calculated by applying the condition of mass balance in the experiment between the fluid, the non-protective layer, and the protective layer. Then, the ratio of the thickness of the protective layer to the thickness of the primary solid being altered can be deduced from Eq. 14.

Identifying the two equations of Table 1 and using equations Eqs. 13 and 14 enables the code’s parameters k m (Eq. 15) and K in (Eq. 16) to be calculated as a function of the model’s parameters.

$$k_{\mathrm{m}} = r_{\mathrm{h}}.$$ (15)

$$K_{\mathrm{{in}}} = \frac{{D\frac{\pi }{2}}}{{r_{\mathrm{h}}}}F_i\frac{{i_{{\mathrm{PS}}}}}{{i_{{\mathrm{PL}}}}}\rho _{\mathrm{{PS}}}{\mathrm{SSA}}\,C_{{\mathrm{PS}}}.$$ (16)

If reaction 2 is a classical dissolution reaction, then reaction 1 occurring at the interface between the primary solid and the protective layer is a reaction that converts one solid into another plus dissolved elements. Therefore, two simple reactions occurring at the same rate are required to produce reaction 1: a congruent dissolution of the primary solid and a backward precipitation of the protective layer.

Composition of the amorphous layer

Prediction of the amorphous layer’s solubility is a priority for understanding the long-term behavior of glass as that determines the stability of the glass. Solubility is measured by fluid analysis under steady-state conditions; “apparent steady-state” refers here to the concentrations required for a reduction of at least three orders of magnitude in the alteration rate.

The former approach using separate simple end-members: SiO 2 , Ca(OH) 2 , AlO(OH), and ZrO 2 21 had the major limitation of not including the chemical interactions between the elements. A richer and more complex model requires new experimental data. Rajmohan et al.57 enumerated a variety of phenomena directly derived from experimental data from simple glasses containing Si, B, Na, Al, Ca, and Zr. These phenomena are useful for improving the choice of end-members stoichiometry and defining the minimum number of end-members required to model those experiments (Table 2).

Table 2 Influence of reaction phenomena on the choice of end-members Full size table

The composition of the amorphous layer is not only a function of the composition of the glass but also a function of the fluid composition, volume, and renewal rate. Therefore, it is meaningful to build a model applicable to any fluid or glass composition. Consequently, a single model is proposed for the six glasses studied by Rajmohan et al.,57 which is then applied to the thirty oxide glasses and the SON68 glass (Table 3) in order to verify whether it is possible to neglect the influence of minor elements.

Table 3 Glass compositions (oxide wt.%/mol%) Full size table

Figure 4 represents the composition of the end-members used to describe the amorphous layer. The polyhedron defined by the end-members must, therefore, include the glass compositions being modeled. In confined media, the Si/Al/Ca/Zr stoichiometry of the amorphous layer tends toward the stoichiometry of the initial glass. This model is not a solid solution model but simply a concatenation of end-members.

Fig. 4 Graphical representation of end-members stoichiometries Full size image

The proposed model is intended to be applied to experimental data in order to check the utility of the end-members to account for glass compositions and fluid pH variations of all available experiments and then to infer the missing parameters, mainly the logK values and diffusion coefficients of the end-members.

Parameters determination

This section presents all the parameters required in the input file of the model, starting with the experimental conditions, the glass composition, or the kinetic parameters already available in the literature. It also describes the method for determining new parameters values from the modeled experiments.

Six glass compositions were studied. The procedures for glass synthesis and preparation are described in ref. 57. Experiments were carried out in a static system at 90 °C with a glass surface area to solution-volume (S/V) ratio of 20 cm−1 (except CJ9: 48 cm−1) and with various imposed pH values (7, 8, 10) at 90 °C for ~150 days. For pH 7 and 8, 0.2 mol of tris (hydroxymethyl) aminomethane (TRIS) buffer solution was prepared. Potassium hydroxide solution was used as a buffer for pH 10 experiments. Experiments in initially pure water at free pH (near 9) and S/V = 80 cm−1 were initiated by Jegou58 and by Gin and Jegou59 and have been running since November 1996. These have been updated by Gin et al.38 Solutions were ultra-filtered to 10,000 Da before analysis. Twenty-four experiments are available for analysis: six glass compositions at four pH values (7, 8, ~9, 10).

The stoichiometry of the glasses is based on the composition of the French SON68 inactive reference glass, following the procedure initiated by Jégou, et al.58 In order to understand the effect of Ca, Al, and Zr, glasses were prepared by adding these oxides to the three major oxides of Si, B, and Na (Table 3). CJ4 is also known as the International Simple Glass ISG.60

r 2 values are measured in pure water based on Si and B almost congruent release. For those glasses in those experimental conditions, protective layer’s thickness is very small compared to altered glass thickness and non-protective layer’s contribution to Si release is negligible. Therefore, Si and B release are almost congruent.

They were measured on each simple glass by Jégou et al.58 at 90 °C, and pH 9.0; n and E a (Eq. 11) are assumed to have the same values as those measured for SON68 glass.61 This is justified here as the forward rates have a negligible effect under the experimental conditions discussed in this paper. The k + values are calculated from r diss pH 9, 90 °C to fit equation (Eq. 11) in order to have k + expressed in m s−1 (Table 4).

Table 4 Dissolution rate parameters determined for the glasses described in Table 3 Full size table

The interdiffusion coefficient varies with temperature and pH according to a simple relation (Eq. 17) proposed by Chave et al.56:

$$D = D_0\;\left[ {OH^ - } \right]^{n\prime }\;e^{ - \frac{{E\prime a}}{{RT}}}$$ (17)

, where n′ is the pH-dependence factor of the interdiffusion coefficient (dimensionless), E a ′ is the activation energy associated with the interdiffusion coefficient (kJ mol−1), and D 0 is the interdiffusion constant (m2 s−1).

In Eq. 17, hydroxide ions are preferred to hydronium ions because the equation is applied in neutral to basic pH conditions. Choosing hydronium ions, the minor species, would introduce the water dissociation constant into the equation, which is also a function of temperature and would change the numerical value of the activation energy, accordingly.

Experimental measurement of the diffusion coefficient is based on the diffusive flux arising in solution once dissolution of the protective layer by its outer surface has stopped. This flux is proportional to D/x PL (Eq. 9). The thickness of the protective layer, x PL , depends on the alteration thickness of the glass, x PS , and the molar density of the protective layer (Eq. 14). Therefore, knowledge of the molar density of the protective layer is required to calculate D. However, only the ratio D/x PL is used in the calculations. As a consequence, the choice is made here to select the molar density of the protective layer so that x PL = x PS (Eq. 14). Collapse of the gel, as identified by small angle X-ray scattering62 and modeled by a Monte Carlo method,54 is very limited in silica-saturated media. Isovolumetric alteration enables D to be calculated by simply applying Eq. 8 to the experimental data. Consequently, the diffusion coefficients given in Table 5 are calculated with this hypothesis of an isovolumetric alteration. Should the protective layer occupy 90% of the altered glass, then the diffusion coefficient would be 0.9 of the value given in Table 5.

Table 5 Diffusion rate parameters: values of the interdiffusion constant, D 0, and pH-dependence factor, n’, for the test glasses listed in Table 3 Full size table

The measurement of the diffusion coefficient can be performed at high surface area to solution-volume ratio where the dissolution thickness due to r 2 can be neglected compared to the diffusion thickness due to r 1 . At low surface area to solution-volume ratio, the formation of non-protective end-members will consume Si, delay saturation of PRI, and sustain glass dissolution, globally increasing the mass fraction of non-protective layer (Fig. 2).

The activation energies of the diffusion coefficients of these glasses are still unknown, though probably close to that of SON68 glass. D 0 values (Table 5) are calculated from 90 °C data57 considering the same activation energy as SON68 glass.56

Uncertainties of the diffusion coefficient values come both from the precision of experimental data and the effectiveness of Eq. 17 in taking pH and temperature variations into account. The value of the diffusion coefficient from Eq. 17 was applied using the parameters of Table 5 to calculate D within a factor of two, i.e., with a precision on the quantities of altered glass equal to a factor of √2. Moreover, under these experimental conditions, a value of 2 × 10−21 m2 s−1 is near the measurement limit and should be considered as a maximum value and not as an accurate measurement of the diffusion coefficient. The measurements determined from experiments at pH 10, all of which are below 2 × 10−21 m2 s−1, were therefore excluded when calculating the regression to determine the parameters indicated in Table 5. The calculated values extrapolated to pH 10 are consistent with the experimental data obtained at pH 10.

Figure 5 to Fig. 6 give representative examples of model results for the case of ISG glass, the glass that contains all the elements. The fit-quality is discussed in the next section.

Fig. 5 Total concentrations in solution at pH 7: comparison between model (blue line) and experiment (red circles) for ISG: pH 7, 90 °C, S/V = 20 cm−1. Real average pH of the experiment (7.2) is slightly above the set-point value (7.0). Al concentrations below-detection limits are not drawn Full size image

Fig. 6 Total concentrations in solution at pH 9: comparison between model and experiment for ISG: unconstrained pH, 90 °C, S/V = 80 cm−1 Full size image

Given the stoichiometry of the elements in the glass and in the end-members, the operation of the model can be summarized as follows:

The solubility products of the SiAlCa and SiAlNa end-members have a first-order effect on the Al activities in soda-lime glass and in soda glass, respectively. Al concentrations are lower in presence of Ca.

The solubility products of the SiZrCa and SiZrNa end-members have a first-order effect on the Zr activities in soda-lime glass and in soda glass, respectively. Given the below-detection limits concentrations of Zr in solution, the log K values indicated for these two end-members were chosen in order for the solution to remain undersaturated, with respect to zirconium oxide.

The Al/Ca and Zr/Ca stoichiometry determines the calcium fraction retained in the amorphous layer at pH ≤8. At pH <8, the SiCa end-member does not form. The solubility of the SiCa end-member determines the calcium activities for more alkaline pH values (e.g., pH 90 °C = 9).

The (Al-2Ca)/Na and (Zr-Ca)/Na stoichiometry determine the sodium fraction retained in the amorphous layer in initially pure water.

The SiAl and Si end-members have a first-order effect on the silicon activities.

The calculations are compared with the results obtained in 18 experiments, i.e., for all six glasses at three different pH values. The pH and elemental concentrations are measured in solution for Si, B, Na, Al, and Ca. The model-experiment comparison therefore covers a hundred curves of concentrations versus time, and a thousand experimental data points. This comparison is used to select the logK values corresponding to each of the end-members indicated in Table 6.

Table 6 End-members: compositions and logK values (90 °C) Full size table

Table 6 indicates the dissolution equations selected for the end-members, as well as the logK 90 °C value for each end-member. These values are determined by applying the model to experiments at pH 7, 8 and pH free57 and using CTDP database.23

The uncertainty presented in Table 6 gives an overview of the ability of the model to reproduce experimental data (discussed in the Results section). It is calculated, so that each available data point (Si, Al, and Ca concentrations), lies within the modeling results whatever the glass, pH, or sampling time.

Geochemical calculation code databases contain many mineral phases, but not all of them are likely to form at temperatures below 100 °C. Some can only be synthesized at high temperatures and cannot be observed in these experimental conditions. They are largely supersaturated in the leaching solutions and taking them into account without any kinetic limitation in the calculations would result in extremely low element concentrations in solution, far below the measured concentrations. Among the precipitating phases, it is important to discriminate between those with fast precipitation kinetics—for which the supply of a stoichiometric element from the glass or from the environment is a growth-limiting factor—and those with slow precipitation kinetics.53,63

Precipitation has been allowed for simple oxides and hydroxides whose precipitation at low temperatures has been established: aluminum hydroxide, portlandite, zirconium oxide, borax, and colemanite. However, all these minerals are undersaturated with respect to fluid experimental compositions and when modeling, precipitation of the amorphous layer’s end-members prevent their formation. Even if they do not form and, as a consequence, have no influence on the modeling results, they do deserve to be taken into account and cannot be removed from the database.

The secondary minerals selected here are those suitable for describing the laboratory experiments performed in initially pure water under oxidizing conditions. In a complex chemical environment, or under the reducing conditions expected in a geological repository, other minerals may have to be considered.

Hydrated calcium silicates and zeolites have been observed in the experiments at pH 10, at least after 150 days. Therefore, experiments at pH 10 could not be used to measure the end-members log K (cf.Table 6). Modeling of these experiments requires the kinetics of zeolite precipitation to be taken into account,34 which is beyond the scope of the paper.

Model application

In this section, the model’s usefulness and ability to account for experimental data at short- and long timescales is presented. An example of the calculation of the amorphous layer composition is given.

Below, we discuss the precision of the fit of the model to the data obtained from the experiments that were used for measuring the interdiffusion coefficient and the end-member solubilities. An examination of the results for ISG glass (Fig. 5 through Fig. 7) prompts the following remarks. These can be generalized to all the glasses and pH below 10:

1. The model accounts precisely for the pH variations with the glass composition (refer to the experiments with unconstrained pH: Fig. 6). The predicted values of pH lie within the uncertainty margin of pH measurement (±0.1 unit). Under these experimental conditions (closed system and high S/V ratio), the concentrations of B, Na (at all pH values), and Ca (at pH 7 and 8) depend mainly on the value of the diffusion coefficient. Variation over time is generally satisfactory, which supports the relevance of a time-squared model. The fit between the model and experimental results is related to the precision of the law describing the variations in the diffusion coefficient versus the pH (the uncertainties are discussed in the dedicated section). For example, at pH 8, for ISG glass, the D = f (pH) law is below the experimental value (Fig. 8); therefore the modeled B concentrations are underestimated for the experiment at pH 8 (Fig. 9). At pH 7 (Fig. 5) and pH 9 (Fig. 7), however, the fit is good. B and Na concentrations are usually known within an uncertainty margin of a factor of 1.4. 2. Si, Al, Zr, and pH 9-Ca concentrations depend on the amorphous layer model and on the pH. Uncertainties given for logK values in Table 6 are given based on the precision with which the concentrations are predicted by the model: the ±0.2 uncertainty in SiAl logK is equivalent to an uncertainty of a factor 1.5 in Si concentration; the ±1 uncertainty in SiAlNa logK is equivalent to an uncertainty of a factor 5 in Al concentration. This also reflects the experimental uncertainty of Al measurement, especially in alkali rich solutions obtained at pH 7.

Fig. 7 Calculated amorphous layer composition for ISG glass: unconstrained pH, 90 °C, S/V = 80 cm−1. a End-members concentration, (b) composition of the amorphous layer, oxygen excepted, in molar percentage Full size image

Fig. 8 pH dependence of interdiffusion coefficient: Interdiffusion coefficient calculated from boron equivalent thicknesses at 90 °C versus pH and glass composition; comparison with curves modeled by equation (Eq. 8) and parameters of Table 5 Full size image

Fig. 9 Total concentrations in solution at pH 8: Comparison between model (blue line) and experiment (red circles) for ISG: at pH 8, 90 °C, S/V = 20 cm−1 Full size image

Figure 7 shows the time-dependent composition of the amorphous layer. A single experiment is shown here as an example, at unconstrained pH. Due to the presence of Ca in sufficient amounts, the SiZrCa and SiAlCa end-members form at the expense of the soda end-members, SiZrNa and SiAlNa. Similarly, the presence of Al in sufficient amounts leads to the formation of the SiAl end-member and the absence of the Si end-member. Finally, as Ca is in excess with respect to the total Al/2 + Zr, and the pH is high enough, the SiCa end-member can form. The relative amounts of each end-member vary over time and thus modify the average composition of the layer (Fig. 7). Si, the most soluble of the elements that form the amorphous layer, accounts for an increasing proportion of the amorphous layer over time, tending toward the initial glass composition ratio.

Experiments at unconstrained pH were initiated several years ago to demonstrate the evolution of the reactive diffusion coefficient over time. Over the long-term, the model accurately accounts for the B concentrations for some glass compositions, but is conservative for others (e.g., ISG, Fig. 10). Experimental data55,64 and Monte Carlo models of glasses17,54,62,65 show that it is not sufficient to reach saturation with respect to Si to observe an inhibition of alteration: the gel must also reorganize. A decrease in the diffusion coefficient over time could therefore be due to a reorganization of the elements within the amorphous layer. Another hypothesis is that more Ca is incorporated in the amorphous layer,12 consistent with the slightly increasing pH. Indeed, an increase in pH causes a strong decrease in the diffusion coefficient for calcium-containing glass (Fig. 8).