Hydrogen‐Bond Structure and Low‐Frequency Dynamics of Electrolyte Solutions: Hydration Numbers from ab Initio Water Reorientation Dynamics and Dielectric Relaxation Spectroscopy

Abstract We present an atomistic simulation scheme for the determination of the hydration number (h) of aqueous electrolyte solutions based on the calculation of the water dipole reorientation dynamics. In this methodology, the time evolution of an aqueous electrolyte solution generated from ab initio molecular dynamics simulations is used to compute the reorientation time of different water subpopulations. The value of h is determined by considering whether the reorientation time of the water subpopulations is retarded with respect to bulk‐like behavior. The application of this computational protocol to magnesium chloride (MgCl2) solutions at different concentrations (0.6–2.8 mol kg−1) gives h values in excellent agreement with experimental hydration numbers obtained using GHz‐to‐THz dielectric relaxation spectroscopy. This methodology is attractive because it is based on a well‐defined criterion for the definition of hydration number and provides a link with the molecular‐level processes responsible for affecting bulk solution behavior. Analysis of the ab initio molecular dynamics trajectories using radial distribution functions, hydrogen bonding statistics, vibrational density of states, water‐water hydrogen bonding lifetimes, and water dipole reorientation reveals that MgCl2 has a considerable influence on the hydrogen bond network compared with bulk water. These effects have been assigned to the specific strong Mg‐water interaction rather than the Cl‐water interaction.

Structural properties of the cation-water hydration shell TABLE S1. Structural properties of the cation-water radial distribution functions obtained from the ab initio MD simulations of the hydrated ions (isolated ion, no counterion). The positions, rmax, and amplitudes, g(rmax), of first peak and the average coordination number (CN) of the cation hydration shell are compared with other ab initio MD studies and available experimental data. Distances in Å.

Mg 2+
This study BOMD PBE-D3 GTH DZVP/PW (1000 Ry) 2.11 13.5 6.0 Di Tommaso [5] CPMD PBE USPP PW (30 Ry) 2.08 12.3 6.0 Callahan et al. [6] XRD 2.0-2.12 6 CaCl2 (6 m) Fulton et al. [8] 2.43 7.2 a) BOMD = Born-Oppenheimer molecular dynamics, CPMD = Car-Parrinello Molecular Dynamics. b) GTH = Goedecker-Teter-Hutter; USPP = Ultra-Soft Pseutopotential; PAW = projector augmented wave; NCPP = Norm-Conserving Pseudopotential. c) DZVP = double-zeta valence polarized; TZV2P = triple-zeta valence doubly polarized; PW = plane wave. Ion pairing of magnesium chloride  Figure S3.   The fraction (f) of water molecules with n number of hydrogen bonds per water molecule, the average number of hydrogen bonds (nHB) per water molecule in bulk water and in aqueous MgCl2 solutions at different concentrations (mol.kg -1 ), obtained from ab initio MD simulations (PBE-D3). Analysis was based on the following configurational criteria: two water molecules are hydrogen-bonded only if their inter-oxygen distance is less than 3.5 Å and, simultaneously, the hydrogen-oxygen distance is less than 2.45 Å and the oxygen-oxygenhydrogen angle is less than 30°. [9] System c (mol.kg -1 ) f0 f1 f2 f3 f4 f5 nHB Mg 2+ -Clion pairing  The fraction (f) of water molecules with n number of hydrogen bonds per water molecule, the average number of hydrogen bonds (nHB) per water molecule in bulk water and in aqueous MgCl2 solutions at different concentrations (mol.kg -1 ), obtained from classical MD simulations using the Lennard-Jones potentials developed by Aqvist [10] and Duboue-Dijon et al. [11] to describe the ion-ion and ion-water interactions together with the SPC/E water model. Analysis was based on the following configurational criteria: two water molecules are hydrogen-bonded only if their inter-oxygen distance is less than 3.5 Å and, simultaneously, the hydrogen-oxygen distance is less than 2.45 Å and the oxygen-oxygen-hydrogen angle is less than 30°. [9] Forcefield c (mol.  S3. Comparison of the average number of hydrogen bonds (HBs) in pure liquid water and aqueous MgCl2 solutions obtained from ab initio MD (PBE-D3) and classical MD (Aqvist [10] and Duboue-Dijon et al. [11] forcefields) simulations.

Distribution of hydrogen bonds in water and electrolyte solutions
In-shell hydrogen bond strength  For Mg 2+ , the first peak is shifted to lower distances and is more intense than Ca 2+ .

Protocol for the calculation of time correlation functions
The methodology below been adopted to compute the hydrogen bonding, , ( ), and reorientation time correlation functions (TCF), 1 ( ), defined as the first-order Legendre polynomials of water dipole, a unit bisector of the H-O-H angle ( ⃗).

Procedure for the categorization of water molecules
The hydrogen and oxygen atoms in each water molecule were labelled OaHbHc where a, b, c = 1, 2, and B in terms of the following criteria: 1 when the atom in the first coordination shell of the nearest ion; 2 when the atom is in the second coordination shell of the nearest ion; B when the above conditions are not met. Assignments were made by comparing the distance between Oa and the nearest magnesium ion with position of the first ( − 1 ) and second ( − 2 ) minima of the Mg-O RDF ( Figure S8), and the distance between Hb (or Hc) and the nearest chlorine ion with first ( − 1 ) and second ( − 2 ) minima of the Cl-H RDF ( Figure S9): Using these geometrical criteria at time t = 0, the water molecules in the MgCl2 solutions were classified in different water subpopulations labelled Wabc ( Figure S10). For example: W111 refers to a subpopulation of water molecules of type O1H1H1 where oxygen is in the first coordination shell of Mg 2+ and both hydrogen atoms are in the first coordination shell of Cl -; W112 refers to a subpopulation of water molecules of type O1H1H2 where oxygen is in the first coordination shell of Mg 2+ , a hydrogen is in the first coordination shell of Cland the other is the second coordination shell of Cl -. The continuous HB and dipole reorientational TCFs were then evaluated for all Wabc subpopulations.

Water exchange between different subpopulations
In this work, we simply tracked the water molecules belonging to a specific water subpopulation at the time origin (t = 0 ps) for calculating the dipole TCF of each water subpopulation: This approximation has been adopted by other groups to compute time correlation functions of dynamical properties of different water layers at solid/water and air/water interfaces. [12][13][14] However, water molecules that are in a particular subpopulation at t = 0 ps, e.g. W12B and W12B in Fig. S16, can exchange to "neighbouring" subpopulations. For example, a water molecule that is W12B initially goes to W1BB, a neighbour subpopulation, but molecules in both subpopulations contribute to the hydration water.
From the graphs in Fig. S16 we can distinguish two kinds of phenomena: short-time fluctuations and long-term exchange. Short term fluctuation (left figure) are fast exchanges of a water molecule between the initial subpopulation (t = 0 ps) and neighboring water subpopulations, which can occur because of rotational or vibrational motion of molecules that are at the margin of the distance criteria distinguishing between subpopulations. A long-term exchange (right figure) is seen when a water molecule exists for a considerable time in different water subpopulations and then returns to initial subpopulation or moves into another one; it could occur because of the translational motion of water between neighboring ions. To analyze the long term exchange in more detail, let us consider a simplified hydration model, Wab, that only considers the positions of the oxygen atom: a, b = 1, 2, or B, where 1 corresponds to O in the 1 st shell of Mg or Cl, 2 corresponds to O in the 2 nd shell of Mg or Cl, and B corresponds to O beyond the 2 nd shell of Mg or Cl. Using this simplified classification method, we can exclude the shortterm fluctuation and focus on the temporal behavior of the long term exchanges. Let us now define the following function: FIGURE S11. Exchange of water molecules between different watr subpopulations. Analysis conducted on the 0.6 mol.kg -1 MgCl2 solution.
∈ where (0) corresponds to the number of water molecules in the subpopulation Wab at time t = 0, and the quantity ( ) is computed as follows: The normalized, time-dependent function ( ) corresponds to the fraction of water molecules initially in Wab that have remained in the same category. Fig. S17 shows that W11 (water molecules in the first hydration shell of Mg and Cl) is quite "rigid" as more than 80% of water molecules remained in the original subpopulation. Also, about 80% of water molecule which is initially bulk water is still belongs to WBB under dynamic exchange equilibrium. In case of water molecules in the second shell either Mg 2+ and Clion, (WB2, W) between 50 and 60 % remains in the same category. In order to minimize the effect originating from the exchange of water molecules from one subpopulation to another, we used the first 8000 steps data (8ps) of the dipole correlation function to fit the biexponential model which give us the averaged relaxation time; this gives more than 70% of water molecules remaining in the initial category. Water molecules of each subpopulation are exchanged with each other, but the relaxation time of dipole correlation function we have obtained is fitted over time intervals that exchange effect is not significant to dipole correlation function because majority of water molecules remain same category.
FIGURE S12. Fraction of water molecules that have remained in the same category at time t. Analysis of 0.6 mol.kg -1 MgCl2 (aq).
Reorientation time correlation function of hydrated Mg 2+ and Cl -FIGURE S14. First-order Legendre reorientational TCFs, P1(t), of the water molecules in the first and second coordination shells of Mg 2+ and Cl -. Results compared to P1(t) of the water molecules in the bulk (beyond the second coordination shell of the ion).
Reorientation time correlation function analysis from classical MD FIGURE S15. Retardation factor for the reorientation relaxation time of the water subpopulations obtained from classical MD simulation of the 0.6 mol.kg -1 solution using the Duboue-Dijon forcefield. [15] Inset: number of water molecules per MgCl2 units in each subpopulation.

Connection between single-water molecule reorientational dynamics and dielectric relaxation spectroscopy measurements
The frequency-dependent dielectric constant (ε) is obtained from the Fourier-Laplace transform of the time-derivative of dielectric decay function (ϕ): where ϕ is given by the normalized autocorrelation function of the total dipole moment ⃗⃗ : [16] ( ) = < ⃗ ⃗ ( ) • ⃗ ⃗ (0) > < (0) 2 > (2) For an ensemble of (water) molecules, each characterised by an unit vector ⃗ defining the orientation of the molecular dipole moment, the total dipole moment is given by: Using this expression for total dipole moment in Eq. 2, the dielectric decay function can be written in terms of autocorrelation and cross-correlation of the molecular dipole: where ( ) denotes the autocorrelation function of the molecular dipole moments, and ( ) denotes the multi-molecular dipole cross-correlation function, Let us assume that the auto-correlation function of the molecular dipole ( ) and the autocorrelation of the total dipole ( ) are governed by one correlation time: where and are the macroscopic decay function and molecular dipole moment decay function, respectively. The expression connecting these two functions has been derived by Kivelson and Madden: [17] τ = [1 + < μ (0) ⃗⃗⃗⃗⃗⃗⃗⃗⃗⃗⃗ • ∑ 1 ( ) ⃗⃗⃗⃗⃗⃗⃗⃗⃗⃗⃗ ≠ > < (0) 2 > >] = (0) (9) Since the steady function (0) corresponds to the Kirkwood dipole orientation correlation factor, , [18] we can simply rewrite Eq. 8 as: = / (10) In our study, we have found that the relaxation time of the molecular dipole correlation function of certain water subpopulation are retarded (2 to 6 times slower) compared to bulk-like water molecules. On the other hand, is the same for all water subpopulation because the macroscopic relaxation time is obtained from the total dipole autocorrelation. it follows from Eq. 9 that water subpopulations with a retarded reorientation dynamics have a lower Kirkwood factor compare to the bulk-like water. Based on the conclusion made by Rinne et al. on single water versus collective water effects on ionspecific reorientation water dynamics, [19] the presence of an ion perturbs the cooperative structure of water in the subpopulations with a retarded reorientation dynamics because ion because of lower Kirkwood dipole orientation correlation factor. The lower dipole correlation factor of the water subpopulation with retarded dipole correlation makes this water subpopulation less contribute to bulk relaxation dynamics. This interpretation is also corroborated by the reduction of multi-molecular dipole cross-correlation contribution in the dielectric loss spectrum calculated, as suggested by Rinne et al. [19] Also, in the dielectric relaxation spectroscopy, the hydration number is calculated from the dielectric loss spectrum of electrolyte solution compare to bulk water. This leads to a decrease in the dielectric loss spectrum (static depolarization) and from static depolarization we can get hydration number.
To summarize, the hydration number computed from the water reorientation dynamic analysis of ab initio and classical MD simulations connects to THz-DR experimental values as follow: the presence of ions in solution can cause the single-molecule reorientation dynamic of some water subpopulations to be slower than bulk water; water subpopulation with retarded reorientational dynamics have lower correlation with other water molecules compared to bulk-like water; in these water subpopulation the multimolecular dipolar correlation contribution to the dielectric loss spectrum is decreased; as a result, the dielectric loss spectrum, which is measured by DRS is lower than that of bulk water.