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

Authors: Yuan Mei; Weihua Liu; A. A. Migdiov; Joël Brugger; A. E. Williams-Jones
Journal: Geofluids
 (2018)

---

## Abstract

We investigated the hydration of the CuCl0 complex in HCl-bearing water vapor at 350°C and a vapor-like fluid density between 0.02 and 0.09 g/cm3 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 [H2O-Cu-Cl]0 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/cm3 (fH2O 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.


---

## Body

## 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 CuCl2- [6–10], and the activity of the chloride ion is the main control on the solubility of Cu minerals at constant temperature (T) and pressure (P). In contrast, in low-density vapors, there is increasing experimental evidence showing that, aside from the fugacity of HCl0, 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–/S3–, [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/cm3 [14], and the speciation and hydration of Au(I)-Cl complex in CO2-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/cm3 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 CuCl0 species. They reported formation constants for hydrated CuCl·nH2O0 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 (fH2O of 1.72–2.05, or fluid density 0.02–0.09 g/cm3; 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–28kg). The PBE exchange correlation-functional was used with a cut-off gradient correction of 5 × 10−5. 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 /cm3 in the simulations and the pressure (P) 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 H2O) was chosen to represent the CuCl0 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 CuCl0 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 CuCl0 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 H2O in the box, and the geometry of the first-shell CuCl(H2O) complex.
Densityg/cm3∗ Pressure(bar) Box size(Å) Simulation time(ps) d C u - C l  (Å) d C u - O  (Å) A n g l e C l - C u - O 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(H2O)0 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 H2O 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 [H2O-CuCl]0 complex.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(H2O)] clusters, second-shell hydration).
(a)
RDF of Cu and O (in all water) (b)
RDF of Cu and O (in solvent water)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(H2O)·mH2O 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.


## 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–28kg). The PBE exchange correlation-functional was used with a cut-off gradient correction of 5 × 10−5. 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 /cm3 in the simulations and the pressure (P) 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 H2O) was chosen to represent the CuCl0 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 CuCl0 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 CuCl0 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 H2O in the box, and the geometry of the first-shell CuCl(H2O) complex.
Densityg/cm3∗ Pressure(bar) Box size(Å) Simulation time(ps) d C u - C l  (Å) d C u - O  (Å) A n g l e C l - C u - O 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(H2O)0 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 H2O 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 [H2O-CuCl]0 complex.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(H2O)] clusters, second-shell hydration).
(a)
RDF of Cu and O (in all water) (b)
RDF of Cu and O (in solvent water)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(H2O)·mH2O 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(H2O)]0 (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/cm3 (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 CuCl2- complex for high-density (0.67–1.03 g/cm3) solutions and higher ionic strength (4 m Cl-) at 300°C.


### 3.2. Hydration of the [CuCl(H2O)]0 Complex

The radial distribution functions of Cu-O at different densities for all simulations times are plotted in Figure2 (left y-axes), together with their integrals representing the coordination number (right y-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 H2O 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(H2O)]0 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/cm3, to 3.8 at 0.03 g/cm3, and to 4 at 0.08 and 0.09 g/cm3. 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 CuCl0 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/cm3 [33].Table 2
Number of H2O molecules and H-bonds from MD simulations for different density and water fugacity at 350°C.
density f H 2 O log ⁡ f H 2 O n(H2O) within 5 Å of CuCl0 Total H-bonds n(H2O) in H-bonds cluster Total H-bondsper 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.25Figure 3
Dynamic hydration number of the CuCl0 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 (Figure4(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/cm3, 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.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(H2O)]0 cluster with O in the solvent water (excluding the water molecule in the [CuCl(H2O)]0 cluster).
(a)
RDF of H (in all water) and O (in all water) (b)
RDF of H (in binding water) and O (in solvent water)In the bulk water (Figure4(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/cm3; 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/cm3.As shown in Figure5(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 (H2O)m, water clusters held together by H-bonds. These clusters are interpreted to make up the second-shell hydration of the [CuCl(H2O)]0 complex and form the species [CuCl(H2O)·(H2O)m]0. The difference in the CN between the 0.02–0.025 g/cm3 density group and the 0.03–0.09 g/cm3 density groups reflects the decrease in hydration with decreasing density.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(H2O)]0 clusters) with O in the solvent water (excluding the water molecule in the [CuCl(H2O)]0 cluster).
(a)
RDF of O (in all water) and O (in all water) (b)
RDF of O (in binding water) and O (in solvent water)


### 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 Å (Figure6). 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 H2O, the CuCl0 complex can be written as [Cu(Cl(H2O)p)((H2O)(H2O)m))]0, where p is the number of water molecules surrounding Cl and m 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(H2O)]0 cluster).


### 3.5. Hydrogen Bonds

An explanation for the strong interaction between solvent water and the [CuCl(H2O)]0 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/cm3. 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(H2O)]0. 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(H2O)n]0 varies from 3 to 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/cm3. N(H2O) 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; N(H2O) in an H-bond cluster is the number of water molecules in the [CuCl(H2O)n]0 cluster linked by H-bonds; the combined N(H2O) value counts all the water molecules within 5 Å of Cu and water molecules beyond 5 Å that are linked by H-bonds.
(a)
6.05 ps; H-bonds_total = 7,N(H2O) within 5 Å of Cu = 2, N(H2O) in H-bond cluster = 3, combined N(H2O) = 3 (b)
16.85 ps; H-bonds_total = 13,N(H2O) within 5 Å of Cu = 5, N(H2O) in H-bond cluster = 5, combined N(H2O) = 8 (c)
19.69 ps; H-bonds_total = 11,N(H2O) within 5 Å of Cu = 4, N(H2O) in H-bond cluster = 4, combined N(H2O) = 4 (d)
23.63 ps; H-bonds_total = 18,N(H2O) within 5 Å of Cu = 2, N(H2O) in H-bond cluster = 8, combined N(H2O) = 8The 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(H2O)]0 cluster over the full time of the simulations are shown in Figures 8(a) and 8(b), respectively; note that the numbers of H2O in Figure 8(b) correspond to the total second-shell hydration (m+p) in the CuCl0 complex, [Cu(Cl(H2O)p)·((H2O)m+1)]0. 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/cm3 yielded smaller clusters than the simulations conducted at higher density (0.03–0.09 g/cm3). Overall, the m+p value is 1–1.5 at a density of 0.02–0.025 g/cm3, and 2.5–3.5 at a density of 0.03–0.09 g/cm3.Figure 8
Running average for the total number of H-bonds (a) and number of water molecules surrounding the [CuCl(H2O)]0 cluster linked by H-bonds (b) at different density.
(a)
Number of total H-bonds (b)
Number of water molecules surrounding CuCl(H2O) cluster linked by H-bondFigure9 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·nH2O0 (n = 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 H2O molecules within 5 Å of CuCl0 (green circles), average number of H2O 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·nH2O0 determined from solubility experiments by Migdisov et al. [3]. For each experimental solution, the average hydration number (Hy) was calculated using the equation Hy=∑i=1nCi∗i/Ctotal, where Ci is the concentration of the hydrated CuCl species with i water molecules,n is the total number of hydrated CuCl species, and Ctotal is the total concentration.


## 3.1. Overall Geometry of the First-Shell CuCl Complex

All MD simulations produced a linear first-shell complex, [CuCl(H2O)]0 (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/cm3 (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 CuCl2- complex for high-density (0.67–1.03 g/cm3) solutions and higher ionic strength (4 m Cl-) at 300°C.


## 3.2. Hydration of the [CuCl(H2O)]0 Complex

The radial distribution functions of Cu-O at different densities for all simulations times are plotted in Figure2 (left y-axes), together with their integrals representing the coordination number (right y-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 H2O 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(H2O)]0 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/cm3, to 3.8 at 0.03 g/cm3, and to 4 at 0.08 and 0.09 g/cm3. 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 CuCl0 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/cm3 [33].Table 2
Number of H2O molecules and H-bonds from MD simulations for different density and water fugacity at 350°C.
density f H 2 O log ⁡ f H 2 O n(H2O) within 5 Å of CuCl0 Total H-bonds n(H2O) in H-bonds cluster Total H-bondsper 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.25Figure 3
Dynamic hydration number of the CuCl0 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 (Figure4(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/cm3, 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.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(H2O)]0 cluster with O in the solvent water (excluding the water molecule in the [CuCl(H2O)]0 cluster).
(a)
RDF of H (in all water) and O (in all water) (b)
RDF of H (in binding water) and O (in solvent water)In the bulk water (Figure4(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/cm3; 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/cm3.As shown in Figure5(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 (H2O)m, water clusters held together by H-bonds. These clusters are interpreted to make up the second-shell hydration of the [CuCl(H2O)]0 complex and form the species [CuCl(H2O)·(H2O)m]0. The difference in the CN between the 0.02–0.025 g/cm3 density group and the 0.03–0.09 g/cm3 density groups reflects the decrease in hydration with decreasing density.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(H2O)]0 clusters) with O in the solvent water (excluding the water molecule in the [CuCl(H2O)]0 cluster).
(a)
RDF of O (in all water) and O (in all water) (b)
RDF of O (in binding water) and O (in solvent water)


## 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 Å (Figure6). 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 H2O, the CuCl0 complex can be written as [Cu(Cl(H2O)p)((H2O)(H2O)m))]0, where p is the number of water molecules surrounding Cl and m 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(H2O)]0 cluster).


## 3.5. Hydrogen Bonds

An explanation for the strong interaction between solvent water and the [CuCl(H2O)]0 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/cm3. 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(H2O)]0. 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(H2O)n]0 varies from 3 to 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/cm3. N(H2O) 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; N(H2O) in an H-bond cluster is the number of water molecules in the [CuCl(H2O)n]0 cluster linked by H-bonds; the combined N(H2O) value counts all the water molecules within 5 Å of Cu and water molecules beyond 5 Å that are linked by H-bonds.
(a)
6.05 ps; H-bonds_total = 7,N(H2O) within 5 Å of Cu = 2, N(H2O) in H-bond cluster = 3, combined N(H2O) = 3 (b)
16.85 ps; H-bonds_total = 13,N(H2O) within 5 Å of Cu = 5, N(H2O) in H-bond cluster = 5, combined N(H2O) = 8 (c)
19.69 ps; H-bonds_total = 11,N(H2O) within 5 Å of Cu = 4, N(H2O) in H-bond cluster = 4, combined N(H2O) = 4 (d)
23.63 ps; H-bonds_total = 18,N(H2O) within 5 Å of Cu = 2, N(H2O) in H-bond cluster = 8, combined N(H2O) = 8The 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(H2O)]0 cluster over the full time of the simulations are shown in Figures 8(a) and 8(b), respectively; note that the numbers of H2O in Figure 8(b) correspond to the total second-shell hydration (m+p) in the CuCl0 complex, [Cu(Cl(H2O)p)·((H2O)m+1)]0. 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/cm3 yielded smaller clusters than the simulations conducted at higher density (0.03–0.09 g/cm3). Overall, the m+p value is 1–1.5 at a density of 0.02–0.025 g/cm3, and 2.5–3.5 at a density of 0.03–0.09 g/cm3.Figure 8
Running average for the total number of H-bonds (a) and number of water molecules surrounding the [CuCl(H2O)]0 cluster linked by H-bonds (b) at different density.
(a)
Number of total H-bonds (b)
Number of water molecules surrounding CuCl(H2O) cluster linked by H-bondFigure9 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·nH2O0 (n = 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 H2O molecules within 5 Å of CuCl0 (green circles), average number of H2O 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·nH2O0 determined from solubility experiments by Migdisov et al. [3]. For each experimental solution, the average hydration number (Hy) was calculated using the equation Hy=∑i=1nCi∗i/Ctotal, where Ci is the concentration of the hydrated CuCl species with i water molecules,n is the total number of hydrated CuCl species, and Ctotal is the total concentration.


## 4. Discussion

The overarching objective of this study was to compare the hydration behavior of the CuCl0 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(H2O)]0 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 H2O 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/cm3).


### 4.1. Nature of Hydrated CuCl0 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 NaCl2- and Na2Cl2(aq). 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 CuCl0 in water vapor, our simulations reveal that the [H2O-CuCl]0 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/cm3 [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 [H2O-CuCl]0 moiety is, by nature, a statistical effect, and time-averaged hydration numbers only reflect the “averaged” nature of the interaction between the CuCl0 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 Figure8(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(H2O)n]0 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/cm3 over the 300–600°C temperature range [43], that is, slightly above the critical density of water (0.32 g/cm3). At fluids density of ~1 g/cm3, 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/cm3 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(H2O)·(H2O)m]0 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 [H2O-CuCl]0 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 [H2O-CuCl]0 cluster by H-bonds as defined by the criteria of Chialvo et al. [31] (distD-H < 3.5 Å and angleD-H-A < 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 (Figure9). The numbers of water molecules contributing to the hydration of the [H2O-CuCl]0 cluster calculated from the RDF and H-bond approach are similar, but the actual H2O molecules may differ, as illustrated in Figure 7. The reason for this is that the H-bond approach also considers the contribution of H2O molecules beyond 5 Å linked to the [H2O-CuCl]0 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(H2O)n]0 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 CuCl0 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·14H2O]0 does not possess a rigid stoichiometry and geometry and instead may represent an average hydration of state of CuCl0 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, H2O (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.


## 4.1. Nature of Hydrated CuCl0 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 NaCl2- and Na2Cl2(aq). 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 CuCl0 in water vapor, our simulations reveal that the [H2O-CuCl]0 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/cm3 [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 [H2O-CuCl]0 moiety is, by nature, a statistical effect, and time-averaged hydration numbers only reflect the “averaged” nature of the interaction between the CuCl0 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 Figure8(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(H2O)n]0 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/cm3 over the 300–600°C temperature range [43], that is, slightly above the critical density of water (0.32 g/cm3). At fluids density of ~1 g/cm3, 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/cm3 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(H2O)·(H2O)m]0 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 [H2O-CuCl]0 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 [H2O-CuCl]0 cluster by H-bonds as defined by the criteria of Chialvo et al. [31] (distD-H < 3.5 Å and angleD-H-A < 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 (Figure9). The numbers of water molecules contributing to the hydration of the [H2O-CuCl]0 cluster calculated from the RDF and H-bond approach are similar, but the actual H2O molecules may differ, as illustrated in Figure 7. The reason for this is that the H-bond approach also considers the contribution of H2O molecules beyond 5 Å linked to the [H2O-CuCl]0 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(H2O)n]0 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 CuCl0 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·14H2O]0 does not possess a rigid stoichiometry and geometry and instead may represent an average hydration of state of CuCl0 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, H2O (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-H2O vapor with density between 0.02 and 0.09 g/cm3 at 350°C, Cu(I) is bonded with one water and one chloride ligand as a linear [H2O-Cu-Cl]0 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.



---

*Source: 4279124-2018-05-02.xml*