Geofluids
Volume 2018 (2018), Article ID 4279124, 12 pages
https://doi.org/10.1155/2018/4279124

Research Article

CuCl Complexation in the Vapor Phase: Insights from Ab Initio Molecular Dynamics Simulations

Yuan Mei,1,2 Weihua Liu,¹ A. A. Migdiov,³ Joël Brugger,² and A. E. Williams-Jones⁴

¹CSIRO Mineral Resources, Clayton, VIC 3168, Australia
²School of Earth, Atmosphere and Environment, Monash University, VIC 3800, Australia
³Earth and Environmental Division, Los Alamos National Laboratory, M.S. J535, P.O. Box 1663, Los Alamos, NM 87545, USA
⁴Earth and Planetary Sciences, McGill University, 3450 University Street, Montreal, QC, Canada H3A 0E8

Correspondence should be addressed to Yuan Mei; yuan.mei@csiro.au

Received 1 December 2017; Accepted 21 March 2018; Published 2 May 2018

Academic Editor: Ferenc Molnar

Copyright © 2018 Yuan Mei 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.

Abstract.

We investigated the hydration of the CuCl⁰ complex in HCl-bearing water vapor at 350°C and a vapor-like fluid density between 0.02 and 0.09 g/cm³ using ab initio molecular dynamics (MD) simulations. The simulations reveal that one water molecule is strongly bonded to Cu(I) (first coordination shell), forming a linear [H₂O-Cu-Cl]⁰ moiety. The second hydration shell is highly dynamic in nature, and individual configurations have short life-spans in such low-density vapors, resulting in large fluctuations in instantaneous hydration numbers over a timescale of picoseconds. The average hydration number in the second shell (m) increased from ~0.5 to ~3.5 and the calculated number of hydrogen bonds per water molecule increased from 0.09 to 0.25 when fluid density (which is correlated to water activity) increased from 0.02 to 0.09 g/cm³ (H₂O 1.72 to 2.05). These changes of hydration number are qualitatively consistent with previous solubility studies under similar conditions, although the absolute hydration numbers from MD were much lower than the values inferred by correlating experimental Cu fugacity with water fugacity. This could be due to the uncertainties in the MD simulations and uncertainty in the estimation of the fugacity coefficients for these highly nonideal “vapors” in the experiments. Our study provides the first theoretical confirmation that beyond-first-shell hydrated metal complexes play an important role in metal transport in low-density hydrothermal fluids, even if it is highly disordered and dynamic in nature.

1. Introduction

Understanding the speciation of metals in hydrothermal fluids is important as it controls metal mobility and mineral solubility in geological and man-made hydrothermal systems, such as those forming ore deposits and geothermal systems and those used in hydrometallurgy [1]. Recent experimental and theoretical studies have shown that metals such as Cu and Mo can be transported in both liquid and vapor phases [2–4]. Copper is transported as chloride species in many hydrothermal fluids [5]. In high-density, reduced-sulfur-poor brines, the dominant copper species is Cu [6–10], and the activity of the chloride ion is the main control on the solubility of Cu minerals at constant temperature and pressure . In contrast, in low-density vapors, there is increasing experimental evidence showing that, aside from the fugacity of HCl⁰, the solubility of Cu minerals is controlled primarily by fluid density, which correlates with the fugacity of water. Specifically, mineral solubility increases exponentially with increasing fluid density in HCl-bearing water vapor [3, 11, 12].

Enabled by advances in computing power and algorithms, molecular dynamics simulations have been applied increasingly to the study of metal speciation in hydrothermal liquids (e.g., Cu(I)-Cl–, [13, 14]; Au(I)-HS–, [15]; Ag(I)-Cl–, [16, 17]; Au(I)-HS–/OH–/, [18, 19]; Cu(I)-HS–/Cl–, [20]; Au(I)-Cl–, [14]; Pd(II)-Cl–/HS–, [21]; Zn(II)-Cl–/HS–, [22, 23]; Ag(I)-HS–/OH–, [24]; Cu(I)-Ac–, [25]). These studies have demonstrated the importance of MD in understanding metal complexation at the molecular-level and, in some cases, have provided complementary information that helped resolve apparent discrepancies in the interpretation of experimental results [13, 22]. Molecular dynamics simulations have also been used to study metal transport in low-density aqueous fluids, namely, to determine the hydration behavior of Cu(I)-Cl and Au(I)-Cl species in fluids with a density as low as 0.1 g/cm³ [14], and the speciation and hydration of Au(I)-Cl complex in CO₂-rich vapors [26]. Concurrently, Migdisov et al. [3] measured the solubility of Cu(I) chloride in HCl-bearing water vapor between 350 and 550°C, with a fluid density down to 0.013 g/cm³ at 350°C. They found that the solubility of CuCl(s) increased exponentially as a function of water fugacity and proposed that this increase reflects an increasing number of hydration waters around the CuCl⁰ species. They reported formation constants for hydrated CuCl·nH₂O⁰ clusters containing one to fourteen water molecules (n = 1–14).

The aim of this study is to conduct MD simulations under conditions similar to those of the experiments of Migdisov et al. [3], especially at conditions with significant change of Cu solubility (fH₂O of 1.72–2.05, or fluid density 0.02–0.09 g/cm³; 350°C) to (1) verify the validity of the assumption that the increasing solubility of Cu with increasing water fugacity reflects a rapid change in hydration number; (2) to determine the nature and dynamics of the hydration sphere; and (3) to characterize the coordination structure of the dominant Cu species in HCl-bearing water vapor.

2. Method and Simulation Details

2.1. Ab initio Molecular Dynamics Simulations

Ab initio molecular dynamics simulations of CuCl complexes in the vapor phase were performed using the Car-Parrinello MD code “CPMD” version 3.17.1 [27]. Details of the calculation procedures are given in Mei et al. [14]; the key elements are repeated here. Vanderbilt ultrasoft pseudopotentials in the CPMD package were used with plane-wave cutoffs of 25 Ry (340.14 eV), a time-step of 4 a.u. (0.097 fs), and a fictitious electron mass of 400 a.u. (3.644 × 10–28 kg). The PBE exchange correlation-functional was used with a cut-off gradient correction of 5 × 10⁻⁵. Temperature was controlled by the Nosé thermostat for both the ions and electrons. All MD simulations were performed in the NVT ensemble at 350°C and run for >50 ps. The density of water was varied from 0.02 to 0.09 g /cm³ in the simulations and the pressure taken from the equations-of-state in the NIST database (Lemmon et al., 2000). The fugacity of water was calculated with HCh package [28] using water model of Kestin et al. [29].

The composition of the simulation box (1 Cu, 1 Cl, and 55 H₂O) was chosen to represent the CuCl⁰ complex in a vapor with low Cu and Cl concentrations. The formal concentrations of Cu and Cl are 1 molal, but since there were only one Cu and one Cl in the simulation box, there were no Cl-Cl or Cu-Cu interactions, which is consistent with the high dilution in vapors. In addition, the CuCl⁰ species (Figure 1) was stable over the full simulation time in all the simulations. Hence, the selected model box reflects the results of experimental studies well, as they [3, 11] have shown that a hydrated CuCl⁰ complex is the dominant species under the simulated conditions. The simulation details are listed in Table 1.

Table 1: Simulation details for a model at 350°C with 1 Cu, 1 Cl, and 55 H₂O in the box, and the geometry of the first-shell CuCl(H₂O) complex.

Density

Pressure
(bar)
Box size
(Å)
Simulation time
(ps)

(Å)

(Å)

0.02
60
43.494
53.1
2.05(7)
1.92(11)
162(11)

0.025
63
40.378
57.8
2.06(7)
1.92(9)
166(8)

0.03
75
37.996
56.1
2.07(7)
1.92(9)
165(8)

0.04
90
34.521
58.0
2.05(6)
1.91(8)
165(9)

0.06
120
30.157
57.3
2.07(6)
1.93(8)
165(9)

0.08
140
27.400
57.2
2.05(7)
1.91(9)
166(8)

0.09
153
26.248
58.0
2.06(7)
1.92(10)
165(7)

Density of water (Cu, Cl not included).

Figure 1: A snapshot of the CuCl(H₂O)⁰ complex and the surrounding water.

2.2. Hydration Number and H-Bonds

To obtain the time average for the ion-pairing and the hydration number, the radial distribution functions (RDF) and their integrals (related to coordination number, CN) for different atom pairs were calculated over the simulation time using VMD [30]. The overall hydration number and the effect of the water molecules were investigated by calculating the RDF of Cu/Cl/O/H with the O atoms in all water molecules. To highlight the effect of the solvent, the RDF of Cu/Cl/O/H with the O of the solvent water (excluding the H₂O molecules which bonded to Cu) was also calculated. The hydration number of Cu was calculated from the integrals of the RDF plots. Based on RDF plots showing a peak in Cu-O distances at ~3–5 Å (Figure 2) and previous MD studies of the hydration of CuCl in low-density fluids [14], a distance cut-off of 5.0 Å was chosen to represent the hydration shell of the [H₂O-CuCl]⁰ complex.

(a) RDF of Cu and O (in all water)

(b) RDF of Cu and O (in solvent water)

Figure 2: The Cu-O radial distribution function (RDF, solid lines, left axes) and the coordination number (CN, dashed lines, right axes) versus the Cu-O distance. (a) The RDF of Cu with all O atoms in the simulation. (b) The RDF of Cu with O atoms in the solvent water (excluding the water molecule in the [CuCl(H₂O)] clusters, second-shell hydration).

Another way to derive hydration numbers from the simulations is to use a geometrical definition of hydrogen bonding (H-bonds) to quantify the number of water molecules (m) in the CuCl(H₂O)·mH₂O clusters at each simulation step. Following the criteria developed by Chialvo et al. [31], a H-bond between the donor (D) and acceptor O or Cl atom (A) would be present, if the D-H distance is less than 3.5 Å and the D-H-A angle is less than 70°. We also quantified the average number of H-bonds among bulk water molecules in the simulation box using the same geometrical approach.

3. Results

3.1. Overall Geometry of the First-Shell CuCl Complex

All MD simulations produced a linear first-shell complex, [CuCl(H₂O)]⁰ (Figure 1). The CuCl distance (2.05–2.07 Å), Cu-O distance (1.91–1.93 Å), and O-CuCl angle (162–166°) at 350°C did not differ significantly over the density range of 0.02–0.09 g/cm³ (Table 1). The bond distance and angles are consistent with results of the in situ X-ray absorption spectroscopy (XAS) studies of Cu(I)-Cl and Cu(I)-HS complexes [6, 12, 32]. The CuCl distance determined in these simulations is also similar to that determined by Mei et al. [14] for the linear complex for high-density (0.67–1.03 g/cm³) solutions and higher ionic strength (4 m Cl⁻) at 300°C.

3.2. Hydration of the [CuCl(H₂O)]⁰ Complex

The radial distribution functions of Cu-O at different densities for all simulations times are plotted in Figure 2 (left -axes), together with their integrals representing the coordination number (right -axes). Figure 2(a) shows the RDF of Cu with all O in the simulation box. The first peaks in Figure 2(a) indicate Cu-O bond distances of 1.91–1.93 Å and correspond to the H₂O bonded directly to Cu(I). There are also weak peaks at ~3.5 Å that indicate the existence of a loosely bonded second water shell around the [CuCl(H₂O)]⁰ complex. This peak is highlighted in Figure 2(b), which shows the RDF of Cu with O calculated after excluding the first-shell water molecules. The hydration numbers of Cu at different densities are listed in Table 2. The number of hydration water molecules based on the integral of the RDF increases with density from 2.7 (1.7, excluding the first-shell oxygen) at a density of 0.02 g/cm³, to 3.8 at 0.03 g/cm³, and to 4 at 0.08 and 0.09 g/cm³. The evolution of the instantaneous hydration number at different values of fluid density is shown in Figure 3. The dynamic hydration number displays a large variation as a function of simulation time, indicating that the stability of the system is lower than is the case in high-density solutions [14]. In other words, the second-shell water molecules of CuCl⁰ in the vapor have short residence times; individual configurations have lifetimes ≤ 2 ps, similar, for example, to the lifetimes of sodium chloride ion pairs at 380°C, 0.55 g/cm³ [33].

Table 2: Number of H₂O molecules and H-bonds from MD simulations for different density and water fugacity at 350°C.

density

n(H₂O) within 5 Å of CuCl⁰
Total H-bonds
n(H₂O) in H-bonds
cluster
Total H-bonds
per water molecule

0.02
53
1.724
2.72
4.78
2.58
0.09

0.025
59
1.769
2.59
4.62
2.24
0.08

0.03
65
1.813
3.83
7.85
4.11
0.14

0.04
75
1.875
3.99
8.35
4.38
0.15

0.06
94
1.973
3.65
9.77
3.36
0.18

0.08
105
2.021
4.04
12.91
3.69
0.23

0.09
112
2.049
4.08
13.57
4.14
0.25

Figure 3: Dynamic hydration number of the CuCl⁰ complex (within 5 Å of Cu) as a function of simulation time at different values of density (the density values are indicated at the top right in each diagram).

3.3. Water-Water Interaction

Here we compare the environment of solvent water molecules (bulk water) to that of the water bound to Cu(I) (bound water). To assess the bulk water, we calculated the RDF of all H with O in the solution over the complete simulation (Figure 4(a)). This RDF shows a very sharp peak at 1 Å, indicating the strong H-O bond within water molecules. The two bands at ~2.0 and 3.4 Å show the interaction between H and O (H-bond) in the solution. In contrast, the RDF plot (Figure 4(b)) of H in the bound water with O in the solvent water shows a peak at ~1.75 Å, and the integral gives 0.3, 0.6, 0.6, and 0.5 oxygen atoms surrounding each hydrogen atom at a distance of 2.5 Å (first minimum in RDF, Figure 4(b)) at a density of 0.02, 0.04, 0.06, and 0.08 g/cm³, respectively. The second peak with a maximum at ~3.2 Å in Figure 4(b) is interpreted to reflect the interaction between O in the solvent water with H in the bond water.

(a) RDF of H (in all water) and O (in all water)

(b) RDF of H (in binding water) and O (in solvent water)

Figure 4: H-O radial distribution function (RDF, solid lines, left axes) and the coordination number (CN, dashed lines, right axes) versus the H-O distance. (a) shows the RDF of all H with all O atoms in the simulation; (b) shows the RDF of H in the [CuCl(H₂O)]⁰ cluster with O in the solvent water (excluding the water molecule in the [CuCl(H₂O)]⁰ cluster).

In the bulk water (Figure 4(a)), the coordination number of H shows systematic changes as a function of solution density. In contrast, the coordination of H belonging to O bound to Cu(I) reflects two groups of solutions with very similar coordination numbers (Figure 4(b)): 0.5–0.6 O at 2.5 Å and ~2 O at 4 Å for simulations at a density of 0.03–0.09 g/cm³; and 0.3 O at 2.5 Å and 1 O at 4 Å for the second group at a density of 0.02–0.025 g/cm³.

As shown in Figure 5(a), the bulk O-O RDF produced a peak at 3 Å, whereas the distribution of O in the bound water with O in the solvent water (Figure 5(b)) produced stronger peaks at the lesser distance of 2.7 Å. The strong interactions between the bound water and the solvent water molecules shown in Figures 4(b) and 5(b) indicate the formation of (H₂O)m, water clusters held together by H-bonds. These clusters are interpreted to make up the second-shell hydration of the [CuCl(H₂O)]⁰ complex and form the species [CuCl(H₂O)·(H₂O)m]⁰. The difference in the CN between the 0.02–0.025 g/cm³ density group and the 0.03–0.09 g/cm³ density groups reflects the decrease in hydration with decreasing density.

(a) RDF of O (in all water) and O (in all water)

(b) RDF of O (in binding water) and O (in solvent water)

Figure 5: O-O radial distribution function (RDF, solid lines, left axes) and the coordination number (CN, dash lines, right axes) versus the O-O distance. (a) All O in the simulation box; (b) O in the binding water (in [CuCl(H₂O)]⁰ clusters) with O in the solvent water (excluding the water molecule in the [CuCl(H₂O)]⁰ cluster).

3.4. Chloride Hydration

The distribution of solvent water O around the Cl bound to Cu(I) is related to peaks at 3.5–3.6 Å (Figure 6). The intensity of these RDF peaks is smaller than that for the O bound to solvent water O (Figure 5(b)), indicating weaker interaction between Cl and solvent water. This suggests that the water cluster surrounding the bound water contributes more to Cu hydration than the water cluster surrounding chloride. If we consider water clusters as acting as ligands in a manner similar to first-shell Cl– and H₂O, the CuCl⁰ complex can be written as [Cu(Cl(H₂O)p)((H₂O)(H₂O)m))]⁰, where is the number of water molecules surrounding Cl and is the number of water molecules surrounding the bound water in the second shell.

Figure 6: Radial distribution function (RDF, solid lines, left axes) and the coordination number (CN, dash lines, right axes) as a function of the Cl-O distance of Cl versus O in solvent water (excluded water molecule in [CuCl(H₂O)]⁰ cluster).

3.5. Hydrogen Bonds

An explanation for the strong interaction between solvent water and the [CuCl(H₂O)]⁰ complex is that water molecules form clusters that attach to metal complexes with the aid of strong hydrogen bonds (H-bonds) [34]. The number of H-bonds per water molecule in the simulation boxes estimated from the bond geometry is listed in Table 2. Figure 7 shows four snapshots of the simulation for a density of 0.08 g/cm³. There H atoms play the roles of donor (D) and Cl or O plays the role of acceptor (A) in the H-bonds. Similar hydrogen bonding structure and the hydration of Cl– in hydrothermal NaCl solutions were also quantitatively addressed by experimental (IR and Raman) and computational (classical MD) approaches [35]. They found that the feature of the NaCl structure has a considerable contribution to strong O-H⋯Cl– bonds, which decreases with increasing temperature due to the ion-pairing. The snapshots taken for different simulation times show clearly the highly dynamic nature of the water clusters bonded to [CuCl(H₂O)]⁰. It is evident in the variable number of water molecules linked by H-bonds (3–8) and the variable number of water within 5 Å of the complex (2–5). The combined number of water molecules around [CuCl(H₂O)n]⁰ varies from 3 to 8.

(a) 6.05 ps; H-bonds_total = 7, (H₂O) within 5 Å of Cu = 2, (H₂O) in H-bond cluster = 3, combined (H₂O) = 3

(b) 16.85 ps; H-bonds_total = 13, (H₂O) within 5 Å of Cu = 5, (H₂O) in H-bond cluster = 5, combined (H₂O) = 8

(c) 19.69 ps; H-bonds_total = 11, (H₂O) within 5 Å of Cu = 4, (H₂O) in H-bond cluster = 4, combined (H₂O) = 4

(d) 23.63 ps; H-bonds_total = 18, (H₂O) within 5 Å of Cu = 2, (H₂O) in H-bond cluster = 8, combined (H₂O) = 8

Figure 7: Hydrogen bonds (red dashed lines) and the sphere of O within 5 Å of Cu (shown by the grey dashed lines) in the simulation for a density of 0.08 g/cm³. (H₂O) within 5 Å of Cu is the hydration number calculated from the RDF and its integral; H-bonds were picked up when the D-H distance was less than 3.5 Å and the D-H-A angle less than 70°; H-bonds_total is the total number of H-bonds in the solution; (H₂O) in an H-bond cluster is the number of water molecules in the [CuCl(H₂O)n]⁰ cluster linked by H-bonds; the combined (H₂O) value counts all the water molecules within 5 Å of Cu and water molecules beyond 5 Å that are linked by H-bonds.

The running averages of the total number of H-bonds in the simulation box and the number of H-bond-linked water molecules surrounding the [CuCl(H₂O)]⁰ cluster over the full time of the simulations are shown in Figures 8(a) and 8(b), respectively; note that the numbers of H₂O in Figure 8(b) correspond to the total second-shell hydration in the CuCl⁰ complex, [Cu(Cl(H₂O)p)·((H₂O)m+1)]⁰. Figures 8(a) and 8(b) display strong oscillations in the size of the clusters at the time scale of 10–15 ps, indicating that the H-bonded clusters have lifetimes in the order of 5–10 ps, as observed by previous studies for H-bonds (the lifetime of H-bonds in supercritical water is 0.2–0.5 ps [36, 37]; see reviews by [38]). The total number of H-bonds increases as a function of density, consistent with the trends observed in Figures 4(a) and 5(a). Considering the water clusters linked by H-bonds, the simulations at a density of 0.02–0.025 g/cm³ yielded smaller clusters than the simulations conducted at higher density (0.03–0.09 g/cm³). Overall, the value is 1–1.5 at a density of 0.02–0.025 g/cm³, and 2.5–3.5 at a density of 0.03–0.09 g/cm³.

(a) Number of total H-bonds

(b) Number of water molecules surrounding CuCl(H2O) cluster linked by H-bond

Figure 8: Running average for the total number of H-bonds (a) and number of water molecules surrounding the [CuCl(H₂O)]⁰ cluster linked by H-bonds (b) at different density.

Figure 9 compares the hydration number derived from the MD simulations (RDF and geometrical H-bonds) with the average hydration number calculated from the thermodynamic properties for hydrated CuCl·nH₂O⁰ ( = 1–14) species derived by Migdisov et al. [3] from their solubility data. It can be seen that the hydration number derived from the MD simulation is significantly lower than that estimated from the data of Migdisov et al. [3]. The reasons for this will be discussed in the next section.

Figure 9: Average number of H₂O molecules within 5 Å of CuCl⁰ (green circles), average number of H₂O in H-bonds clusters (pink diamonds), and average hydration numbers based on the solubility study (blue hatches) of Migdisov et al. [3]. Also shown is the hydration number calculated using the thermodynamic properties of hydrated CuCl·H₂O⁰ determined from solubility experiments by Migdisov et al. [3]. For each experimental solution, the average hydration number (Hy) was calculated using the equation , where is the concentration of the hydrated CuCl species with water molecules, n is the total number of hydrated CuCl species, and is the total concentration.

4. Discussion

The overarching objective of this study was to compare the hydration behavior of the CuCl⁰ complex derived from the MD with that determined experimentally by Migdisov et al. [3] for a similar temperature (350°C) and fluid density and to explain any differences in this behavior predicted by the two approaches. Hydration numbers for the [CuCl(H₂O)]⁰ moiety were derived from MD simulations by calculating the integral of the Cu-O paired distribution function and by counting the H-bonded water molecules using a geometrical definition of H-bonds (Table 2). As illustrated in Figure 9, the simulations confirm a trend of increasing hydration with increasing solution density, which tallies with the solubility measurements. However, the hydration numbers predicted using both methods are much lower in the simulations than those predicted by the experimentally based thermodynamic model.

Significantly, the MD simulations reveal that the H₂O clusters are short-lived (a few picoseconds) and the second hydration shell is very dynamic. As a result, the hydration number varies greatly over short periods of time (Figures 3 and 7). This indicates that the hydration water in the second shell is only weakly bonded, even compared to similar simulations at higher density for CuCl and Au-Cl complexes [14, 26]. This raises an important question about the nature of the thermodynamic models that we use to model such solutions with density between those of (near) ideal gases and those of condensed solutions (density  >  ~0.6 g/cm³).

4.1. Nature of Hydrated CuCl⁰ Species

Speciation-based thermodynamic models, as opposed, for example, to models that assume full dissociation of complexes (e.g., Pitzer, 1973), allow for successful extrapolations beyond experimental conditions because they include much important intrinsic information about the chemistry of the solutions [39]. The concept of “speciation,” however, can become ambiguous in some situations, in which case, it is necessary to reconsider the significance and validity of the assumptions linking physical reality and thermodynamic models. For example, Driesner et al. [33] illustrated the formation of a large variety of NaCl clusters in concentrated solutions. In such solutions, the speciation was defined by the time-averaged proportions of clusters such as Na and Na₂. Owing to the short residence times and variable geometry of different NaCl clusters, it is difficult to derive a speciation model with a precise geometry and stoichiometry for accurate extrapolation.

In the case of CuCl⁰ in water vapor, our simulations reveal that the [H₂O-CuCl]⁰ moiety is very stable (residence time > 50 ps) with a well-defined linear structure, confirming the information obtained via solubility, in situ UV-Vis, and in situ XAS experiments at a density > 0.19 g/cm³ [6–9, 12, 40–42]. In contrast, the second hydration shell revealed in this study is highly disordered and dynamic (Figures 3 and 7). The hydration of the [H₂O-CuCl]⁰ moiety is, by nature, a statistical effect, and time-averaged hydration numbers only reflect the “averaged” nature of the interaction between the CuCl⁰ cluster and the solvent molecules within certain time scales.

4.2. Hydrogen Bonding: Vapor versus High-Density Liquid

To evaluate the effect of water-water interaction, we counted the total number of H-bonds in the solutions with different fluid density (or water fugacity). As shown in Figure 8(a), the running average of the number of H-bonds shows good convergence over the whole simulation time and increases with increasing fluid density (Table 2). The short lifetime of the [CuCl(H₂O)n]⁰ clusters in low-density fluids reflects the weaker interaction between water molecules in the solvent at low-density compared to high-density, because changing configurations in low-density fluids requires breaking fewer H-bonds and hence a lower activation energy in low-density than in high-density fluids [43].

Following the approach of Kalinichev [43], the total H-bonds in the simulation box were expressed as “H-bond total per water molecule” by dividing total H-bonds with 55, the total number of water molecules in the simulation box (Table 2). Kalinichev [43] (using Monte Carlo and classical MD data presented in [38]) calculated the number of H-bonds per molecule in bulk water for a wide range in density over a temperature range of 25–1000°C. Their calculations showed a near-linear increase in the number of H-bonds as a function of density for bulk water, and a small decrease in their number with increasing temperature for a given density. The degree of H-bonding in the bulk fluid is an important parameter to consider when using a geometrical approach to define the hydration of a cluster; in high-density solutions, all water molecules interact via H-bonding, so that hydration defined in this manner becomes infinite. This happens when the number of H-bond per water molecule ~1 at a density of ~0.4 g/cm³ over the 300–600°C temperature range [43], that is, slightly above the critical density of water (0.32 g/cm³). At fluids density of ~1 g/cm³, the theoretical percolation threshold for the network of hydrogen bonds in water is known to be ~1.4–1.7, which has since been confirmed by multiple molecular simulations [44–46] and experiments [47, 48]. Our ab initio calculations at very low-density range (0.02–0.09 g/cm³ at 350°C) indicated 0.09–0.25 H-bonds per water molecule. These values are consistent with the trends reported by Kalinichev [43].

4.3. Defining Microscopic Hydration Numbers and Comparison with Experiments

Owing to the short lifetime and lack of definite geometry for second-shell complexes in low-density vapor, defining individual hydrated complexes [CuCl(H₂O)·(H₂O)m]⁰ is not a trivial task and will depend to some degree on the definition of what constitutes a cluster of particular stoichiometry (m). In this study, we used two approaches to quantify the hydration of the [H₂O-CuCl]⁰ cluster (Figure 1). (i) In the first approach (RDF approach), the hydration shell is obtained by calculating the number of water molecules in a defined sphere of fixed radius (5 Å) around the Cu(I) atom, based on the RDF. (ii) The second approach (H-bond approach) is to count the number of water molecules in the whole simulation box that are linked to the [H₂O-CuCl]⁰ cluster by H-bonds as defined by the criteria of Chialvo et al. [31] ( < 3.5 Å and < 70°).

Overall, these two approaches produce very similar trends and absolute values for the change of hydration as function of fluids density; however, both yielded smaller values than the numbers inferred from the experimental data (Figure 9). The numbers of water molecules contributing to the hydration of the [H₂O-CuCl]⁰ cluster calculated from the RDF and H-bond approach are similar, but the actual H₂O molecules may differ, as illustrated in Figure 7. The reason for this is that the H-bond approach also considers the contribution of H₂O molecules beyond 5 Å linked to the [H₂O-CuCl]⁰ cluster by H-bonds (Figures 7(a) and 7(d)). However, in some cases, the H-bond approach ignores the potential contribution of water molecules that are bonded to the Cu(I) ion via Van der Waals forces. For example, in Figure 7(b), some water molecules are linked to Cu by Van der Waals forces (grey dashed lines); these water molecules may be part of larger hydrated cluster and not considered by the H-bond approach. Consequently, both approaches provide a lower estimate of the actual hydration number (see four examples in Figure 7). Note that the combined number, n, of water molecules around [CuCl(H₂O)n]⁰ in snapshots can total 8 (Figure 7), and this value is also very dynamic.

Although our analysis of MD data underestimates the hydration number of the CuCl complex predicted by the experimental solubility study [3], the MD results, nonetheless, show the same trend as those estimated from the experimental study. There are two likely reasons for the discrepancy between the hydration numbers derived from our ab initio MD simulations and the solubility experiments. One of these reasons, given the nature of ab initio MD (limited in simulation time and number of molecules), is that our calculations may have underestimated the formation of large, fairly stable clusters. The limitation in simulation size of ab initio MD means that our model does not reproduce the spatial heterogeneity in density, and thus it is difficult to evaluate all possible long-range effects arising from interactions of the CuCl⁰ complex with solvent molecules having a large radius. For example, in the solubility study, the results of experiments at the highest density predicted a hydration number of 14 [3], whereas the MD analysis predicted a number of 4. This could be because a species such as [CuCl·14H₂O]⁰ does not possess a rigid stoichiometry and geometry and instead may represent an average hydration of state of CuCl⁰ in a solution with molecules having a large radius (i.e., outside the hydration sphere of 5 Å calculated by RDF approach, Figure 7). This is undefined in MD analysis, as it represents a solution that is intermediate between a gas (mostly noninteracting molecules) and a solution (mostly interacting water molecules). Mixed classical and ab initio MD simulations, ground-proofed using the current work, may help to provide further insight into the nature of hydration and the long-range effect of the solvent in low-density fluids at these intermediate conditions. Another reason for the discrepancy between the hydration number predicted by our simulations and that predicted from solubility experiments is that the latter is sensitive to the choice of fugacity model used. In treating the experimental solubility data, Migdisov et al. [3] employed the Lewis-Randall model for gas mixtures (ideal mixture of real gases), in which the behavior of Cu-bearing clusters was assumed to obey the Ideal Gas Law. The authors assumed that most of the nonideality in the behavior of these species can be addressed by accounting for the formation of hydrogen bonds between gaseous CuCl and gas-solvent, H₂O (cluster formation). The validity of this assumption, however, has not been verified, and we do not exclude the possibility that, in addition to hydration, other processes may contribute to the deviation of the formed moieties from ideal behavior.

5. Conclusions

The ab initio molecular dynamics (MD) simulations show that in a HCl-H₂O vapor with density between 0.02 and 0.09 g/cm³ at 350°C, Cu(I) is bonded with one water and one chloride ligand as a linear [H₂O-Cu-Cl]⁰ moiety, with a highly dynamic second-shell hydration shell that has short life-span configurations. The changes of hydration number as a function of water fugacity are qualitatively consistent with previous solubility under similar conditions, with lower hydration numbers than the experimentally inferred values [3]. The discrepancy could be due to (i) our MD analysis underestimating the “thermodynamic” hydration number (e.g., ignoring the effects of density heterogeneity) in high temperature fluids; and (ii) uncertainty in the estimation of the fugacity coefficients for these highly nonideal “vapors” in the experiments. The calculated number of hydrogen bond per water molecule decreases with decreasing water density, indicating weak interaction between water molecules in such low-density fluids. This study provides the first theoretical confirmation that beyond-first-shell hydrated metal complexes play an important role in metal transport in low-density hydrothermal fluids, even if it is highly disordered and dynamic in nature.

Conflicts of Interest

There are no conflicts of interest regarding the publication of this study.

Acknowledgments

Research funding was provided by the CSIRO Office of the Chief Executive (OCE) fellowship to Yuan Mei and the Australian Research Council (ARC) to Weihua Liu (FT130100510). The MD calculations in this study were supported by resources provided by the Pawsey Supercomputing Centre with funding from the Australian Government and the Government of Western Australia and the high-performance computers at CSIRO.

References

J. Brugger, W. Liu, B. Etschmann, Y. Mei, D. M. Sherman, and D. Testemale, “A review of the coordination chemistry of hydrothermal systems, or do coordination changes make ore deposits?” Chemical Geology, vol. 447, pp. 219–253, 2016.

N. C. Hurtig and A. E. Williams-Jones, “An experimental study of the transport of gold through hydration of AuCl in aqueous vapour and vapour-like fluids,” Geochimica et Cosmochimica Acta, vol. 127, pp. 305–325, 2014.

A. A. Migdisov, A. Y. Bychkov, A. E. Williams-Jones, and V. J. Van Hinsberg, “A predictive model for the transport of copper by HCl-bearing water vapour in ore-forming magmatic-hydrothermal systems: Implications for copper porphyry ore formation,” Geochimica et Cosmochimica Acta, vol. 129, pp. 33–53, 2014.

A. E. Williams-Jones and C. A. Heinrich, “100th Anniversary Special Paper: Vapor transport of metals and the formation of magmatic-hydrothermal ore deposits,” Economic Geology, vol. 100, no. 7, pp. 1287–1312, 2005.

T. M. Seward and H. L. Barnes, “Metal transport by hydrothermal ore fluids,” in Geochemistry of Hydrothermal Ore Deposits, H. L. Barnes, Ed., vol. 3rd, pp. 435–486, Wiley, NY, USA, 1997.

J. Brugger, B. Etschmann, W. Liu et al., “An XAS study of the structure and thermodynamics of Cu(I) chloride complexes in brines up to high temperature (400°C, 600 bar),” Geochimica et Cosmochimica Acta, vol. 71, no. 20, pp. 4920–4941, 2007.

J. L. Fulton, M. M. Hoffmann, and J. G. Darab, “An X-ray absorption fine structure study of copper(I) chloride coordination structure in water up to 325°C,” Chemical Physics Letters, vol. 330, no. 3-4, pp. 300–308, 2000.

W. Liu, D. C. McPhail, and J. Brugger, “An experimental study of copper(I)-chloride and copper(I)-acetate complexing in hydrothermal solutions between 50°C and 250° and vapor-saturated pressure,” Geochimica et Cosmochimica Acta, vol. 65, no. 17, pp. 2937–2948, 2001.

M. Louvel, A. Bordage, B. Tripoli, D. Testemale, J.-L. Hazemann, and J. Mavrogenes, “Effect of S on the aqueous and gaseous transport of Cu in porphyry and epithermal systems: Constraints from in situ XAS measurements up to 600°C and 300 bars,” Chemical Geology, vol. 466, pp. 500–511, 2017.

R. Zhong, J. Brugger, Y. Chen, and W. Li, “Contrasting regimes of Cu, Zn and Pb transport in ore-forming hydrothermal fluids,” Chemical Geology, vol. 395, pp. 154–164, 2015.

S. M. Archibald, A. A. Migdisov, and A. E. Williams-Jones, “An experimental study of the stability of copper chloride complexes in water vapor at elevated temperatures and pressures,” Geochimica et Cosmochimica Acta, vol. 66, no. 9, pp. 1611–1619, 2002.

W. Liu, J. Brugger, B. Etschmann, D. Testemale, and J.-L. Hazemann, “The solubility of nantokite (CuCl(s)) and Cu speciation in low-density fluids near the critical isochore: An in-situ XAS study,” Geochimica et Cosmochimica Acta, vol. 72, no. 16, pp. 4094–4106, 2008.

D. M. Sherman, “Complexation of Cu⁺ in Hydrothermal NaCl Brines: Ab initio molecular dynamics and energetics,” Geochimica et Cosmochimica Acta, vol. 71, no. 3, pp. 714–722, 2007.

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-5000bar),” Geochimica et Cosmochimica Acta, vol. 131, pp. 196–212, 2014.

X. Liu, X. Lu, R. Wang, H. Zhou, and S. Xu, “Speciation of gold in hydrosulphide-rich ore-forming fluids: Insights from first-principles molecular dynamics simulations,” Geochimica et Cosmochimica Acta, vol. 75, no. 1, pp. 185–194, 2011.

X. Liu, X. Lu, R. Wang, and H. Zhou, “Silver speciation in chloride-containing hydrothermal solutions from first principles molecular dynamics simulations,” Chemical Geology, vol. 294-295, pp. 103–112, 2012.

G. S. Pokrovski, J. Roux, G. Ferlat et al., “Silver in geological fluids from in situ X-ray absorption spectroscopy and first-principles molecular dynamics,” Geochimica et Cosmochimica Acta, vol. 106, pp. 501–523, 2013.

Y. Mei, D. M. Sherman, W. Liu, and J. Brugger, “Complexation of gold in S3--rich hydrothermal fluids: Evidence from ab-initio molecular dynamics simulations,” Chemical Geology, vol. 347, pp. 34–42, 2013.

G. S. Pokrovski, M. A. Kokh, D. Guillaume et al., “Sulfur radical species form gold deposits on Earth,” Proceedings of the National Acadamy of Sciences of the United States of America, vol. 112, no. 44, pp. 13484–13489, 2015.

Y. Mei, D. M. Sherman, W. Liu, and J. Brugger, “Ab initio molecular dynamics simulation and free energy exploration of copper(I) complexation by chloride and bisulfide in hydrothermal fluids,” Geochimica et Cosmochimica Acta, vol. 102, pp. 45–64, 2013.

Y. Mei, B. Etschmann, W. Liu et al., “Palladium complexation in chloride- and bisulfide-rich fluids: Insights from ab initio molecular dynamics simulations and X-ray absorption spectroscopy,” Geochimica et Cosmochimica Acta, vol. 161, pp. 128–145, 2015.

Y. Mei, D. M. Sherman, W. Liu, B. Etschmann, D. Testemale, and J. Brugger, “Zinc complexation in chloride-rich hydrothermal fluids (25-600°C): A thermodynamic model derived from ab initio molecular dynamics,” Geochimica et Cosmochimica Acta, vol. 150, pp. 265–284, 2015.

Y. Mei, B. Etschmann, W. Liu, D. M. Sherman, D. Testemale, and J. Brugger, “Speciation and thermodynamic properties of zinc in sulfur-rich hydrothermal fluids: Insights from ab initio molecular dynamics simulations and X-ray absorption spectroscopy,” Geochimica et Cosmochimica Acta, vol. 179, pp. 32–52, 2016.

M. He, X. Liu, X. Lu, C. Zhang, and R. Wang, “Structures and Acidity Constants of Silver-Sulfide Complexes in Hydrothermal Fluids: A First-Principles Molecular Dynamics Study,” The Journal of Physical Chemistry A, vol. 120, no. 42, pp. 8435–8443, 2016.

F. Lai, L. Liu, and W. Cao, “Complexation of copper in acetate-rich low-temperature hydrothermal fluids: Evidence from ab initio molecular dynamics simulations,” Chemical Geology, vol. 476, pp. 100–118, 2018.

Y. Mei, W. Liu, J. Brugger, A. A. Migdisov, and A. E. Williams-Jones, “Hydration Is the Key for Gold Transport in CO₂–HCl–H₂O Vapor,” ACS Earth and Space Chemistry, vol. 1, no. 7, pp. 368–375, 2017.

R. Car and M. Parrinello, “Unified approach for molecular dynamics and density-functional theory,” Physical Review Letters, vol. 55, no. 22, pp. 2471–2474, 1985.

Y. V. Shvarov, “HCh: New potentialities for the thermodynamic simulation of geochemical systems offered by windows,” Geochemistry International, vol. 46, no. 8, pp. 834–839, 2008.

J. Kestin, J. V. Sengers, B. Kamgar-Parsi, and J. M. H. Levelt Sengers, “Thermophysical properties of fluid H₂O,” Journal of Physical and Chemical Reference Data, vol. 13, no. 1, pp. 175–183, 1984.

W. Humphrey, A. Dalke, and K. Schulten, “VMD: visual molecular dynamics,” Journal of Molecular Graphics, vol. 14, no. 1, pp. 33–38, 1996.

A. A. Chialvo, P. T. Cummings, and H. D. Cochran, “Solvation structure, hydrogen bonding, and ion pairing in dilute supercritical aqueous NaCl mixtures,” International Journal of Thermophysics, vol. 17, no. 1, pp. 147–156, 1996.

B. E. Etschmann, W. Liu, D. Testemale et al., “An in situ XAS study of copper(I) transport as hydrosulfide complexes in hydrothermal solutions (25-592°C, 180-600bar): Speciation and solubility in vapor and liquid phases,” Geochimica et Cosmochimica Acta, vol. 74, no. 16, pp. 4723–4739, 2010.

T. Driesner, T. M. Seward, and I. G. Tironi, “Molecular dynamics simulation study of ionic hydration and ion association in dilute and 1 molal aqueous sodium chloride solutions from ambient to supercritical conditions,” Geochimica et Cosmochimica Acta, vol. 62, no. 18, pp. 3095–3107, 1998.

D. C. Lacy, J. Mukherjee, R. L. Lucas, V. W. Day, and A. S. Borovik, “Metal complexes with varying intramolecular hydrogen bonding networks,” Polyhedron, vol. 52, pp. 261–267, 2013.

G. V. Bondarenko, Y. E. Gorbaty, A. V. Okhulkov, and A. G. Kalinichev, “Structure and hydrogen bonding in liquid and supercritical aqueous NaCl solutions at a pressure of 1000 bar and temperatures up to 500°C: A comprehensive experimental and computational study,” The Journal of Physical Chemistry A, vol. 110, no. 11, pp. 4042–4052, 2006.

R. D. Mountain, “Comparison of a fixed-charge and a polarizable water model,” The Journal of Chemical Physics, vol. 103, no. 8, pp. 3084–3090, 1995.

T. I. Mizan, P. E. Savage, and R. M. Ziff, “Temperature dependence of hydrogen bonding in supercritical water,” The Journal of Physical Chemistry C, vol. 100, no. 1, pp. 403–408, 1996.

A. G. Kalinichev, “Molecular simulations of liquid and supercritical water: Thermodynamics, structure, and hydrogen bonding,” Reviews in Mineralogy and Geochemistry, vol. 42, pp. 83–129, 2001.

T. M. Seward and T. Driesner, “Hydrothermal solution structure: experiments and computer simulations,” in Aqueous Systems at Elevated Temperatures and Pressures: Physical Chemistry in Water, Steam and Hydrothermal Solutions, D. A. Palmer, R. Fernández-Prini, and A. H. Havey, Eds., pp. 149–182, Elsevier, 2004.

W. Liu, J. Brugger, D. C. McPhail, and L. Spiccia, “A spectrophotometric study of aqueous copper(I)-chloride complexes in LiCl solutions between 100°C and 250°C,” Geochimica et Cosmochimica Acta, vol. 66, no. 20, pp. 3615–3633, 2002.

W. Liu and D. C. McPhail, “Thermodynamic properties of copper chloride complexes and copper transport in magmatic-hydrothermal solutions,” Chemical Geology, vol. 221, no. 1-2, pp. 21–39, 2005.

Z. Xiao, C. H. Gammons, and A. E. Williams-Jones, “Experimental study of copper(I) chloride complexing in hydrothermal solutions at 40 to 300°C and saturated water vapor pressure,” Geochimica et Cosmochimica Acta, vol. 62, no. 17, pp. 2949–2964, 1998.

A. G. Kalinichev, “Universality of hydrogen bond distributions in liquid and supercritical water,” Journal of Molecular Liquids, vol. 241, pp. 1038–1043, 2017.

H. E. Stanley and J. Teixeira, “Interpretation of the unusual behavior of H₂O and D₂O at low temperatures: Tests of a percolation model,” The Journal of Chemical Physics, vol. 73, no. 7, pp. 3404–3422, 1980.

A. G. Kalinichev and S. V. Churakov, “Thermodynamics and structure of molecular clusters in supercritical water,” Fluid Phase Equilibria, vol. 183-184, pp. 271–278, 2001.

L. B. Pártay, P. Jedlovszky, I. Brovchenko, and A. Oleinikova, “Percolation transition in supercritical water: A Monte Carlo simulation study,” The Journal of Physical Chemistry B, vol. 111, no. 26, pp. 7603–7609, 2007.

M. Bernabei, A. Botti, F. Bruni, M. A. Ricci, and A. K. Soper, “Percolation and three-dimensional structure of supercritical water,” Physical Review E: Statistical, Nonlinear, and Soft Matter Physics, vol. 78, no. 2, Article ID 021505, 2008.

Y. Gorbaty and G. V. Bondarenko, “Transition of liquid water to the supercritical state,” Journal of Molecular Liquids, vol. 239, pp. 5–9, 2017.