The dissociation mechanism and thermodynamic properties of HCl_(aq) in hydrothermal fluids (to 700 °C, 60 kbar) by ab initio molecular dynamics simulations

Abstract

HCl is one of the most significant volatiles in the Earth’s crust. It is well established that chloride activity and acidity (pH) play important roles in controlling the solubility of metals in aqueous hydrothermal fluids. Thus, quantifying the dissociation of HCl in aqueous solutions over a wide range of temperature and pressure is crucial for the understanding and numerical modeling of element mobility in hydrothermal fluids. Here we have conducted ab initio molecular dynamics (MD) simulations to investigate the mechanism of HCl_(aq) dissociation and to calculate the thermodynamic properties for the dissociation reaction at 25–700 °C, 1 bar to 60 kbar, i.e. including high temperature and pressure conditions that are geologically important, but difficult to investigate via experiments.

Our results predict that HCl_(aq) tends to associate with increasing temperature, and dissociate with increasing pressure. In particular, HCl_(aq) is highly dissociated at extremely high pressures, even at high temperatures (e.g., 60 kbar, 600–700 °C). At 25 °C, the calculated logK_d values (6.79 ± 0.81) are close to the value (7.0) recommended by IUPAC (International Union of Pure and Applied Chemistry) and some previous experimental and theoretical studies (Simonson et al.., 1990; Sulpizi and Sprik, 2008, 2010). The MD simulations indicate full dissociation of HCl at low temperature; in contrast, some experiments were interpreted assuming significant association at high HCl concentrations (≥1 m HCl_tot) even at room T (logK_d ∼0.7; e.g., Ruaya and Seward, 1987; Sretenskaya, 1992; review in Tagirov et al., 1997). This discrepancy is most likely the result of difficulties in the experimental determination of minor (if any) concentration of associated HCl_(aq) under ambient conditions, and thus reflects differences in the activity models used for the interpretation of the experiments. With increasing temperature, the discrepancy between our MD results and previous experimental studies, and between different studies, becomes smaller as the degree of HCl association increases. The MD simulations and available experimental studies show consistent results at hydrothermal conditions (300–700 °C, up to 5 kbar). The new thermodynamic properties based on the MD results provide an independent check of the dissociation constants for HCl_(aq), and the first dataset on HCl dissociation in high P-T fluids (up to 60 kbar, 700 °C) beyond available experimental conditions. Our results will enable prediction of the role of HCl in controlling element mobility in deep earth hydrothermal systems, including fluids associated with ultra-high pressure metasomatism in subduction zones.

Introduction

Hydrogen chloride (HCl) is one of the most significant volatiles in the Earth’s crust. In aqueous solutions, the neutral HCl_(aq) ion pair dissociates to the hydrogen ion (H⁺) and chloride anion (Cl⁻) through a reversible chemical reaction:

$${\text{HCl}}_{(\text{aq})}\rightleftharpoons {\text{H}}_{(\text{aq})}^{+}+{\text{Cl}}_{(\text{aq})}^{-}$$

Both H⁺ and Cl⁻ are important in controlling mineral solubility and element mobility in natural and man-made hydrothermal fluids. The chloride anion (Cl⁻) is one of the most abundant anions in hydrothermal fluids, and plays an important role in promoting element transport in hydrothermal fluids by forming stable chloro-complexes with many metals (e.g., Cu⁺, Co²⁺, Zn²⁺, Pb²⁺, Pd²⁺, etc.), whereas pH is an essential control of mineral solubility and element speciation in aqueous solutions (as reviewed in Seward and Barnes, 1997, Brugger et al., 2016). A recent study also reveals the importance of the pH of subduction zone fluids in controlling metal and volatile cycles (Galvez et al., 2016). Hence, the knowledge of the dissociation mechanism and thermodynamic properties of HCl_(aq) over a wide range of temperature and pressure is crucial for the understanding of element transport and deposition in hydrothermal fluids.

The logarithm of the dissociation constant of HCl_(aq) (Eq. (1), logK_d) can be expressed as;

$$\log K_d(P\text{,}T)=\log a_{H^{+}}+\log a_{{\mathit{Cl}}^{-}}-\log a_{{\mathit{HCl}}_{(\mathit{aq})}}$$

where a refers to the activities of the subscripted species. As summarized in Table 1, experiments up to 700 °C/4 kbar (one experiment to 4 kbar, and two to 2 kbar) and thermodynamic extrapolations have quantified the dissociation of HCl in aqueous fluids. At low temperature, significant discrepancies in the dissociation constant (up to 6 orders of magnitude) exist between different studies (see Table 1); as a result, some studies assume that HCl is significantly associated at room-temperature in concentrated HCl solutions, whereas other see no significant association. For instance, a few experimental studies obtained the dissociation constant of HCl_(aq) by measuring the partial vapor pressure to determine the activity of HCl_(aq) ($a_{{\mathit{HCl}}_{(\mathit{aq})}}$) (e.g., Robinson and Bates, 1971, Marsh and McElroy, 1985). In these studies, the equilibrium of HCl_(aq) ⇌ HCl_(g) was calculated based on Raoult’s Law, yielding a value of logK_d (HCl_(aq)) around 6.2 at 25 °C. Textbooks (e.g., p. 952 of Greenwood and Earnshaw, 1984) and the IUPAC dataset (p. 46 of Perrin, 1982) list a similar logK_d value of 7, making HCl a strong acid with no association at room temperature. However, Clegg and Brimblecombe (1986) pointed out that the approach based on Raoult’s Law is inappropriate for strong acids like HCl, and the Henry’s Law constant is a better approach to describe equilibrium between the vapor pressure of HCl solution and dissociated H⁺ and Cl⁻ ions.

In contrast to the above, several studies based on potentiometric and conductance data (Table 1) gave values as low as 0.9 (e.g., Johnson and Pytkowicz, 1978, Sretenskaya, 1992), implying significant ion pairing even at room temperature and suggesting that HCl is a relatively weak acid. The extrapolation of Ruaya and Seward (1987)’s AgCl_(s) solubility measurements at 100–350 °C also indicated a logK_d value of 0.7 at room temperature, and Sverjensky et al. (1991) recommended a logK_d value at 0.86 based on a fit of the experimental solubility data for alkali minerals.

Simonson et al. (1990) showed that the discrepancies in logK_d of HCl_(aq) can be due to the different chemical models employed to calculate the dissociation constants based on the available experimental data. Based on Holmes et al. (1987)’s excess thermodynamic data, Simonson et al. (1990) calculated two sets of the dissociation constants of HCl_(aq) using two types of data treatment (as listed in Table 1), and obtained a logK_d at 25 °C of 7.6 based on the ion association-interaction model (Pitzer’s model), and 0.2 based on the activity expansion-chemical equilibrium model.

The logK_d from different studies generally converge with increasing temperature to mostly overlap at temperatures ≥ 300 °C (Tagirov et al., 1997, Pokrovskii, 1999, Ho et al., 2001). This reflects the increased association of HCl making measurements more achievable. There are still some small discrepancies, in particular between mineral solubility and conductance measurements, and in low-density (<0.6 g cm⁻³) fluids. For example, there is up to 1 order of magnitude difference among the logK_d values reported by Franck, 1956, Frantz and Marshall, 1984 and Sverjensky et al. (1991).

The Helgeson-Kirkham-Flowers (HKF) model has been widely used to fit logK_d(HCl_(aq)) at high temperature (e.g., Sverjensky et al., 1991, Tagirov et al., 1997, Pokrovskii, 1999). The logK_d values extrapolated by the HKF model, based on mineral equilibrium data for K-feldspar-muscovite-quartz (Sverjensky et al., 1991) and AgCl(s) solubility experiments (Tagirov et al., 1997), are within a range of half a log unit at high temperature (300–600 °C, 0.5–2 kbar).

Over the past few decades, ab initio MD has been playing an increasingly significant role in understanding the molecular structures and thermodynamics of aqueous systems over a wide range of T, P, and solution compositions. In particular, the dissociation of acids at room temperature has been investigated quantitatively using such ab initio MD approaches. Trout and Parrinello (1998) first calculated the free energy of H₂O dissociation using ab initio distance-constrained thermodynamic integration. Although their results are 16.7 kJ/mol less than the experimental values due to the choice of the density functional and finite size (32 water molecules) of the system due to limited computing power at the time, they first demonstrated the capability of quantitative prediction of thermodynamic properties for aqueous species. Later on, Sprik (2000) calculated the logK_d of liquid water at ambient conditions using ab initio thermodynamic integration based on the coordination number of a given oxygen by hydrogen atoms and obtained values that are within one log unit of the experimental value.

In recent years, new techniques have been used to provide a molecular-level understanding of the association of HCl_(aq). Chialvo et al., 2002, Chialvo et al., 2003, Chialvo and Simonson, 2007 conducted potential of mean force (PMF) calculations based on classical force-field molecular dynamics simulation to predict the dissociation constants as a function of solution density. They used Lennard-Jones potentials between the oxygen in rigid H₂O and H₃O⁺ molecules and the Cl⁻ ion to model intermolecular interactions, and obtained reasonable representations of properties such as dielectric constant of water at near critical conditions and the logK_d values as a function of density changes. However, those studies did not provide direct association/dissociation constants as a function of temperature and pressure. Murakhtina et al. (2006) used density functional theory (DFT) based ab initio MD to calculate the H NMR chemical shift of HCl_(aq), and showed good agreement with experimental results. Sulpizi and Sprik, 2008, Sulpizi and Sprik, 2010 applied ab initio thermodynamic integration by constraining a dummy proton to calculate the free energy of dissociation of a series of acids; they obtained excellent agreement for logK_d (57 °C) values to within 0.4 log units for several small acid molecules (HCl, H₂S, formic acid), and within 2 log units for the wide range of acids investigated. They obtained logK_d values of 6.7–7.1 for HCl_(aq) at 57 °C (330 K), in excellent agreement with the IUPAC value of 7 at 25 °C (Perrin, 1982). Tummanapelli and Vasudevan (2014) calculated pKa values of weak organic acids at a temperature of 57 °C (330 K) using ab initio metadynamics simulations and achieved accurate predictions (within 0.2 log units).

Recently ab initio MD has also been widely used in molecular-level understanding of the geometry of metal complexes and the energetics of ligand exchange reactions relevant for metal transport in hydrothermal fluids (e.g., Cu(I)–Cl⁻, Sherman, 2007, Mei et al., 2014; U(VI)-Cl⁻, Bühl et al., 2008; Au(I)–HS⁻, Liu et al., 2011; Cu(I)–HS⁻–Cl⁻, Mei et al., 2013a; Ag(I)–Cl⁻, Liu et al., 2012, Pokrovski et al., 2013; Ag(I)-HS⁻/OH⁻, He et al., 2016; Au(I)–HS⁻/OH⁻/S3⁻, Mei et al., 2013b, Pokrovski et al., 2015; Au(I)–Cl⁻, Mei et al., 2014, Mei et al., 2017). Using distance-constrained thermodynamic integration, our previous studies successfully predicted the quantitative thermodynamic properties for the formation constants for Cu(I)–Cl⁻–HS⁻ (Mei et al., 2013a), Zn(II)–Cl⁻ (Mei et al., 2015b), Zn(II)–HS⁻ (Mei et al., 2016) and Pd(II)–Cl⁻–HS⁻ (Mei et al., 2015a). Two studies have successfully measured acid dissociation constants using ab initio MD simulations under hydrothermal conditions. Liu et al. (2013) obtained agreement to within 2 pKa units for molybdic acid to 300 °C, and Liu et al. (2015) obtained agreements of within 1 unit for pKa₁ of As(OH)_3(aq) to 300 °C, and for pKa_{1 to 3} of As(HS)_3(aq) at 25 °C (no experimental values exist at T ≥ 25 °C).

Despite all of the experimental and theoretical studies dedicated to HCl association, there is no ab initio MD study of HCl_(aq) dissociation at elevated temperature, and the mechanism of the increasing ion pairing/association of HCl_(aq) with increasing temperatures has not been clarified. There is also no data for HCl dissociation constants at ultra-high P-T relevant to mantle metasomatism and ultra-high pressure subduction environments (up to 60 kbar, 1000 °C). Under such conditions the role of HCl in controlling pH and metal complexation is unknown.

The present study aims to address these knowledge gaps in order to improve our understanding of HCl_(aq) dissociation/association, provide an independent check of the dissociation constants, and provide the first dataset for HCl dissociation at ultra-high P-T fluids. Hence, we conducted ab initio molecular dynamics (MD) simulations to investigate the association mechanism of HCl_(aq) as a function of temperature and pressure (25–700 °C, 1 bar to 60 kbar), and then calculated dissociation constants of HCl_(aq) based on distance-constraint thermodynamic integration. The dissociation constants obtained from ab initio MD simulations were then used to fit thermodynamic properties using the HKF model and DEW model (Sverjensky et al., 2014). These properties allow estimation of the role of HCl in controlling fluid pH and element mobility in deep earth fluids.