Thermodynamic Properties of Aqueous Species Calculated Using the HKF Model: How Do Different Thermodynamic and Electrostatic Models for Solvent Water Affect Calculated Aqueous Properties?

Miron, George D.; Leal, Allan M. M.; Yapparova, Alina

doi:https://doi.org/10.1155/2019/5750390

Geofluids

On this page

Abstract Introduction Methods Results Discussion Conclusions Data Availability Conflicts of Interest Acknowledgments Supplementary Materials References Copyright Related Articles

Special Issue

Advances in Numerical Simulations of Hydrothermal Ore Forming Processes

View this Special Issue

Research Article | Open Access

Volume 2019 | Article ID 5750390 | https://doi.org/10.1155/2019/5750390

Thermodynamic Properties of Aqueous Species Calculated Using the HKF Model: How Do Different Thermodynamic and Electrostatic Models for Solvent Water Affect Calculated Aqueous Properties?

George D. Miron,¹Allan M. M. Leal,²and Alina Yapparova³

Academic Editor: Yuan Mei

Received27 May 2018

Accepted13 Sept 2018

Published07 Feb 2019

Abstract

Thermodynamic properties of aqueous species are essential for modeling of fluid-rock interaction processes. The Helgeson-Kirkham-Flowers (HKF) model is widely used for calculating standard state thermodynamic properties of ions and complexes over a wide range of temperatures and pressures. To do this, the HKF model requires thermodynamic and electrostatic models of water solvent. In this study, we investigate and quantify the impact of choosing different models for calculating water solvent volumetric and dielectric properties, on the properties of aqueous species calculated using the HKF model. We identify temperature and pressure conditions at which the choice of different models can have a considerable effect on the properties of aqueous species and on fluid mineral equilibrium calculations. The investigated temperature and pressure intervals are 25–1000°C and 1–5 kbar, representative of upper to middle crustal levels, and of interest for modeling ore-forming processes. The thermodynamic and electrostatic models for water solvent considered are: Haar, Gallagher and Kell (1984), Wagner and Pruß (2002), and Zhang and Duan (2005), to calculate water volumetric properties, and Johnson and Norton (1991), Fernandez and others (1997), and Sverjensky and others (2014), to calculate water dielectric properties. We observe only small discrepancies in the calculated standard partial molal properties of aqueous species resulting from using different water thermodynamic models. However, large differences in the properties of charged species can be observed at higher temperatures (above 500°C) as a result of using different electrostatic models. Depending on the aqueous speciation and the reactions that control the chemical composition, the observed differences can vary. The discrepancy between various electrostatic models is attributed to the scarcity of experimental data at high temperatures. These discrepancies restrict the reliability of the geochemical modeling of hydrothermal and ore formation processes, and the retrieval of thermodynamic parameters from experimental data at elevated temperatures and pressures.

1. Introduction

Aqueous fluids play a major role in geochemical processes that occur in the Earth’s crust [1–9]. Field observations as well as experimental and theoretical studies have shown that these fluids contain a variety of ions and aqueous complexes that are responsible for mineral dissolution-precipitation and the mobility of elements [7, 10–15]. Fluid-rock interactions are ubiquitously present in ore deposition, metasomatic, and hydrothermal processes (e.g., [16–19]).

At elevated temperature and pressure conditions, at which direct observations are not possible and experiments are difficult to perform, thermodynamic modeling is key to gain a better understanding of the chemistry of aqueous fluids involved in fluid-rock interactions and for quantifying elemental mobility (e.g., [11, 20–24]). The reliability of geochemical and reactive transport modeling, however, among other things, depends on the accuracy of calculated thermodynamic properties of aqueous species (ions and complexes) at the temperature and pressure conditions of interest [25–30].

The most widely used model for calculating thermodynamic properties of aqueous ions and complexes in geochemical applications is the Helgeson-Kirkham-Flowers (HKF) equation of state. The HKF model has been originally formulated for the calculation of standard partial molal thermodynamic properties of aqueous species up to 1000°C and 5 kbar, and water density [31–41]. These conventional limitations are related to the accuracy of the model and the accuracy of the calculated water density and dielectric constant. Sverjensky et al. [9] developed an empirical model for calculating the dielectric constant of water allowing the HKF model to be applied at deep crustal-level conditions with pressures up to 60 kbar and temperatures up to 1200°C [4, 9, 42–45].

Although a versatile model, the HKF model is not accurate enough when applied to thermodynamic conditions close to the critical point of water, for calculating properties of neutral species, and in low-density solutions [46–48]. Over the past years, however, a variety of thermodynamic models (or equations of state) for calculating the properties of aqueous solvent and solutes have been formulated with superior calculation accuracy at those conditions for which the HKF model is not sufficiently reliable. A review of these models is presented in Dolejš [47]. Despite the recent efforts in the formulation of more accurate models, the HKF model remains dominant in the geochemical modeling community, and this will most likely continue in the near future given its widespread use in many geochemical software.

The success of the HKF model is notably due to the extensive parameter dataset ([40]; SUPCRT datasets, e.g., slop98.dat and slop16.dat) that covers a large part of the periodic table and due to the broad set of correlations that were developed to estimate standard state properties and equation of state coefficients for systems with little or no experimental data [38, 49]. The standard state thermodynamic properties from the SUPCRT dataset together with the HKF model are being used for geochemical modeling in several Gibbs energy minimization (GEM) code packages such as GEM-Selektor [50], Reaktoro [51], Perplex [52], and ChemApp [53]. Furthermore, the SUPCRT dataset, together with the HKF model, was applied for generating the bulk of the existing law of mass action (LMA) thermodynamic databases used in popular geochemical modeling codes such as PHREEQC [54], Geochemist’s Workbench [55], EQ3/6 [56], SOLVEQ, and CHILLER [57]. The SUPCRT dataset is distributed with the HKF model [38, 39] and is still being further extended [9, 43, 58–60].

The HKF model considers the standard partial molal properties of aqueous species as the sum of nonsolvation and solvation contributions. The nonsolvation contributions are calculated using semi-empirical HKF parameters. In order to calculate the solvation contributions, the HKF model uses the Born theory [61], which defines the solvation contributions as functions of the water dielectric constant. In turn, models that calculate the dielectric constant are functions of the water density. Consequently, in order to obtain the properties of aqueous species using the HKF model, one needs to choose an electrostatic model for calculating the dielectric constant of water and a thermodynamic model to compute the water density.

Both electrostatic and thermodynamic models rely on data obtained from experiments. Experimental data on the dielectric constant of water are only available up to 550°C and 5 kbar, and no new experiments at elevated temperatures and pressures have so far been reported [2, 9, 46, 62]. Several models have been developed to predict the dielectric constant of water at temperatures up to 1200°C and pressures of 60 kbar [9, 32, 63–65]. However, some discrepancies are known between the existing electrostatic models, especially at temperature and pressure conditions for which no experimental data are available [9, 47, 66].

In contrast to dielectric constant, there exists a significant amount of experimental data on water density at temperature and pressure conditions relevant to hydrothermal processes [67]. Thus, models for calculating thermodynamic properties of water are more reliable than those for calculation of electrostatic properties. The International Association for the Properties of Water and Steam (IAPWS) has recommended accurate thermodynamic models (equations of state) for pure water calibrated using reliable, high-quality experimental data [67]. Other equations of state for water have been developed for high-pressure conditions that extend the recommended applicability ranges of the IAPWS models [68, 69].

The focus of this study is to establish how the choice of thermodynamic and electrostatic models for water solvent can affect the calculated properties of aqueous species and to assess the discrepancies in the calculated properties between commonly used thermodynamic and electrostatic models in the framework of the Helgeson-Kirkham-Flowers model. The selected thermodynamic models are Haar, Gallagher, and Kell [70]—HGK84 (also known as IAPS84), Wagner and Pruß [67]—IAPWS95, and Zhang and Duan [68]—ZD05. The electrostatic models are Johnson and Norton [46]—JN91, Fernandez and others [63]—FE97, and Sverjensky and others [9]—SV14.

The selection made for this study is by no means exhaustive. These specific models were chosen because they are widely used in geochemical modeling of hydrothermal processes [4, 9, 40]. The HGK84 thermodynamic and JN91 electrostatic models are most frequently used with the HKF model being part of the SUPCRT92 software. The IAPWS95 thermodynamic and FE97 electrostatic models are the recommended standard from the International Association for the Properties of Water and Steam and can replace HGK84 and JN91. The ZD05 thermodynamic and SV14 electrostatic models are used for extrapolations of the HKF model above the previous 5 kbar limit, up to 60 kbar.

Our comparison study is restricted to temperatures up to 1000°C and pressures up to 5 kbar, and water density , which are within the conventional ranges of applicability of all considered models and are representative of upper to middle crustal conditions and thus relevant for hydrothermal and ore formation processes.

In the first part, we briefly present the chosen models and their applicability intervals (Table 1). In the following part, we compare the volumetric and dielectric properties of pure water solvent, calculated using different thermodynamic (HGK84, IAPWS95, and ZD05) and electrostatic (JN91, FE97, and SV14) models and their effect on the calculated thermodynamic properties of aqueous species and reactions (, , , , and ). Finally, we compare results from equilibrium calculations concerning mineral solubilities, mineral assemblages, and fluid-rock interaction. Some recommendations on the model applicability and interchangeability are provided in the conclusions. This study serves as practical information for many modelers that use the HKF model and makes them aware of the differences that may arise at certain temperatures and pressures related to conditions of ore deposition and hydrothermal processes.

2. Thermodynamic Models

2.1. Helgeson-Kirkham-Flowers Model

The initial development of the HKF model is described in detail in a series of papers [31–34]. It was later revised by Tanger and Helgeson [35] and Oelkers and Helgeson [71] to allow the calculation of standard partial molal properties of aqueous species at temperatures up to 1000°C and pressures up to 5 kbar. The model is also limited to , , or for bar. This conventional limitation is related to accuracy of the pressure- and temperature-dependent solvent function and its partial derivatives needed for calculating the properties of charged aqueous species [40, 72].

In the HKF model, the contributions to the molal properties are separated into nonsolvation () and solvation contributions ():

The nonsolvation contributions are determined using semi-empirical heat capacity and volume equations and the solvation contributions using the dielectric constant of water. At high temperatures and low pressures, the solvation contributions predominate over the nonsolvation contributions [35]. The HKF model uses an additivity relationship [73] between the standard partial molal thermodynamic properties of the electrolyte () with properties of the aqueous ions (): where represents the number of moles of the ion in one mole of the electrolyte (Tanger and Helgeson [35]; (1)). This permits the use of experimental data for aqueous electrolytes to infer the standard partial molal properties of their associated aqueous ions. Tanger and Helgeson [35] used experimental data for partial molal volumes, compressibilities, and isobaric heat capacities of electrolytes, to regress the parameters in the semi-empirical expressions for the nonsolvation contributions for electrolytes (parameters , , , , , and in the HKF model). From the resulting parameters for electrolytes, using the additivity principle and the convention that all properties of the H⁺ ion are zero at any temperature and pressure, Tanger and Helgeson [35] calculated the HKF parameters of several aqueous ions and complexes. These form the core HKF dataset from which different correlations were later developed and used to generate HKF parameters for many other ions and complexes [38, 39, 49]. In some cases, to improve the accuracy of calculations, the HKF parameters are optimized against experimental data [27, 29, 58, 74].

The solvation contributions are formulated using the Born theory [61] as functions of the water dielectric constant. The standard partial molal Gibbs-free energy of solvation for the ion is expressed by where represents the dielectric constant of water and is the conventional Born coefficient of the ion, which is expressed as (see Tanger and Helgeson [35], their eq. (13)) where represents an HKF model constant (1.66027∙10⁵ Å∙cal∙mol⁻¹), is the charge of the species, is a function of temperature and solvent density [72], and is the effective electrostatic radius of the species at reference temperature and pressure : where is the conventional Born coefficient at reference temperature and pressure.

Other partial molal solvation properties such as volume, , entropy, , and isobaric heat capacity, , are calculated as pressure and temperature derivatives of the free Gibbs energy of solvation using the formulae where , , and are the Born functions (see Tanger and Helgeson [35], their appendix and ):

The temperature and pressure derivatives of the dielectric water constant, which are necessary for the , , and evaluation, are obtained from [46] which relate the dielectric constant derivatives to the water density, the water isothermal compressibility: the isobaric expansivity and the isobaric temperature derivative of the isobaric expansivity

Equations (3), (4), (5), (6), (7), (8), (9), (10), (11), (12), (13), and (14) show the dependence of the HKF model on both the water density and the water dielectric constant.

2.2. Haar, Gallagher, and Kell Water Equation of State (HGK84)

The Haar, Gallagher, and Kell [70] water equation of state (also known as IAPS84) is represented as a Helmholtz free energy function, , of temperature, , and density, . This equation of state is assumed to be accurate up to 1000°C and 15 kbar, except within a region near the critical point, where it is not accurate enough for calculating derivative properties like the isothermal compressibility, , and the isobaric expansivity, [46, 67].

The equation of state outside the critical region can be written as where represents a conversion factor and is the molecular weight of water. Other thermodynamic properties can be calculated by combining the derivatives of the Helmholtz free energy function ([67]; their Table 6.3). This equation of state is used in SUPCRT92 [40] to calculate the thermodynamic properties of water, except for the region close to the critical point where the model described by Sengers et al. [75] is used.

2.3. Wagner and Pruß Water Equation of State (IAPWS95)

The IAPWS95 [67] formulation was adopted as the scientific standard for calculation of the thermodynamic properties of ordinary water substance for general and scientific use, by the International Association for the Properties of Water and Steam (IAPWS). The IAPWS95 equation of state showed significant improvement in accuracy of the prediction of thermodynamic properties of water close to the critical point of water, in metastable low temperature fields and for extrapolations to very high temperature and pressure. The equation of state for water is expressed as a relation for the Helmholtz free energy of water as a function of temperature and water density. The model was parameterized by simultaneously fitting measurements of different properties of water over a wide range of conditions (multiproperty fitting). This formulation is recommended for calculating water fluid properties up to 1000°C and 10 kbar and has a smooth extrapolation to 100 kbar.

2.4. Zhang and Duan Water Equation of State (ZD05)

Zhang and Duan [68] developed an empirical pressure-volume-temperature equation of state for water with 14 parameters. The parameters were optimized against selected experimental data covering a wide range of temperatures and pressures, and data from molecular dynamics and Monte Carlo simulations. Compared to the IAPWS95, high-pressure density data was used in the regression analysis, beyond the 10 kbar recommended limit of calculations of the IAPWS95. Below 10 kbar pressure, a similar data selection as for IAPWS95 was used. The authors suggested that the equation can be applied for conditions between 0 and 2000°C and between 1 kbar and the water ionization limit (200–300 kbar).

3. Electrostatic Models

3.1. Johnson and Norton Water Electrostatic Model (JN91)

Johnson and Norton [46] developed an empirical equation for calculating the dielectric constant of water in a recommended range of conditions from 0 to 1000°C and 1 to 5000 bar as a function of temperature, pressure, and water density. They modified the equation proposed by Uematsu and Frank [64] to accurately represent the values from the model of Helgeson and Kirkham [32] for temperatures ≤350°C and the semi-empirical equation developed by Pitzer [76] for predictions at temperatures from 550°C to 1000°C. The formulation developed by Johnson and Norton [46] has a good representation of the experimental data, the approximations of Helgeson and Kirkham [32] equations, and the approximations and predictions of the Pitzer [76] equation [46]. When evaluating the experimental and predicted data for the JN91 model refinement, Johnson and Norton [46] calculated the density of water and its derivatives using the HGK84 model.

3.2. Fernandez et al. Water Electrostatic Model (FE97)

The dielectric constant model proposed by Fernández et al. [63] uses the Kirkwood equation. Fernández et al. [63] reviewed the existing data on the dielectric constant of water and used it for regressing a 12-parameter equation. In the process of evaluating and regressing the experimental data for the dielectric constant, the water density from the dielectric constant formulation of Fernández et al. [63] was calculated at the temperature and pressure of the experiments using the IAPWS95 equation of state. The recommended range of applicability of this model is from −35°C to 600°C and up to a pressure of 12 kbar. The International Association for the Properties of Water and Steam recommends the use of the electrostatic model of Fernández et al. [63] combined with the IAPWS95 thermodynamic model.

3.3. Sverjensky et al. Water Electrostatic Model (SV14)

In order to extend the application of the HKF model for calculating properties of aqueous species to upper mantle conditions, Sverjensky et al. [9] used the new ab initio molecular dynamics simulation data and provided the first empirical relation for the dielectric constant of water up to 60 kbar and 1200°C. They first used the equation derived by Franck et al. [65] and recalibrated it using the experimental data of Heger et al. [77] up to 550°C, data from ab initio molecular dynamics at 727°C from Pan et al. [78], and estimates from the formulation of Fernández et al. [63] up to 927°C. Sverjensky et al. [9] suggest that the equation of Franck et al. [65] can be applied up to a density of 1100 kg/m³. Based on the previous observations that the dielectric constant of water is linearly dependent on the natural logarithm of the water density [77], Sverjensky et al. [9] formulated an empirical equation for as a function of at higher densities. For regressing the parameters in this equation, values for calculated with the recalibrated equation of Franck et al. [65] were used. They suggested that the linear relation holds up to high pressures (60 kbar) by comparing their results to the results of ab initio molecular dynamic simulations [78]. The density of water was calculated from the equation of Zhang and Duan [68]. The authors recommended the use of the model for temperatures of 100–1200°C and pressures of 1–60 kbar.

4. Setup and Methods

We investigate how different thermodynamic and electrostatic models (Table 1) affect the calculation of standard molal properties of aqueous species using the HKF model. As underlined before, these models are used to calculate the solvation contributions to the standard molal properties. In the following comparison, we focus on differences in the calculated standard state properties that arise from the solvation contributions ((3), (6), (7), and (8)) calculated using different solvent property models. No differences result from the nonsolvation contributions. Unless stated otherwise, we use the same species-specific set of HKF parameters (, , , , , and ). The temperature, pressure, and density conditions we considered were 25–1000°C, 1–5 kbar, and , respectively (data points below the 350 kg/m³ HKF density limit were excluded from the plots). The calculated solvent properties refer to that of pure water. No solute effects on the thermodynamic and electrostatic properties of the solvent are considered. The thermodynamic properties of solutes refer to standard partial molal properties at infinite dilution and the apparent Gibbs energy in the Benson-Helgeson convention [79, 80]. In order to model real systems, an activity coefficient model that corrects for the nonideal mixing contribution (excess properties of the solution) is necessary. In the current study, when calculating fluid-mineral equilibria, the activity coefficients () of aqueous species are calculated by the extended Debye-Hückel equation [34], which is often used together with the HKF model where and represent the Debye-Hückel solvent parameters, the ion size parameter, a semi-empirical extended term parameter, and the effective ionic strength (corrected for ion pairing and complexing). For the aqueous fluid-mineral equilibrium calculations, the parameters suitable for a NaCl major background electrolyte were used with and at 25°C, both with temperature and pressure dependence [34, 81]. This equation is applicable up to a 3-molal background electrolyte concentration.

The volumetric and dielectric properties of water solvent and the standard thermodynamic properties of ions, complexes, and reactions were computed using the ThermoFun library (Miron and others, in preparation, https://bitbucket.org/gems4/thermofun). The library calculates standard thermodynamic properties of substances and reactions at the temperature and pressure of interest, using an extensive collection of equations of state and empirical models for substances and reactions. The Reaktoro ([51]; http://reaktoro.org, [82]) geochemical modeling tool was used for the aqueous fluid-mineral equilibrium calculations. Reaktoro is an open-source framework in C++ and Python that uses state-of-the-art numerical methods and that can be used for simulating chemical systems governed by equilibrium, kinetics, or a combination of both.

5. Results

In the following examples, we performed calculations with different combinations of thermodynamic and electrostatic models together with the HKF model and looked at the differences in the calculated properties compared to the baseline, where we used the IAPWS95 thermodynamic model and the FE97 electrostatic model, which currently are the recommended standard.

5.1. Comparison between Calculated Water Solvent Properties

In this section, we compare values of water properties (density, isobaric expansivity, isothermal compressibility, isobaric temperature derivative of isobaric expansivity, dielectric constant, and values of Born coefficients) calculated using different combinations of thermodynamic and electrostatic models.

First, we assess the effect of using different thermodynamic models on the volumetric properties of water: density, isobaric expansivity, isothermal compressibility, and isobaric temperature derivative of isobaric expansivity.

The IAPWS95, HGK84, and ZD05 thermodynamic models produce consistent water density values, , in most of the considered temperature and pressure interval (25 to 1000°C, 1 to 5 kbar), with deviations of HGK84 and ZD05 from IAPWS95 predominantly less than 1%. However, the discrepancies increase with increasing temperature and pressure, as seen in Figure 1(a). The largest discrepancies between IAPWS95 and ZD05 can be seen at 1 kbar above 550°C, with deviations abruptly increasing to 2.5% at 600°C (calculated densities 374 kg/m³ and 365 kg/m³ for IAPWS95 and ZD05). The biggest mismatch between IAPWS95 and HGK84 occurs at and , where the deviations increase to 1.5%.

(a)

(b)

(c)

(d)

The values of isobaric expansivity, , and isothermal compressibility, , calculated using the three different models (Figures 1(b) and 1(c)) are overall in good agreement between 50 to 1000°C and 1 to 5 kbar. For , the deviations of HGK84 and ZD05 from the IAPWS95 model do not go beyond 5%, with the exceptions at (with increasing pressure deviations go up to 18%) and at , 1 kbar for ZD05 (deviations go up to 11%). For , the deviations remain below 5.2%, except at 1 kbar, , for ZD05 compared to IAPWS95.

For the isobaric temperature derivative of the isobaric expansivity, , the values obtained using HGK84 and ZD05 have discrepancies from the IAPWS95 model smaller than 25% only at conditions between 100 and 600°C (Figure 1(d)). Below 100°C () and in the region between 600 and 800°C, the differences go up to 140%. At these temperature and pressure conditions, the isobaric expansivity changes from very close to zero positive to negative values, resulting in large percent deviations for very small values. Overall, the ZD05 model always shows increasing discrepancies for all calculated water volumetric properties at 1 kbar pressure with increasing temperature.

Secondly, we look at the effect that different combinations of thermodynamic and electrostatic models have on the values of the dielectric constant of water. The electrostatic models used to calculate the water dielectric properties are functions of the water density and its derivatives. Unlike the volumetric properties of water, which are only affected by the choice of the thermodynamic model, the dielectric constant of water () is influenced both by the choice of the electrostatic and thermodynamic models.

Exchanging the thermodynamic models produce in general small differences. The use of the FE97 electrostatic model (recommended standard by the IAPWS) for the dielectric constant of water, combined with the HGK84, and ZD05 thermodynamic models, results in deviations from the IAPWS95 model not larger than 2.1% (except for ZD05 at 1 kbar), as seen in Figure 2(a). As expected ( being a function of ), the percent deviations are similar to the ones calculated for water density (Figure 1(a)), with a slightly larger spread overall.

(a)

(b)

Large discrepancies result when exchanging the electrostatic models. The differences of the JN91 and SV14 electrostatic models from FE97 for the calculated dielectric constant, , using the IAPWS95 thermodynamic model to calculate increase with increasing temperature above 500°C (Figure 2(b)). Below 500°C, the discrepancies stay below 5%, while above 500°C, they go up to 19%. Overall, the differences decrease with increasing pressure, with SV14 showing a 9.5% difference and JN91 a 2% difference from FE97 at 1000°C and 5 kbar.

Thirdly, we look at the differences in Born functions, calculated using different electrostatic models. , , and Born coefficients all depend on (see (9), (10), and (11)), but as observed in Figure 2(a), the dependence of on the choice of a thermodynamic model is weak; therefore, we only compare the impact of using different electrostatic models for calculating Born coefficients. We use FE97 as baseline and the IAPWS95 thermodynamic model to calculate the water density and its derivatives.

The differences for the Born functions (, , and ) calculated using the JN91 and SV14 models from the FE97 model show a similar pattern as in the case of (Figure 2(b)). For the Born function (Figure 3(a)), the discrepancies from FE97 are under 10% below 500°C with the exception of SV14 at temperatures lower than 100°C (which is the recommended lower limit for the SV14 model, Sverjensky et al. [9]). Above 500°C, the discrepancies increase with temperature but decrease with pressure, with a maximum difference of 35% at 950°C and 2 kbar. For the Born function (Figure 3(b)), the maximum disagreement between JN91 and FE97 is 48% and that between SV14 and FE97 is 35%. Above 600°C, the differences start to increase with temperature and pressure between FE97 and JN91 and increase with temperature but decrease with pressure between FE97 and SV14. The Born function (Figure 3(c)) shows increasing discrepancies below 200°C for SV14, while the overall differences increase with temperature and decrease with pressure above 600°C for SV14 and JN91 compared to FE97.

(a)

(b)

(c)

Figures EA1, EA2, and EA3 in the electronic supplementary material are analogous to Figures 1–3, but show the absolute differences in property values instead.

5.2. Comparison of the Calculated Standard Properties of the Al³⁺ and OH⁻ Ions and of the Muscovite Forming and Water Dissociation Reactions

The use of different combinations of thermodynamic and electrostatic models will affect the values of standard state properties of aqueous ions and complexes, as well as equilibrium constants of reactions, obtained from the HKF model. Discrepancies in the calculated values of , , , and propagate in the calculated values of , , , and . The dependence of the standard state properties on the water dielectric constant and its derivatives can be seen in (3), (6), (7), (8) (9), (10), and (11). Furthermore, the electrostatic models are dependent on the thermodynamic model used to calculate the water density and its derivatives ((12), (13), and (14)). The dependence is also visible when looking at similar patterns of the comparisons between the calculated properties of aqueous species, presented below, and those of water as seen in Figures 1–3 and Figures EA1, EA2, and EA3 from the electronic supplementary material.

In the first part of this section, we chose the Al³⁺ aqueous species to investigate the effect of different thermodynamic and electrostatic models on its partial molal standard state properties: Gibbs energy, entropy, heat capacity, and molar volume. Aluminum is a major rock-forming element commonly showing low solubility in aqueous solutions for a wide range of crustal conditions. The high charge and high value for the conventional Born coefficient of Al³⁺ amplify the effect of using different water properties models (see (4) and (5)). The results showing differences in the calculated thermodynamic properties of Al³⁺ are presented in Figure 4, comparing thermodynamic models, and Figure 5, comparing electrostatic models.

(a)

(b)

(c)

(d)

(a)

(b)

(c)

(d)

First, we compared how different thermodynamic models affect the standard molar properties of Al³⁺, using the FE97 electrostatic model for the calculation of the dielectric properties of water. Despite the large visual differences observed in Figure 4, the standard molal properties of Al³⁺ are in good agreement for all thermodynamic models, since the deviations are actually relatively small. For most of the investigated temperatures and pressures (25 to 1000°C, 1 to 5 kbar), except at 1 kbar and temperatures above 800°C, the differences in , , , and do not exceed 1.5 kJ·mol⁻¹, 7 J·mol⁻¹·K⁻¹, 20 J·mol⁻¹·K⁻¹, and 1 cm³, respectively. The largest deviation is observed between the IAPWS95 and ZD05 models at 1 kbar above 500°C and between IAPWS95 and HGK84 above 800°C. Differences in between IAPWS95 and HGK84 increase with temperature and pressure with a maximum of 3.2 kJ·mol⁻¹ at 1000°C and 5 kbar.

Second, we compared how different electrostatic models affect the standard molar properties of Al³⁺. The thermodynamic model IAPWS95 is used for the calculation of the volumetric properties of water. As Figure 5 shows, inconsistency between the selected electrostatic models (JN91, SV14, and FE97) for temperatures above 500°C strongly affects the calculated standard properties of Al³⁺. However, with increasing pressure at temperatures above 500°C, although large, the mismatch decreases. For the values, there is good agreement (differences <1.5 kJ·mol⁻¹) between the investigated models up to 500°C (Figure 5(a)). Above this temperature, the discrepancies between JN91 and SV14 from FE97 strongly increase with increasing temperature but decrease with increasing pressure with maximum differences of 58 and 41 kJ·mol⁻¹ between FE97 and JN91, and SV14 and FE97, respectively. For the calculated and , large differences between FE97 and SV14 (Figures 5(b) and 5(c)) can be seen below 100°C (the lower limit of the SV14 model), while the values calculated with the JN91 and FE97 models are in good agreement at these conditions. Above 400°C, the differences between FE97 and SV14, and between FE97 and JN91 for and , increase with increasing temperature, but decrease with increasing pressure. There is good agreement in the calculated (Figure 5(d)) with noticeable differences increasing above 600°C (except at 1 kbar above 500°C) but decreasing with pressure.

Figures EA4 and EA5 (electronic supplementary material) show differences in the standard state properties for the OH⁻ species, calculated using different combinations of the thermodynamic and electrostatic models. The observed trends are similar to those from Figures 4 and 5 for the Al³⁺ ion. The discrepancies between the models are less pronounced due to the lower absolute value of the charge of OH⁻ (−1 compared to +3).

In the second part of this section, we compare the equilibrium constant values of the muscovite-forming and water dissociation reactions

The muscovite forming reaction (17) is not characteristic for natural systems, but it is a complex reaction that includes several ions and it is used in many LMA thermodynamic databases to calculate geochemical equilibria (e.g., SOLTHERM-2011.XPT and SOLVEQ database). Therefore, it was chosen together with the water dissociation reaction to illustrate and quantify the effect of using different thermodynamic and electrostatic models on values of reaction equilibrium constants.

When comparing IAPWS95 to HGK84 and ZD05 formulations for the volumetric properties of water, the differences in the values do not exceed 0.5 log units for the muscovite reaction (Figure 6(a)), and 0.1 log units for the water reaction (Figure 6(b)). This is similar to the case of individual standard state properties of ions. When different formulations for the electrostatic properties of water are used, the discrepancies between FE97 compared to JN91 and SV14 remain below 0.5 log units for the muscovite reaction (Figure 6(c)) and below 0.1 log units for the water reaction (Figure 6(d)) between 25 and 400°C. The discrepancies start to increase sharply above 400°C, as it was also observed for the standard state properties of individual ions (Figure 5). The differences for the muscovite reaction exceed 7 log units, while the maximum difference for the water reaction is 1.45 log units at 950°C, 2 kbar, between calculations using FE97 and JN91.

(a)

(b)

(c)

(d)

Figure 6

Differences in values of muscovite precipitation and water dissociation reactions, obtained from different thermodynamic models: (a and b) full lines—difference between IAPWS95 and HGK84, and dashed lines—between IAPWS95 and ZD05 models for calculating water volumetric properties, using FE97 for water dielectric properties; (c and d) full lines—difference between FE97 and JN91, and dashed lines—between FE97 and SV14 models for calculating water dielectric properties, using IAPWS95 for water volumetric properties. A vertical dashed line at 500°C is shown on (c and d).

To see if we can determine the use of which electrostatic model produces values that are more accurate, we compare them with an independent equation for calculating the water dissociation constant. Bandura and Lvov [83] made an exhaustive literature review concerning the dissociation of water and proposed an empirical equation (BL06), with seven adjustable parameters, recommended to be valid for temperatures between 0–800°C and for densities between 0 and 1.2 g/cm³. We compare the BL06 model with the water dissociation constants calculated using the JN91, FE97, and SV14 electrostatic models (Figure 7). For computing the reaction constant ((18)), the Gibbs energy of water was calculated using the IAPWS95 thermodynamic model, while the Gibbs energy of OH⁻ was calculated using the HKF model. In addition to exchanging different electrostatic models (used for calculating the solvation contributions), one comparison was done using a different set of HKF parameters for OH⁻ (used for calculating the nonsolvation contributions). First, we used the original OH⁻ HKF parameters from Tanger and Helgeson [35] and tested each electrostatic model (SV14, FE97, and JN91; Figure 7 dashed, dotted, and dashed-dotted lines). Second, we used the OH⁻ HKF parameters reported in the DEW model (Deep Earth Water model, v. 11.0.1) dataset together with the SV14 electrostatic model (Figure 7, full lines). In the DEW dataset, it is stated that the HKF parameters for OH⁻ were optimized against values calculated with the BL06 model. When using the HKF parameters for OH⁻ from Tanger and Helgeson [35], at temperatures smaller than 600°C, differences from BL06 remain below 0.5 log units except at 1 kbar when using SV14 (Figure 7, dashed lines). If the HKF parameters for OH⁻ are taken from the DEW dataset and the SV14 electrostatic model is used, differences from BL06 remain below 0.2 log units (expect at 1 kbar, Figure 7, full lines). For all cases, above 600°C discrepancies increase with increasing temperature, up to 4 log units (between BL06 and SV14 at 950°C and 2 kbar), and decrease with increasing pressure.

Figure 7

Differences in calculated water dissociation constants from values obtained with the empirical model of Bandura and Lvov [83] (BL06) and values obtained from using SV14, FE97, and JN91 electrostatic models for calculating the standard state properties of OH⁻. The Gibbs energy of H₂O was calculated using the IAPWS95 thermodynamic model. Full lines—differences between BL06 and the SV14 electrostatic model using the nonsolvation HKF parameters reported by Sverjensky et al. [9]. The remaining differences are obtained using the nonsolvation HKF parameters reported by Tanger and Helgeson [35]. Dashed lines—differences between BL06 and SV14. Dotted lines—between BL06 and FE97. Dashed-dotted lines—between BL06 and JN91.

5.3. Effect of Using Different Electrostatic Models in the Modeling of Fluid Mineral Equilibria and Fluid-Rock Interaction

To assess the impact of choosing different models for calculating the water dielectric properties on mineral equilibria and fluid-rock interaction, we preformed several equilibrium calculations between a single mineral, mineral assemblages, simple granite rock composition on the one hand, and aqueous solutions of varying salinity on the other hand at different temperature and pressure conditions. We compare JN91, SV14, and FE97 electrostatic models; the IAPWS95 thermodynamic model was used for calculating water density and its derivatives.

First, we simulated the solubility of a single mineral (calcite and muscovite) and of a simple mineral assemblage (wollastonite plus quartz) in pure water and of another assemblage (enstatite plus forsterite) in the NaCl water solution (Figure 8). The solubility of calcite, muscovite, and wollastonite plus quartz was calculated from 200 until 800°C at 2 kbar. The enstatite plus forsterite assemblage was equilibrated at 700°C and 4 kbar with chlorine concentration increasing from 0.001 to 2 mol/kg H₂O.

(a)

(b)

(c)

(d)

Figure 8

Mineral solubility in water as a function of temperature at 2 kbar, calculated using the JN91, SV14, and FE97 models for water dielectric properties and IAPWS95 for water volumetric properties: (a) calcite solubility, (b) muscovite solubility, and (c) wollastonite + quartz solubility; (a, b, c) in pure water at 2 kbar with varying temperature; (d) enstatite + forsterite solubility in NaCl solution of varying molality at 700°C and 4 kbar. Total amounts of dissolved elements as well as the species concentrations and negative pH are shown. Symbols represent element concentrations measured in experiments.

In the single mineral solubility example (Figures 8(a) and 8(b)), differences in the total element concentrations as well as in the pH and in the amounts of aqueous species start to appear above 500 to 600°C. Using the SV14 and FE97 models results in increasingly lower mineral solubilities with increasing temperature, compared to the JN91 model. In the calcite solubility calculation, the discrepancies increase with the increasing charge of an ion (with concentration of Ca²⁺ having the largest and that of CaCO₃⁰ having the smallest discrepancy at 800°C, Figure 8(a)). The calculated total Ca concentration is in good agreement with the experimental data (400 to 600°C) of Fein and Walther [84]. For the muscovite solubility (Figure 8(b)), a large mismatch is observed in the case of Al(OH)₃⁰, which despite being a neutral aqueous complex has a rather large conventional Born coefficient () [27]. This leads to a difference of about 0.28 log units in the dissolved Al and K at 800°C. The pH starts to show differences only above 700°C, with 0.1 log units smaller at 800°C when using the SV14 and FE97 models, compared to the JN91 model.

In the case of wollastonite plus quartz solubility (Figure 8(c)), there is not much disagreement between the JN9, SV14, and FE97 models, with the only major difference being visible for the pH above 600°C, with a 0.5-unit difference at 800°C. The calculated silica and calcium concentrations are in good agreement with the experimental data (300 to 660°C) of Xie and Walther [85].

In the case of enstatite plus forsterite solubility (Figure 8(d)), for the total amount of dissolved Mg and pH, differences larger than 0.1 log units exist only at chlorine concentrations below 0.1 Cl mol/kg H₂O.

Secondly, the composition of a fluid in equilibrium with albite, K-feldspar, andalusite, and quartz mineral assemblage, with increasing chlorine concentration, was calculated (Figure 9). The system was equilibrated at 650°C and 2 kbar with Cl concentration ranging from 0.01 to 2 mol/kg H₂O. For the total amounts of dissolved elements in equilibrium with albite, K-feldspar, andalusite, and quartz mineral assemblage, no significant differences are visible when using the different electrostatic models (JN91, SV14, and FE97) for calculating the water dielectric properties (Figure 9(a)). Small discrepancies are apparent in the total dissolved aluminum at low Cl concentrations. The modeled data agrees well with the experimental points from Pak et al. [86]. Even though there is a scatter in the experimental data, the Si, Al, Na, and K concentrations calculated using the SV14 and JN91 models are less than 0.1 log units different. Figure 9(b) shows the calculated speciation for the same system (molality of individual species); differences in the concentrations of aqueous ions (e.g., Na⁺, K⁺) up to 0.25 log units can be observed; similar differences of opposite sign can be observed in the concentration of neutral complexes (e.g., NaCl⁰, KCl⁰).

(a)

(b)

Finally, we simulated the reaction of granite with a simple composition (37 wt% albite, 30 wt% quartz, 21 wt% K-feldspar, and 12 wt% muscovite) with a 1 NaCl mol/kg H₂O solution (Figure 10). The equilibrium was calculated at 600°C and 2 kbar (consitions at which significant dicrepanices exist between electrostatic models), from 0.1 to 600 water to solid ratios (W/S). The concentrations of total dissolved Na, K, Si, and Al, as well as the pH show consistent results when using the FE97, SV14, and JN91 models.

6. Discussion

6.1. Discrepancies between Calculated Water Solvent and Aqueous Species Properties

In the HKF model, the solvation contribution to the partial molal standard state properties of aqueous ions and complexes is calculated using the dielectric constant of the water, which is calculated using models that are functions of water density. Therefore, differences in the calculated water density, dielectric constant, and their derivatives with respect to temperature and pressure propagate in the values of the aqueous ions and species properties (notice similarities between Figures 1–5, and EA1–EA5 in the electronic supplementary material).

As presented in Results, overall, only small differences arise in the calculated standard state properties of aqueous ions and complexes when comparing HGK84, IAPWS95, and ZD05 formulations for the water density calculation. A considerable agreement exists between various models from the early HGK84 formulation up to present models [2, 67]. This is due to a large amount of available experimental data on the volumetric properties of water and a similarity in the functional form of the HGK84 and IAPWS95 thermodynamic models. Some notable differences can be seen at high and low water densities for the heat capacity and at low water density for the volume. Although the recommended lower pressure limit for ZD05 is set at 1 kbar, due to significant disagreements with the other two models for the calculated and , a more appropriate lower limit would be 2 kbar. The disagreements in the calculated density further translate into differences in the derivatives of the dielectric constant (Born functions) and in the standard state properties of aqueous species such as , , and (Figures 1 and 4). Differences between 1 and 3 kJ·mol⁻¹ in the calculated Gibbs energy between IAPWS95, ZD05, and HGK84, although relatively small, can be significant when estimating reference temperature and pressure values from high-temperature and -pressure measurements.

The investigated formulations for calculating the dielectric properties of water are semi-empirical fits to selected experimental data [46] with some additional theoretical considerations and to molecular dynamics simulation data beyond the experimental interval [9, 63]. Although the models extrapolate smoothly outside their recommended interval, they start to disagree beyond the upper temperature limit of the available experimental data (~500°C, experimental data of Heger et al. [77] up to 550°C and 5 kbar) [47]. The electrostatic models are functions of the water density and its derivatives and were developed using a supporting water thermodynamic model. Strictly speaking, one has to use the same model to calculate the water density that was used when fitting a dielectric properties model. However, as shown in the results (Figures 1, 2(a), 4, 6(a), and 6(b) and from the electronic suplementary material Figures EA1, EA2A, and EA4), there are no major differences when using different models for calculating the water density and its derivatives below 800°C. Thus, the impact of combining thermodynamic models for water, as discussed above, in the investigated temperature and pressure range, with any water electrostatic model, is minimal. For example, if water volumetric properties are calculated using HGK84, instead of IAPWS95, at 1000°C and 5 kbar, the difference for between JN91 and SV14 would be 0.7% smaller (see Figure 3(a)). Outside the investigated temperature pressure range (25 to 1000°C and 1 to 5 kbar) and outside the reported recommended conventional model bounds (Table 1), the differences may become significant [47].

The solvation contributions to the standard partial molal properties such as , , and ((6), (7), (8)) are dependent on the , , and Born functions ((9), (10), and (11)). Any differences in calculated values of the Born functions will propagate to the respective thermodynamic properties. Large differences can be seen in the calculated values for , , and at high temperatures, because the investigated electrostatic models disagree above the highest available temperature experimental point (500°C). Below 150°C, the disagreements between the calculated values for and using SV14 in comparison with FE97 and JN91, propagate to the values of and . The , , and Born functions are dependent on the first- and second-order derivatives of with respect to pressure and temperature. The models for calculating having different mathematical form were fitted to experimental values of and not the derivatives. It is recommended that a good agreement between the experimental data and the model should be maintained not only for but also for its derivatives. A model might appear to visually reproduce the experimental data, but due to its mathematical formulation and not a robust physical background can introduce larger errors in the property derivatives.

The standard state Gibbs energy of ions calculated using the HKF model is dependent on the inverse dielectric constant ((3)). The differences in the dielectric constant, obtained from different models that appear above 500°C, are translated into differences in the calculated standard state Gibbs energy. In the solvation contribution formulation ((3)), the inverse of the dielectric constant is multiplied with the conventional Born coefficient . This leads to an increase in the differences of the calculated values with increasing charge of the ion (absolute value) and with increasing value of the reference conventional Born coefficient () in the HKF model ((4) and (5)). The neutral species have zero charge and is usually set or is close, to zero, thus their solvation contribution is close to zero and the differences of the calculated values are only about a few hundreds of J/mol. By contrast, for charged species, the solvation contribution is not close to zero and increases with the absolute value of the charge and the Born term parameter (). This translates in differences of the calculated values that are on the order of thousands of J/mol. The aqueous speciation can be seen in the aqueous mineral equilibrium calculations. Where the concentration of an element is dominated by a neutral aqueous species (e.g., Figures 8(c), 8(d), and 9 for total Si with SiO₂⁰ dominant species), the discrepancies are small. In the systems where the concentration is dominated by charged ions or complexes, the discrepancies are visible (e.g., Figures 8(a) and 8(b) for total Ca, Al, and K).

Differences in the calculated water dissociation constant from the BL06 model can have two main sources. The first source is the dielectric constant model (JN91, FE97, and SV14) that is used to calculate the solvation contribution to the Gibbs energy of OH⁻ (e.g., Figure EA5, electronic suplementary material; as discussed in the previous paragraph). At high temperatures, the solvation contribution to the Gibbs energy has the largest impact [35]. If we consider that the BL06 mode is accurate, the large discrepancies above 600°C (Figure 7) suggest that the main source of mismatch should be the electrostatic models. But above 400°C, only the reported data of Quist [87] was used for parameterizing the BL06 equation with stated deviations of up to 0.82 log units from experimental data [83], still below the maximum discrepancy from different electrostatic models. The second source for the observed differences can be the HKF parameter set used for calculating the nonsolvation contribution to the Gibbs energy of OH⁻. [35]. The OH⁻ HKF parameters from the core dataset [35] were derived from electrolyte experimental data using the additivity principle ((2)) and are independent from measurements of water dissociation constant. To resolve the discrepancies with the BL06 model, the developers of the DEW dataset (Deep Earth Water model, v. 11.0.1) optimized the OH⁻ HKF parameters. They refined them against values for the water dissociation constant calculated with the BL06 model, leading to an improved agreement up to 700°C (Figure 7). One side effect of adjusting the HKF parameters of OH⁻ (a core species) is that it breaks the consistency with the electrolyte data as presented in Tanger and Helgeson [35]. For example, from (2), the sum of the parameter for Na⁺ and OH⁻, 18.18 and 4.15, is equal to 22.33. This is the value of for NaOH electrolyte obtained by Tanger and Helgeson [35] from heat capacity measurements. In the DEW dataset, the values of for Na⁺ and OH⁻ are 18.18 and 12.12. Their sum is 30.3, which is not anymore consistent with the original value for the NaOH electrolyte of 22.33.

In the geochemical literature, aqueous species are presented in their hydrated or dehydrated form. Conventionally, the dehydrated form of an aqueous species can be transformed to its hydrated counterpart by the reaction with water:

The thermodynamic properties of such a reaction are set to 0 by convention. Using different models for calculating the thermodynamic properties for water will result in differences in the calculated thermodynamic properties of the hydrated species. Calculating these properties using one water model and then doing equilibrium calculations with another one will break the consistency of reaction (19). If determined from a hydration reaction, the properties of hydrated species get a contribution from the water properties and their temperature and pressure dependence. Furthermore, if properties calculated using a reaction similar to 19 are used to fit HKF parameters, the water properties and any discrepancies in them are further propagated in the optimized HKF parameters [88, 89].

For a large number of aqueous species, the HKF parameters were generated from correlations [38, 49]. Using solubility and speciation data, the discrepancies between the modeled and measured data can be minimized by adjusting the reference standard state properties and the HKF parameters (parameters , , , , , and ) [27, 29, 58, 89]. If the solubility data are measured at conditions where the water dielectric constant models are in disagreement (e.g., above 500°C), then the discrepancies that result from the dielectric model used to calculate the solvation contributions will be fitted into the HKF parameters of the nonsolvation contribution to the standard state properties.

6.2. Discrepancies in Fluid Mineral Equilibria and Fluid-Rock Interaction Modeling

When modeling fluid-rock interaction processes, depending on the aqueous speciation, large or small discrepancies in the total calculated amounts of dissolved elements can be observed.

In the case of calcite and muscovite solubility calculations, the total concentrations of Ca, Al, and K are controlled by charged species such as Ca²⁺, Ca(OH)⁺, K⁺, and Al(OH)₄⁻ and species with relatively large conventional Born coefficients, e.g., Al(OH)₃⁰ [27]; therefore, increasing discrepancies are observed above 500°C (Figure 8(a) and 8(b)). Conversely, in the case of wollastonite plus quartz solubility (Figure 8(c)), above 500°C, the concentration of Ca is controlled by CaSiO₃⁰ species and only a small discrepancy is visible (similar to the concentration of Si). Below 500°C, CaOH⁺ controls the concentration of Ca, but at these temperatures, there is good agreement between different electrostatic models. The enstatite plus forsterite solubility with increasing chlorine concentration shows how the change in speciation affects the discrepancy between the calculations (Figure 8(d)). At low chlorine concentration where MgOH⁺ controls the Mg concentration, a ~0.3 log unit difference can be seen. With increasing chlorine concentration, the MgCl₂⁰ neutral species becomes dominant, resulting in almost identical Mg concentrations calculated with the three electrostatic models.

In the case of the mineral assemblage (low albite, sanidine K-feldspar, andalusite, and quartz, Figure 9) in equilibrium with a NaCl aqueous solution, the ratio of Na to K in the dissolved fluid is constrained by the following reaction:

The major species governing the concentration of dissolved Na, K, and Si are neutral, and this results in small differences for the calculated total element concentrations. Additionally, the K and Na species appear on both sides of the exchange reaction and differences in the values of aqueous species properties calculated by different models cancel out. This has an effect on mineral exchange and partitioning reactions, which involve elements with similar speciation, resulting in a reduced, sometimes even vanishing discrepancy between different water dielectric models (at ). Even though there is no visible difference in the total dissolved Na, K, and Si, there is a visible difference in the concentration of different species when using different models (Figure 9(b)). A difference in the speciation of major elements can have an impact on pH, which in turn will affect the aqueous composition of minor and trace elements such as ore metals. For example, if a species such as CaOH⁺ is more stable relative to Ca²⁺, more OH⁻ will be consumed in order to form the neutral species, thus reducing the pH (e.g., Figure 8(a), where above 600°C CaOH⁺ is more stable when using JN91, resulting in a smaller pH, when compared with SV14 and FE97). The calculated total Ca concentration in equilibrium with wollastonite and quartz is not affected by changing the water dielectric properties models, but due to the impact that different models have on speciation, a discrepancy in the calculated pH is observed (Figure 8(c)). Part of this is due to HSiO₃⁺ which is more stable when using the JN91 model. SiO₂⁰ reacts with OH⁻ forming HSiO₃⁺, resulting in a lower pH (dashed black line). In contrast, for the solubility of muscovite, the large discrepancy in calculated total Al and K concentrations is not visible in the case of the calculated pH. A difference in pH is only visible above 700°C.

The salinity of metamorphic and hydrothermal fluids can vary from very dilute to hypersaline [22, 90, 91], with the majority being around 3.2–4 wt% NaCl [91]. This has the effect that chloride species play a major role in the composition of aqueous fluids (e.g., NaCl⁰, KCl⁰, MgCl⁺, MgCl₂⁰, CaCl⁺, and CaCl₂⁰). With increasing temperature, complexation is favored and species such as NaCl⁰ and KCl⁰ become more stable compared to cations such as Na⁺ and K⁺. This was observed in the simulation modeling a granite reacting with a 1 mol/kg H₂O NaCl solution at different water-to-solid ratios (Figure 10). No significant differences resulted when using the JN91, SV14, or FE97 model, because the species controlling the elemental concentration were SiO₂⁰, NaCl⁰, KCl⁰, and NaAl(OH)₄⁰.

6.3. Impact of Using Different Models on the Activity Coefficients

Besides being used for calculating standard state properties of ions and complexes based on the Born functions (i.e., HKF model), the dielectric constant of water is also used to calculate the Debye-Hückel limiting law coefficients [34, 63]. These are necessary for calculating the activity coefficients, describing the asymptotic concentration dependence of the thermodynamic properties at infinite dilution [34, 63]. In order to obtain the Debye-Hückel coefficients, values of the first- and the second-order derivatives of with respect to temperature and pressure are needed, as well as the first- and the second-order derivatives of . Differences in the values of these derivatives that arise from applying different models propagate to the calculated activity coefficients (excess thermodynamic properties). However, in this study, we did not investigate the impact of using different models for the water properties calculation on the Debye-Hückel limiting law coefficients (i.e., the same activity model has been used in all cases).

6.4. HKF Drawbacks and Outlook

The main disadvantage of the HKF model is that it becomes increasingly inaccurate in predicting the properties of solutes close to the critical point of water and at low water solvent densities [47]. This makes HKF not applicable to model solute properties in vapor like water densities. The model has also been criticized for inaccurate prediction of the properties of nonelectrolyte solutes in the near critical and supercritical regions of water [48, 92, 93]. Additionally, the effect of the solute on the dielectric constant is not considered; i.e., the calculated properties are for pure water solvent. This prevents modeling of mixed solvent electrolyte solutions (e.g., CO₂-rich fluids; Galvez et al. [44]). As shown above, another disadavantage results from the discrepancies in the water electrostatic models, necessary in the HKF model to be able to calculate the standard state properties of aqueous ions and complexes.

One way to resolve the existing discrepancies between different models for calculating the dielectric constant of water is producing new experimental data. New data retrieved from experiments or from molecular dynamics simulations on the dielectric constant of water, especially at elevated temperatures and pressures, is highly desirable. The experimental data needs to be fitted to a model with robust physical background, avoiding data overfitting, and taking care that the derivatives of the function representing the property are correctly represented. In addition, the dielectric properties of water can be indirectly retrieved from solubility measurements at high-pressure and -temperature conditions. This can be done in simple systems such as quartz solubility (see McKenzie and Helgeson [94]) accompanied by spectroscopic studies on the aqueous speciation for an accurate description of the fluid chemistry. Furthermore, one could use the HKF parameters for OH⁻ that were independently derived from electrolyte data, at temperatures below 500°C (where the nonsolvation contributions are dominant) and the measured water dissociation constants at temperatures above 500°C (where the solvation electrostatic contributions are dominant) for back calculating the dielectric constant of water. If the HKF parameters derived from electrolyte data and the measured water dissociation constants are accurate enough, the retrieved values for the dielectric constant can be used to constrain electrostatic models at elevated temperatures.

Another way to avoid the existing disagreements introduced by different water dielectric constant models in the values of the aqueous species properties is the use of the so-called “density” models [2, 95]. These models are functions which describe the properties of reactions involving aqueous species and minerals at elevated temperatures and pressures using the density of water and several adjustable parameters, without requiring the dielectric constant of water, and better reproducing the standard state properties close to the critical point of water. However, to be able to have a comprehensive thermodynamic database, for the numerous species needed when modeling complex geochemical systems, a large number of experimental data points are needed to be able to fit the empirical parameters of the density models. In comparison, the HKF model parameters are largely correlated. This gives the possibility to generate sets of parameters for species with little or no experimental data available [9], but how well these parameters, resulting from correlations, reproduce experimental, and natural data, is still a matter of debate [29, 48].

7. Conclusions

We have studied how different models for evaluating thermodynamic and electrostatic properties of water affect the calculated standard state thermodynamic properties of aqueous species using the HKF model.

Overall, there is good agreement between the three investigated thermodynamic models for water (HGK84, IAPWS95, and ZD05) in a wide range of temperature and pressure conditions. However, increasing discrepancies were identified for the standard state properties at (for Gibbs energy, entropy, heat capacity, and volume) and water (for volume). A more appropriate pressure lower limit for the ZD05 model would be 2 kbar, while HGK84 and IAPWS95 are in good agreement between 1 and 5 kbar and 25 and 800°C. The general good match between these models is due to the large amount of high-quality experimental data available that were used in their parameterization and the similarity of the governing equations.

Large disagreement is observed in the calculated standard state properties of aqueous ions and complexes (using the HKF model) when different electrostatic models (JN91, FE97, and SV14) are used for calculating the dielectric constant of water and its temperature and pressure derivatives. The deviations are evident at temperatures above 500°C. These discrepancies arise from the lack of experimental data for the dielectric constant of water at higher temperatures and pressures. To resolve the existing differences, new data on the dielectric constant of water from direct or indirect measurements (e.g., mineral solubility, water dissociation) and from ab initio molecular dynamics are necessary.

The discrepancies in the calculated properties of aqueous ions and complexes increase with increasing absolute charge and the standard conventional Born coefficient of the species. Therefore, the use of different models has little impact on neutral species, which have zero charge and conventional Born coefficient close to zero (e.g., SiO₂⁰). For mineral exchange reactions and element partitioning (e.g., exchange between ions with similar hydrated ionic radius), the differences will be minimized, appearing with similar values on both sides of the reactions. When modeling aqueous mineral equilibrium, variation can be small or large depending on the major species that control the elemental concentration. An aqueous phase in equilibrium with minerals dominated by ions and charged species shows increasing deviation in the total concentration of elements above 500°C. Systems with aqueous fluids containing a variable amount of chlorine (e.g., NaCl, typical for many hydrothermal and metamorphic fluids) will be dominated at elevated temperatures by neutral species such as NaCl⁰, KCl⁰, CaCl₂⁰, and MgCl₂⁰ and therefore show very small discrepancies. Even if there is no difference in the calculated concentration of major elements because of the impact that different models have on the speciation, a difference in the calculated pH can result. This will have an impact on the calculated concentration of minor and trace elements such as metals and rare earth elements.

The shortage of thermodynamic data at higher temperatures inhibits the modeling of high-temperature magmatic-hydrothermal systems, demanding new mineral solubility experiments at such conditions. However, different electrostatic models will produce different values for the dielectric properties of water in the temperature range of interest (above 500°C), and these differences will be propagated to the nonsolvation HKF model parameters, evaluated from the solubility experiments. Therefore, care has to be taken when using mineral solubility experiments (done above 500°C) to retrieve standard state thermodynamic properties, such as solubility constants or Gibbs free energies of formation for aqueous species, and when modeling these data using different electrostatic models for water together with the HKF model.

The IAPWS95 and FE97 models for calculating the water thermodynamic properties and dielectric constant are the recommended standards. They should replace the HGK84 and JN91 models for modeling calculations at upper to middle crustal conditions. The core HKF dataset containing standard state properties and HKF parameters for ions and complexes [35] was built using experimental data retrieved at temperatures below 500°C, and the parameters for many species were subsequently derived from correlations built upon these data. To be fully consistent, a general update of the HKF parameter set due to exchanging water solvent models would seem necessary, but the differences in the models below 500°C are small. Unless the parameters were derived from experiments done at higher temperatures, the significant effort would only bring minor improvements. Efforts should be focused instead on resolving the discrepancies in the water dielectric constant above 500°C, or providing parameter sets (covering a rather large number of species) for models that do not require calculating the dielectric constant of water.

Data Availability

The data used to support the findings of this study are included and are referenced within the article.

Conflicts of Interest

The authors declare that there is no conflict of interest regarding the publication of this paper. Data supporting the results are available in the manuscript and at the cited references.

Acknowledgments

Support to G. D. Miron by the German Federal Ministry for Education and Research (BMBF), ThermAc project (02NUK039A), and to A. Yapparova by the Swiss National Science Foundation through grants (CRSII2160757) (COTHERM II - COmbined hydrological, geochemical and geophysical modelling of geotTHERMal systems II) and (200020_172851) (Key problems of heat and mass transfer in the Earth’s crust) is gratefully acknowledged. The authors would like to thank D. A. Kulik for his constant support and discussions on the topic. Comments from D. Dolejš and four anonymous reviewers have greatly contributed to improving the manuscript.

Supplementary Materials

Supplementary material contains analogous figures to the ones present in the manuscript: Figures S1, S2, and S3 are analogous to Figures 1, 2, and 3, but show the absolute differences in property values instead. Figures S4 and S5 are analogous to Figures 4 and 5, but show differences in the standard state properties for the OH species, calculated using different combinations of the thermodynamic and electrostatic models. (Supplementary Materials)

References

E. H. Oelkers, P. Benezeth, and G. S. Pokrovski, “Thermodynamic databases for water-rock interaction,” Reviews in Mineralogy and Geochemistry, vol. 70, no. 1, pp. 1–46, 2009.
View at: Publisher Site | Google Scholar
D. Dolejš and C. E. Manning, “Thermodynamic model for mineral solubility in aqueous fluids: theory, calibration and application to model fluid-flow systems,” Geofluids, vol. 10, no. 1-2, 40 pages, 2010.
View at: Publisher Site | Google Scholar
M. E. Galvez, O. Beyssac, I. Martinez et al., “Graphite formation by carbonate reduction during subduction,” Nature Geoscience, vol. 6, no. 6, pp. 473–477, 2013.
View at: Publisher Site | Google Scholar
M. E. Galvez, J. A. D. Connolly, and C. E. Manning, “Implications for metal and volatile cycles from the pH of subduction zone fluids,” Nature, vol. 539, no. 7629, pp. 420–424, 2016.
View at: Publisher Site | Google Scholar
C. E. Manning, “Mobilizing aluminum in crustal and mantle fluids,” Journal of Geochemical Exploration, vol. 89, no. 1–3, pp. 251–253, 2006.
View at: Publisher Site | Google Scholar
K. I. Shmulovich, B. W. D. Yardley, and G. G. Gonchar, Fluids in the Crust: Equilibrium and Transport Properties, Springer Netherlands, 1995.
W. S. Fyfe, N. J. Price, and A. B. Thompson, Fluids in the Earth's Crust: Their Significance in Metamorphic, Tectonic, and Chemical Transport Processes, Elsevier Scientific Pub. Co., 1978.
K. Kyser and E. E. Hiatt, “Fluids in sedimentary basins: an introduction,” Journal of Geochemical Exploration, vol. 80, no. 2-3, pp. 139–149, 2003.
View at: Publisher Site | Google Scholar
D. A. Sverjensky, B. Harrison, and D. Azzolini, “Water in the deep Earth: the dielectric constant and the solubilities of quartz and corundum to 60 kb and 1200 °C,” Geochimica et Cosmochimica Acta, vol. 129, pp. 125–145, 2014.
View at: Publisher Site | Google Scholar
A. Stefánsson, “Dissolution of primary minerals of basalt in natural waters: I. Calculation of mineral solubilities from 0°C to 350°C,” Chemical Geology, vol. 172, no. 3-4, pp. 225–250, 2001.
View at: Publisher Site | Google Scholar
A. P. Gysi, “Numerical simulations of CO2 sequestration in basaltic rock formations: challenges for optimizing mineral-fluid reactions,” Pure and Applied Chemistry, vol. 89, no. 5, pp. 581–596, 2017.
View at: Publisher Site | Google Scholar
J. C. Ayers and E. B. Watson, “Solubility of apatite, monazite, zircon, and rutile in supercritical aqueous fluids with implications for subduction zone geochemistry,” Philosophical Transactions of the Royal Society of London. Series A: Physical and Engineering Sciences, vol. 335, no. 1638, pp. 365–375, 1991.
View at: Publisher Site | Google Scholar
P. Y. Azimov and S. A. Bushmin, “Solubility of minerals of metamorphic and metasomatic rocks in hydrothermal solutions of varying acidity: thermodynamic modeling at 400–800°C and 1–5 kbar,” Geochemistry International, vol. 45, no. 12, pp. 1210–1234, 2007.
View at: Publisher Site | Google Scholar
P. Bénézeth, A. Stefánsson, Q. Gautier, and J. Schott, “Mineral solubility and aqueous speciation under hydrothermal conditions to 300 °C - the carbonate system as an example,” Reviews in Mineralogy and Geochemistry, vol. 76, no. 1, pp. 81–133, 2013.
View at: Publisher Site | Google Scholar
S. E. Drummond, Boiling and Mixing of Hydrothermal Fluids: Chemical Effects on Mineral Precipitation, PhD Thesis. Pennsylvania State University, 1981.
L. P. Baumgartner and J. M. Ferry, “A model for coupled fluid-flow and mixed-volatile mineral reactions with applications to regional metamorphism,” Contributions to Mineralogy and Petrology, vol. 106, no. 3, pp. 273–285, 1991.
View at: Publisher Site | Google Scholar
T. Baker, E. V. Achterberg, C. G. Ryan, and J. R. Lang, “Composition and evolution of ore fluids in a magmatic-hydrothermal skarn deposit,” Geology, vol. 32, no. 2, pp. 117–120, 2004.
View at: Publisher Site | Google Scholar
H. A. Alfredsson, B. S. Hardarson, H. Franzson, and S. R. Gislason, “CO₂ sequestration in basaltic rock at the Hellisheidi site in SW Iceland: stratigraphy and chemical composition of the rocks at the injection site,” Mineralogical Magazine, vol. 72, no. 1, pp. 1–5, 2008.
View at: Publisher Site | Google Scholar
S. Arnórsson, E. Gunnlaugsson, and H. Svavarsson, “The chemistry of geothermal waters in Iceland. II. Mineral equilibria and independent variables controlling water compositions,” Geochimica et Cosmochimica Acta, vol. 47, no. 3, pp. 547–566, 1983.
View at: Publisher Site | Google Scholar
D. Dolejš and T. Wagner, “Thermodynamic modeling of non-ideal mineral–fluid equilibria in the system Si–Al–Fe–Mg–Ca–Na–K–H–O–Cl at elevated temperatures and pressures: implications for hydrothermal mass transfer in granitic rocks,” Geochimica et Cosmochimica Acta, vol. 72, no. 2, pp. 526–553, 2008.
View at: Publisher Site | Google Scholar
L. Andre, N. Spycher, T. Xu, F. D. Vuataz, and K. Pruess, Modeling Brine-Rock Interactions in an Enhanced Geothermal System Deep Fractured Reservoir at Soultz-Sous-Forets (France): a Joint Approach Using Two Geochemical Codes: Frachem and ToughReact, p. 71, 2006, https://escholarship.org/uc/item/17b2550x.
G. D. Miron, T. Wagner, M. Wälle, and C. A. Heinrich, “Major and trace-element composition and pressure–temperature evolution of rock-buffered fluids in low-grade accretionary-wedge metasediments, Central Alps,” Contributions to Mineralogy and Petrology, vol. 165, no. 5, pp. 981–1008, 2013.
View at: Publisher Site | Google Scholar
T. J. Wolery, “On the thermodynamic framework of solutions (with special reference to aqueous electrolyte solutions),” American Journal of Science, vol. 290, no. 3, pp. 296–320, 1990.
View at: Publisher Site | Google Scholar
A. Yapparova, G. D. Miron, D. A. Kulik, G. Kosakowski, and T. Driesner, “An advanced reactive transport simulation scheme for hydrothermal systems modelling,” Geothermics, vol. 78, pp. 138–153, 2019.
View at: Publisher Site | Google Scholar
M. Engi, “Thermodynamic data for minerals: a critical assessment,” in The Stability of Minerals, G. Price and N. Ross, Eds., pp. 267–328, Kluwer Academic Publishers, 1992.
View at: Google Scholar
T. J. Wolery and C. F. Jové Colón, “Chemical thermodynamic data. 1. The concept of links to the chemical elements and the historical development of key thermodynamic data,” Geochimica et Cosmochimica Acta, vol. 213, pp. 635–676, 2017.
View at: Publisher Site | Google Scholar
G. D. Miron, T. Wagner, D. A. Kulik, and C. A. Heinrich, “Internally consistent thermodynamic data for aqueous species in the system Na–K–Al–Si–O–H–Cl,” Geochimica et Cosmochimica Acta, vol. 187, pp. 41–78, 2016.
View at: Publisher Site | Google Scholar
B. M. Tutolo, X. Z. Kong, W. E. Seyfried Jr, and M. O. Saar, “Internal consistency in aqueous geochemical data revisited: applications to the aluminum system,” Geochimica et Cosmochimica Acta, vol. 133, pp. 216–234, 2014.
View at: Publisher Site | Google Scholar
G. D. Miron, T. Wagner, D. A. Kulik, and B. Lothenbach, “An internally consistent thermodynamic dataset for aqueous species in the system Ca-Mg-Na-K-Al-Si-O-H-C-Cl to 800 °C and 5 kbar,” American Journal of Science, vol. 317, no. 7, pp. 755–806, 2017.
View at: Publisher Site | Google Scholar
C. A. J. Appelo, “Principles, caveats and improvements in databases for calculating hydrogeochemical reactions in saline waters from 0 to 200 °C and 1 to 1000 atm,” Applied Geochemistry, vol. 55, pp. 62–71, 2015.
View at: Publisher Site | Google Scholar
H. C. Helgeson and D. H. Kirkham, “Theoretical prediction of the thermodynamic properties of aqueous electrolytes at high pressures and temperatures. III. Equation of state for aqueous species at infinite dilution,” American Journal of Science, vol. 276, no. 2, pp. 97–240, 1976.
View at: Publisher Site | Google Scholar
H. C. Helgeson and D. H. Kirkham, “Theoretical prediction of the thermodynamic behavior of aqueous electrolytes at high pressures and temperatures; I, summary of the thermodynamic/electrostatic properties of the solvent,” American Journal of Science, vol. 274, no. 10, pp. 1089–1198, 1974.
View at: Publisher Site | Google Scholar
H. C. Helgeson and D. H. Kirkham, “Theoretical prediction of the thermodynamic behavior of aqueous electrolytes at high pressures and temperatures; II, Debye-Huckel parameters for activity coefficients and relative partial molal properties,” American Journal of Science, vol. 274, no. 10, pp. 1199–1261, 1974.
View at: Publisher Site | Google Scholar
H. C. Helgeson, D. H. Kirkham, and G. C. Flowers, “Theoretical prediction of the thermodynamic behavior of aqueous electrolytes by high pressures and temperatures; IV, calculation of activity coefficients, osmotic coefficients, and apparent molal and standard and relative partial molal properties to 600 degrees C and 5 kb,” American Journal of Science, vol. 281, no. 10, pp. 1249–1516, 1981.
View at: Publisher Site | Google Scholar
J. C. Tanger and H. C. Helgeson, “Calculation of the thermodynamic and transport properties of aqueous species at high pressures and temperatures; revised equations of state for the standard partial molal properties of ions and electrolytes,” American Journal of Science, vol. 288, no. 1, pp. 19–98, 1988.
View at: Publisher Site | Google Scholar
E. H. Oelkers, H. C. Helgeson, E. L. Shock, D. A. Sverjensky, J. W. Johnson, and V. A. Pokrovskii, “Summary of the apparent standard partial molal Gibbs free energies of formation of aqueous species, minerals, and gases at pressures 1 to 5000 bars and temperatures 25 to 1000 °C,” Journal of Physical and Chemical Reference Data, vol. 24, no. 4, pp. 1401–1560, 1995.
View at: Publisher Site | Google Scholar
E. L. Shock, H. C. Helgeson, and D. A. Sverjensky, “Calculation of the thermodynamic and transport properties of aqueous species at high pressures and temperatures: standard partial molal properties of inorganic neutral species,” Geochimica et Cosmochimica Acta, vol. 53, no. 9, pp. 2157–2183, 1989.
View at: Publisher Site | Google Scholar
E. L. Shock, D. C. Sassani, M. Willis, and D. A. Sverjensky, “Inorganic species in geologic fluids: correlations among standard molal thermodynamic properties of aqueous ions and hydroxide complexes,” Geochimica et Cosmochimica Acta, vol. 61, no. 5, pp. 907–950, 1997.
View at: Publisher Site | Google Scholar
D. A. Sverjensky, E. L. Shock, and H. C. Helgeson, “Prediction of the thermodynamic properties of aqueous metal complexes to 1000°C and 5 kb,” Geochimica et Cosmochimica Acta, vol. 61, no. 7, pp. 1359–1412, 1997.
View at: Publisher Site | Google Scholar
J. W. Johnson, E. H. Oelkers, and H. C. Helgeson, “SUPCRT92: a software package for calculating the standard molal thermodynamic properties of minerals, gases, aqueous species, and reactions from 1 to 5000 bar and 0 to 1000°C,” Computers & Geosciences, vol. 18, no. 7, pp. 899–947, 1992.
View at: Publisher Site | Google Scholar
M. D. Schulte, E. L. Shock, and R. H. Wood, “The temperature dependence of the standard-state thermodynamic properties of aqueous nonelectrolytes,” Geochimica et Cosmochimica Acta, vol. 65, no. 21, pp. 3919–3930, 2001.
View at: Publisher Site | Google Scholar
S. Facq, I. Daniel, G. Montagnac, H. Cardon, and D. A. Sverjensky, “In situ Raman study and thermodynamic model of aqueous carbonate speciation in equilibrium with aragonite under subduction zone conditions,” Geochimica et Cosmochimica Acta, vol. 132, pp. 375–390, 2014.
View at: Publisher Site | Google Scholar
S. Facq, I. Daniel, G. Montagnac, H. Cardon, and D. A. Sverjensky, “Carbon speciation in saline solutions in equilibrium with aragonite at high pressure,” Chemical Geology, vol. 431, pp. 44–53, 2016.
View at: Publisher Site | Google Scholar
M. E. Galvez, C. E. Manning, J. A. D. Connolly, and D. Rumble, “The solubility of rocks in metamorphic fluids: a model for rock-dominated conditions to upper mantle pressure and temperature,” Earth and Planetary Science Letters, vol. 430, pp. 486–498, 2015.
View at: Publisher Site | Google Scholar
B. Debret and D. A. Sverjensky, “Highly oxidising fluids generated during serpentinite breakdown in subduction zones,” Scientific Reports, vol. 7, no. 1, p. 10351, 2017.
View at: Publisher Site | Google Scholar
J. W. Johnson and D. Norton, “Critical phenomena in hydrothermal systems; state, thermodynamic, electrostatic, and transport properties of H₂O in the critical region,” American Journal of Science, vol. 291, no. 6, pp. 541–648, 1991.
View at: Publisher Site | Google Scholar
D. Dolejš, “Thermodynamics of aqueous species at high temperatures and pressures: equations of state and transport theory,” Reviews in Mineralogy and Geochemistry, vol. 76, no. 1, pp. 35–79, 2013.
View at: Publisher Site | Google Scholar
A. V. Plyasunov and E. L. Shock, “Correlation strategy for determining the parameters of the revised Helgeson-Kirkham-Flowers model for aqueous nonelectrolytes,” Geochimica et Cosmochimica Acta, vol. 65, no. 21, pp. 3879–3900, 2001.
View at: Publisher Site | Google Scholar
E. L. Shock and H. C. Helgeson, “Calculation of the thermodynamic and transport properties of aqueous species at high pressures and temperatures: correlation algorithms for ionic species and equation of state predictions to 5 kb and 1000 °C,” Geochimica et Cosmochimica Acta, vol. 52, no. 8, pp. 2009–2036, 1988.
View at: Publisher Site | Google Scholar
D. A. Kulik, T. Wagner, S. V. Dmytrieva et al., “GEM-Selektor geochemical modeling package: revised algorithm and GEMS3K numerical kernel for coupled simulation codes,” Computational Geosciences, vol. 17, no. 1, pp. 1–24, 2013.
View at: Publisher Site | Google Scholar
A. M. M. Leal, D. A. Kulik, and G. Kosakowski, “Computational methods for reactive transport modeling: a Gibbs energy minimization approach for multiphase equilibrium calculations,” Advances in Water Resources, vol. 88, pp. 231–240, 2016.
View at: Publisher Site | Google Scholar
J. A. D. Connolly, “The geodynamic equation of state: what and how,” Geochemistry, Geophysics, Geosystems, vol. 10, no. 10, 2009.
View at: Publisher Site | Google Scholar
G. Eriksson, K. Hack, and S. Petersen, “ChemApp - a programmable thermodynamic calculation interface,” in Werkstoffwoche '96, Symposium 8 Simulation, Modellierung, Informationssysteme, DGM Informationsgesellschaft mbH, Dresden, Germany, 1997.
View at: Google Scholar
D. L. Parkhurst and C. A. J. Appelo, “Description of input and examples for PHREEQC version 3—a computer program for speciation, batch-reaction, one-dimensional transport, and inverse geochemical calculations,” U.S. Geological Survey Techniques and Methods, vol. 6, no. A43, p. 497, 2013.
View at: Google Scholar
C. M. Bethke, Geochemical and Biogeochemical Reaction Modeling, 2008, Cambridge University Press, p. 547.
T. J. Wolery, EQ3/6, a Software Package for Geochemical Modeling of Aqueous Systems: Package Overview and Installation Guide (Version 70), Lawrence Livermore National Laboratory, 1992.
M. H. Reed, “Calculation of simultaneous chemical equilibria in aqueous-mineral-gas systems and Its application to modeling hydrothermal processes,” in Techniques in Hydrothermal Ore Deposits Geology, J. P. Richards and P. B. Larson, Eds., pp. 109–124, Society of Economic Geologists, 1998.
View at: Google Scholar
B. Tagirov and J. Schott, “Aluminum speciation in crustal fluids revisited,” Geochimica et Cosmochimica Acta, vol. 65, no. 21, pp. 3965–3992, 2001.
View at: Publisher Site | Google Scholar
Y. Mei, W. Liu, J. Brugger, D. M. Sherman, and J. D. Gale, “The dissociation mechanism and thermodynamic properties of HCl(aq) in hydrothermal fluids (to 700 °C, 60 kbar) by ab initio molecular dynamics simulations,” Geochimica et Cosmochimica Acta, vol. 226, pp. 84–106, 2018.
View at: Publisher Site | Google Scholar
Y. Mei, W. Liu, D. M. Sherman, and J. Brugger, “Metal complexation and ion hydration in low density hydrothermal fluids: ab initio molecular dynamics simulation of Cu(I) and Au(I) in chloride solutions (25–1000 °C, 1–5000 bar),” Geochimica et Cosmochimica Acta, vol. 131, pp. 196–212, 2014.
View at: Publisher Site | Google Scholar
M. Born, “Volumen und Hydratationswärme der Ionen,” Zeitschrift für Physik, vol. 1, no. 1, pp. 45–48, 1920.
View at: Publisher Site | Google Scholar
D. P. Fernández, Y. Mulev, A. R. H. Goodwin, and J. M. H. L. Sengers, “A database for the static dielectric constant of water and steam,” Journal of Physical and Chemical Reference Data, vol. 24, no. 1, pp. 33–70, 1995.
View at: Publisher Site | Google Scholar
D. P. Fernández, A. R. H. Goodwin, E. W. Lemmon, J. M. H. Levelt Sengers, and R. C. Williams, “A formulation for the static permittivity of water and steam at temperatures from 238 K to 873 K at pressures up to 1200 MPa, including derivatives and Debye–Hückel coefficients,” Journal of Physical and Chemical Reference Data, vol. 26, no. 4, pp. 1125–1166, 1997.
View at: Publisher Site | Google Scholar
M. Uematsu and E. U. Frank, “Static dielectric constant of water and steam,” Journal of Physical and Chemical Reference Data, vol. 9, no. 4, pp. 1291–1306, 1980.
View at: Publisher Site | Google Scholar
E. U. Franck, S. Rosenzweig, and M. Christoforakos, “Calculation of the dielectric constant of water to 1000°C and very high pressures,” Berichte der Bunsengesellschaft für physikalische Chemie, vol. 94, no. 2, pp. 199–203, 1990.
View at: Publisher Site | Google Scholar
M. L. Marshall, “Dielectric constant of water discovered to be simple function of density over extreme ranges from −35 to +600 C and to 1200 MPa (12000 Atm.), believed universal,” Nature Precedings, pp. 1–40, 2008, http://precedings.nature.com/documents/2472/version/1.
View at: Google Scholar
W. Wagner and A. Pruß, “The IAPWS formulation 1995 for the thermodynamic properties of ordinary water substance for general and scientific use,” Journal of Physical and Chemical Reference Data, vol. 31, no. 2, pp. 387–535, 2002.
View at: Publisher Site | Google Scholar
Z. Zhang and Z. Duan, “Prediction of the PVT properties of water over wide range of temperatures and pressures from molecular dynamics simulation,” Physics of the Earth and Planetary Interiors, vol. 149, no. 3-4, pp. 335–354, 2005.
View at: Publisher Site | Google Scholar
C. Zhang and Z. Duan, “A model for C–O–H fluid in the Earth’s mantle,” Geochimica et Cosmochimica Acta, vol. 73, no. 7, pp. 2089–2102, 2009.
View at: Publisher Site | Google Scholar
L. Haar, J. G. Gallagher, and G. S. Kell, NBS/NRC Steam Tables, Hemisphere pub., Washington, D. C., 1984.
E. H. Oelkers and H. C. Helgeson, “Calculation of the thermodynamic and transport properties of aqueous species at high pressures and temperatures: aqueous tracer diffusion coefficients of ions to 1000 °C and 5 kb,” Geochimica et Cosmochimica Acta, vol. 52, no. 1, pp. 63–85, 1988.
View at: Publisher Site | Google Scholar
E. L. Shock, E. H. Oelkers, J. W. Johnson, D. A. Sverjensky, and H. C. Helgeson, “Calculation of the thermodynamic properties of aqueous species at high pressures and temperatures. Effective electrostatic radii, dissociation constants and standard partial molal properties to 1000 °C and 5 kbar,” Journal of the Chemical Society, Faraday Transactions, vol. 88, no. 6, pp. 803–826, 1992.
View at: Publisher Site | Google Scholar
F. Millero, “The partial molar volumes of electrolytes in aqueous solutions,” Water and Aqueous Solution, pp. 519–564, 1972.
View at: Google Scholar
B. R. Tagirov, A. V. Zotov, and N. N. Akinfiev, “Experimental study of dissociation of HCl from 350 to 500°C and from 500 to 2500 bars: thermodynamic properties of HCl°(aq),” Geochimica et Cosmochimica Acta, vol. 61, no. 20, pp. 4267–4280, 1997.
View at: Publisher Site | Google Scholar
J. M. H. L. Sengers, B. Kamgar-Parsi, F. W. Balfour, and J. V. Sengers, “Thermodynamic properties of steam in the critical region,” Journal of Physical and Chemical Reference Data, vol. 12, no. 1, pp. 1–28, 1983.
View at: Publisher Site | Google Scholar
K. S. Pitzer, “Dielectric constant of water at very high temperature and pressure,” Proceedings of the National Academy of Sciences, vol. 80, no. 14, pp. 4575-4576, 1983.
View at: Publisher Site | Google Scholar
K. Heger, M. Uematsu, and E. U. Franck, “The static dielectric constant of water at high pressures and temperatures to 500 MPa and 550°C,” Berichte der Bunsengesellschaft für Physikalische Chemie, vol. 84, no. 8, pp. 758–762, 1980.
View at: Publisher Site | Google Scholar
D. Pan, L. Spanu, B. Harrison, D. A. Sverjensky, and G. Galli, “Dielectric properties of water under extreme conditions and transport of carbonates in the deep earth,” Proceedings of the National Academy of Sciences, vol. 110, no. 17, pp. 6646–6650, 2013.
View at: Publisher Site | Google Scholar
S. W. Benson, Thermochemical Kinetics. Methods for the Estimation of Thermochemical Data and Rate Parameters, John Wiley & Sons, New York, 1968.
H. C. Helgeson, J. M. Delaney, H. W. Nesbitt, and D. K. Bird, “Summary and critique of the thermodynamic properties of rock-forming minerals,” American Journal of Science, vol. 278A, p. 230, 1978.
View at: Google Scholar
V. A. Pokrovskii and H. C. Helgeson, “Thermodynamic properties of aqueous species and the solubilities of minerals at high pressures and temperatures; the system Al₂O₃-H₂O-NaCl,” American Journal of Science, vol. 295, no. 10, pp. 1255–1342, 1995.
View at: Publisher Site | Google Scholar
A. M. M. Leal, M. J. Blunt, and T. C. LaForce, “Efficient chemical equilibrium calculations for geochemical speciation and reactive transport modelling,” Geochimica et Cosmochimica Acta, vol. 131, pp. 301–322, 2014.
View at: Publisher Site | Google Scholar
A. V. Bandura and S. N. Lvov, “The ionization constant of water over wide ranges of temperature and density,” Journal of Physical and Chemical Reference Data, vol. 35, no. 1, pp. 15–30, 2006.
View at: Publisher Site | Google Scholar
J. B. Fein and J. V. Walther, “Calcite solubility in supercritical CO₂-H₂O fluids,” Geochimica et Cosmochimica Acta, vol. 51, no. 6, pp. 1665–1673, 1987.
View at: Publisher Site | Google Scholar
Z. Xie and J. V. Walter, “Wollastonite + quartz solubility in supercritical NaCl aqueous solutions,” American Journal of Science, vol. 293, no. 3, pp. 235–255, 1993.
View at: Publisher Site | Google Scholar
T. M. Pak, C. A. Hauzenberger, and L. P. Baumgartner, “Solubility of the assemblage albite+K-feldspar+andalusite+quartz in supercritical aqueous chloride solutions at 650 °C and 2 kbar,” Chemical Geology, vol. 200, no. 3-4, pp. 377–393, 2003.
View at: Publisher Site | Google Scholar
A. S. Quist, “Ionization constant of water to 800.deg. and 4000 bars,” The Journal of Physical Chemistry, vol. 74, no. 18, pp. 3396–3402, 1970.
View at: Publisher Site | Google Scholar
G. D. Miron, D. A. Kulik, S. V. Dmytrieva, and T. Wagner, “GEMSFITS: code package for optimization of geochemical model parameters and inverse modeling,” Applied Geochemistry, vol. 55, pp. 28–45, 2015.
View at: Publisher Site | Google Scholar
Y. Shvarov, “A suite of programs, OptimA, OptimB, OptimC, and OptimS compatible with the Unitherm database, for deriving the thermodynamic properties of aqueous species from solubility, potentiometry and spectroscopy measurements,” Applied Geochemistry, vol. 55, pp. 17–27, 2015.
View at: Publisher Site | Google Scholar
B. W. D. Yardley and J. T. Graham, “The origins of salinity in metamorphic fluids,” Geofluids, vol. 2, no. 4, 256 pages, 2002.
View at: Publisher Site | Google Scholar
P. Nehlig, “Salinity of oceanic hydrothermal fluids: a fluid inclusion study,” Earth and Planetary Science Letters, vol. 102, no. 3-4, pp. 310–325, 1991.
View at: Publisher Site | Google Scholar
N. N. Akinfiev and L. W. Diamond, “Thermodynamic description of aqueous nonelectrolytes at infinite dilution over a wide range of state parameters,” Geochimica et Cosmochimica Acta, vol. 67, no. 4, pp. 613–629, 2003.
View at: Publisher Site | Google Scholar
J. Sedlbauer, J. P. O'Connell, and R. H. Wood, “A new equation of state for correlation and prediction of standard molal thermodynamic properties of aqueous species at high temperatures and pressures,” Chemical Geology, vol. 163, no. 1–4, pp. 43–63, 2000.
View at: Publisher Site | Google Scholar
W. F. McKenzie and H. C. Helgeson, “Estimation of the dielectric constant of H₂O from experimental solubilities of quartz, and calculation of the thermodynamic properties of aqueous species to 900 °C at 2 kb,” Geochimica et Cosmochimica Acta, vol. 48, no. 11, pp. 2167–2177, 1984.
View at: Publisher Site | Google Scholar
G. M. Anderson, S. Castet, J. Schott, and R. E. Mesmer, “The density model for estimation of thermodynamic parameters of reactions at high temperatures and pressures,” Geochimica et Cosmochimica Acta, vol. 55, no. 7, pp. 1769–1779, 1991.
View at: Publisher Site | Google Scholar

Copyright

Copyright © 2019 George D. Miron et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

PDF Download Citation

Download other formats

Order printed copies

Views

5400

Downloads

2622

Citations

Geofluids

Advances in Numerical Simulations of Hydrothermal Ore Forming Processes

Thermodynamic Properties of Aqueous Species Calculated Using the HKF Model: How Do Different Thermodynamic and Electrostatic Models for Solvent Water Affect Calculated Aqueous Properties?

Abstract

1. Introduction

2. Thermodynamic Models

2.1. Helgeson-Kirkham-Flowers Model

2.2. Haar, Gallagher, and Kell Water Equation of State (HGK84)

2.3. Wagner and Pruß Water Equation of State (IAPWS95)

2.4. Zhang and Duan Water Equation of State (ZD05)

3. Electrostatic Models

3.1. Johnson and Norton Water Electrostatic Model (JN91)

3.2. Fernandez et al. Water Electrostatic Model (FE97)

3.3. Sverjensky et al. Water Electrostatic Model (SV14)

4. Setup and Methods

5. Results

5.1. Comparison between Calculated Water Solvent Properties

5.2. Comparison of the Calculated Standard Properties of the Al3+ and OH− Ions and of the Muscovite Forming and Water Dissociation Reactions

5.3. Effect of Using Different Electrostatic Models in the Modeling of Fluid Mineral Equilibria and Fluid-Rock Interaction

6. Discussion

6.1. Discrepancies between Calculated Water Solvent and Aqueous Species Properties

6.2. Discrepancies in Fluid Mineral Equilibria and Fluid-Rock Interaction Modeling

6.3. Impact of Using Different Models on the Activity Coefficients

6.4. HKF Drawbacks and Outlook

7. Conclusions

Data Availability

Conflicts of Interest

Acknowledgments

Supplementary Materials

References

Copyright

5.2. Comparison of the Calculated Standard Properties of the Al³⁺ and OH⁻ Ions and of the Muscovite Forming and Water Dissociation Reactions