Molecular Dynamics Simulation of Na+–Cl– Ion-Pair in Water

Nov 19, 2015 - *Telephone: +91-22-2576-4199. Fax: +91-22-2576-7152 E-mail: ... Local structures around the ions are studied using the radial distribut...
1 downloads 7 Views 2MB Size
Article pubs.acs.org/JPCB

Molecular Dynamics Simulation of Na+−Cl− Ion-Pair in Water− Methanol Mixtures under Supercritical and Ambient Conditions Sonanki Keshri, Atanu Sarkar, and B. L. Tembe* Department of Chemistry, Indian Institute of Technology, Bombay, Powai, Mumbai, 400076, India S Supporting Information *

ABSTRACT: Constrained molecular dynamics simulations have been performed to investigate the structure and thermodynamics of Na+−Cl− ion-pair association in water− methanol mixtures under supercritical and ambient conditions in dilute solutions. From the computed potentials of mean force (PMFs) we find that contact ion pairs (CIPs) are more stable than all other associated states of the ion pairs in both ambient and supercritical conditions. Stabilities of CIPs increase with increase in the mole fraction of methanol. In supercritical conditions, major changes in PMFs occur as we go from xmethanol = 0.00 to xmethanol = 0.50. The stable solvent shared ion pair (SShIP) which occurs in xmethanol = 0.00 and 0.25, vanishes when xmethanol is 0.50 or greater. The stabilities of these ion pairs increase with increasing temperature. Local structures around the ions are studied using the radial distribution functions, density profiles, angular distribution functions, running coordination numbers and excess coordination numbers. Preferential solvation analysis shows that both Na+ and Cl− ions are preferentially solvated by water. From the calculation of enthalpies and entropies, we find that Na+−Cl− ion-pair association in water−methanol binary mixtures is endothermic and driven by entropy both in ambient as well as under supercritical conditions.

I. INTRODUCTION At temperatures and pressures above its critical point, a substance exists in a state that exhibits the properties of both liquid and gas.1−3 The properties of supercritical fluids (SCFs) such as compressibility, diffusivity, density, and solvation of different solutes can be modulated through minor variations of temperature and pressure when the system is very close to the its critical conditions.4−7 These fascinating features of SCFs lead to their important potential as environmentally benign green alternatives to toxic organic solvents and there have been many industrial and technological applications using SCFs. Industrial separations using SCFs are well-known with the most famous example being decaffeination.8 In conditions of the high dilution that prevail in different processes, solute’s thermodynamic properties scale as the solvent isothermal compressibility. Ion solvation under different temperatures and pressures govern many important biological, geological and chemical processes. In particular, the thermodynamic properties of ions in supercritical solvents have crucial roles in many geological processes.9 Therefore, it is important to simulate the solvation structure and dynamics of ions under both ambient and supercritical conditions at the microscopic level for understanding these physical and chemical processes in such hydrothermal systems. Alcohol−water mixtures are interesting and important for many fields ranging from fundamental molecular science to extensive industrial applications. Mixtures of water and alcohols are also interesting because of their complex dynamics, which results from the presence of hydrophobic groups and hydrogen© XXXX American Chemical Society

bonding. These mixtures are also important because of their extensive use as solvents and are thus in the center of longstanding experimental and theoretical investigations.10−17 Simple alcohols such as methanol and water mix completely from a macroscopic viewpoint, but their solutions are known to be microscopically inhomogeneous and show nonideal properties at ambient conditions and also at higher temperatures and pressures. This heterogeneity is induced by hydrogen bond interactions existing in solution.18 Methanol as well as water are highly associated liquids, but there is a considerable difference in the nature of their hydrogen bond network. Water molecules form a three-dimensional, tetrahedrally coordinated structure, via hydrogen bonded network.19,20 Methanol, on the other hand, forms a zigzag polymer chain through hydrogen bonds.21−24 In the case of water−methanol mixtures, both the solvents are connected by a common hydrogen-bonded network.25 The mixture is still heterogeneous at a molecular level, as the bulkier methyl group of methanol remains unable to replace the hydrogen from the tetrahedrally coordinated structure of liquid water. Recent studies have also shown that the structure and dynamics of hydrogen-bonded networks of water−methanol mixtures depend strongly on their composition.26 In recent years, methanol and water−methanol mixtures under supercritical conditions have attracted the attention of different research groups due to their industrial and Received: June 6, 2015 Revised: November 17, 2015

A

DOI: 10.1021/acs.jpcb.5b05401 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B

simulations was studied by Bujnicka et al.47 They have found that the chloride ion is preferentially solvated by water and calcium ion is preferentially solvated by methanol molecules. Rybicki et al. studied the solvation of Mg2+ in water−methanol binary mixtures.48 Their study indicated preferential solvation of Mg2+ ion by water molecules. Experimental solvation studies in binary mixtures deficient in water, probed using NMR relaxation techniques, revealed that sodium ion was preferentially hydrated whereas chloride ion did not show specific solvation by either components.49−51 In methanol rich solvents, self-diffusion measurements52 supported preferential hydration of alkali ions, but for halide ions, the preferential solvation was by methanol. Solvation structure of Na+−Cl− ion pair in supercritical water has been extensively studied by Cummings et al.53−61 wherein they have studied the thermodynamics and kinetics of ion pair association in supercritical water. They have also studied the hydrogen bonding in supercritical water.61 In our earlier study, we have performed MD simulations of Na+−Cl− ion pair in supercritical methanol as a function of temperature.62 Thus, extensive computer simulation studies of the ion pairs in the water−methanol mixtures can be devised to give reasonable structural predictions for the mixtures under study and can be extended to other hydrogen bonded liquids. Although water−methanol binary mixtures have been extensively investigated at room temperature both theoretically and experimentally, very few theoretical and experimental studies of this system have been made at supercritical conditions. To the best of our knowledge, studies of the solvation structure of the ion-pairs in binary mixtures of water−methanol under supercritical conditions have not yet been reported. The aim of this work is to determine the solvation structure around the Na+− Cl− ion pair in water−methanol mixtures under supercritical conditions with the goal of developing a molecular understanding of the association structure of ion pair under supercritical conditions. Since experimental density data for the mixtures are available at 623 K and 350 bar, all our simulations correspond to these conditions. At 623 K, pure water is below the critical temperature of water which is 647 K but we have performed simulations at 623 K to compare these PMFs with the PMFs of other compositions of water− methanol mixtures at 623 K. Our main purpose here is to compute the PMFs at low concentrations (the molarity of Na+−Cl− is ∼0.025 M in all the mixtures studied), i.e., for dilute solutions, wherein the ion−ion effects are small and the structural properties may be different in more concentrated solutions. Since corresponding PMFs at the same molarity were not available in ambient conditions, we have also computed them for comparison. The methodology and computational details are described in section II. Results and discussion are given in section III and our conclusions are summarized in section IV.

technological applications. From a technological point of view, such mixtures have served as useful industrial solvents for a variety of separation processes. Recently they are also used in solar thermal systems.27 Physical chemists have been attracted by their eccentric, unusual nonideal behavior, especially in the low concentration range.28 Water and methanol have been proposed as green solvents for conversion of biomass,29 biodiesel production,30 and chemical recycle of monomers via depolymerization.31 Some of these processes use supercritical water and others use supercritical methanol. On the basis of this, it is likely that water−methanol mixtures can be advantageously used as solvent medium since solvent properties can be varied around those of methanol and water continuously at appropriate temperature and pressure conditions. Dixit et al.32 observed a microscopically inhomogeneous evidence of water−methanol mixtures at room temperature via molecular dynamics simulations and attributed this behavior to small but non-negligible hydrophobic interactions in aqueous solutions. Thus, solvent properties at higher temperature will be greatly affected by hydrophobic groups because at higher temperature, the hydrogen bond is considerably weakened as shown by Raman measurements for water by Ikushima et al.33 For supercritical water−methanol mixtures, hydrogen bond interactions will be completely different from those at ambient conditions. Thus, a deeper understanding of the microscopic structural properties of supercritical aqueous alcoholic mixtures will lead to marked improvements in practical applications for environmental, mechanical, chemical, biological, and geothermal industries. Variation of temperature and pressure as well as addition of electrolytes affects the solvation structure.34 For aqueous methanolic solutions, many experimental and theoretical studies have been carried out in ambient conditions which provide an insight into the structural, thermodynamic and hydrogen bonding properties. Hawlicka and Swiatla-Wojcik35−40 employed the flexible Bopp−Jancso−Heinzinger (BJH)41,42 and Pálinkás−Hawlicka−Heinzinger (PHH)43 models of water and methanol respectively, in their MD simulations of water−methanol mixtures containing an Na+− Cl− ion-pair. They observed that in methanol deficit mixtures, cations and anions are preferentially solvated by methanol molecules.35 In methanol rich systems, no preferential solvation of cations has been found; but the anions are solvated only by methanol molecules.35 In their MD simulations of ions in binary mixtures, Hawlicka et al.37 observed that in water rich mixtures and in equimolar mixtures, the cations are preferentially solvated by methanol. They have found a weak preferential hydration of Na+ in methanol-rich solvents. In their more recent study40,44 they have observed that preferential solvation of anions by methanol molecules becomes less pronounced with decreasing charge density of a solute and it vanishes for the discharged chloride ion. In contrast to preferential hydration of Na+ in water deficit solvents, its uncharged counterpart has been found to be preferentially solvated by methanol molecules over the whole range of solvent composition. However, Day et al.45 have found that both Na+ and Cl− ions are preferentially solvated by water molecules in water−methanol binary mixtures. Feakins and Watson have studied the preferential solvation of ions in nonaqueous solvents and their aqueous mixtures.46 Owczarek et al. studied the influence of CaCl2 on micro heterogeneity of water−methanol mixtures by molecular dynamics simulations.34 Solvation of Ca2+ in aqueous methanol via MD

II. METHODOLOGY AND COMPUTATIONAL DETAILS We have performed MD simulations using GROMACS package (version 4.5.4).63 In the present study, the TIP3P model of Jorgensen et al. has been used for water.64 The united atom optimized potentials for liquid simulation (OPLS)65 force fields have been used for methanol. For Na+−Cl− ion pair, the potential model of Smith and Dang has been used. The initial configurations of the system were generated using Packmol.66 The force fields and geometrical parameters for water, methanol and the Na+−Cl− ion pair are given in Table S1 B

DOI: 10.1021/acs.jpcb.5b05401 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B Table 1. Details of the Simulation Cellsa xMeOH

NMeOH

NWater

Ntotal

NNa+

NCl‑

ρ (kg/m3) (exptl.)

0.00 0.25 0.50 0.75 1.00

− 234 341 413 462

1411 704 341 138 −

1411 938 682 551 462

1 1 1 1 1

1 1 1 1 1

659.41 524.02 443.22 408.18 384.52

ρ (kg/m3) (simulations) 653.62 521.12 440.21 398.05 381.49

± ± ± ± ±

3.23 3.39 3.15 3.50 3.10

Key: xMeOH = mole fraction of methanol and xWater = (1 − xMeOH), NMeOH and NWater = number of molecules of methanol and water, NNa+ and NCl‑ = number of Na+ and Cl− ions respectively, Ntotal = total number of solvent molecules, L = simulation box length, and ρ = density at 623 K and 350 bar pressure. a

Table 2. Details of the Simulation Cellsa xMeOH

NMeOH

NWater

NTotal

NNa+

NCl‑

ρ (kg/m3) (exptl.)

0.00 0.25 0.50 0.75 1.00

− 375 681 954 952

2137 1127 681 318 −

2137 1502 1362 1272 952

1 1 1 1 1

1 1 1 1 1

998.99 838.13 884.50 842.24 791.02

ρ (kg/m3) (simulations) 997.59 836.42 880.46 841.02 790.21

± ± ± ± ±

3.46 5.97 3.22 3.30 3.13

Key: xMeOH = mole fraction of methanol and xWater = (1 − xMeOH), NMeOH and NWater = number of molecules of methanol and water, NNa+ and NCl‑ = number of Na+ and Cl− ions respectively, Ntotal = total number of solvent molecules, L = simulation box length, and ρ = density at 298 K and 1 bar pressure. a

Half of the box length is used as a cutoff for Lennard-Jones forces. The electrostatic interactions are treated by PME method,70 with a Coulomb cutoff of 1.5 nm and an interpolation order of 4. For nonbonded van der Waals interactions, a 1.5 nm cutoff is used. The potential energy of our system is minimized using the steepest decent algorithm with a tolerance for force of 1000 kJ mol−1 nm−1 and convergence is obtained in all the cases. Subsequent to energy minimization, we have performed equilibration for 2 ns in the NPT ensemble. Following this equilibration procedure, MD simulations are initiated. Periodic boundary conditions (PBCs) are used along with the minimum image criterion.71 The neighbor list is updated every 10 steps. SHAKE algorithm72 is used to maintain the constant bond lengths and bond angles of the solvent molecules during simulations. The temperature of the system is maintained at 623 K using a velocity rescaling thermostat73 with a relaxation time of 0.1 ps. Pressure of the system is fixed at 350 bar using the Berendsen barostat74 with a relaxation time of 0.5 ps. The equations of motion are integrated using the leapfrog algorithm75 with a time step of 2 fs. At the beginning of the simulation, the velocities of the atoms were assigned from a Maxwell distribution at a desired temperature. We have performed 81 simulations ranging from 0.2 to 1.0 nm with a constant distance interval of 0.01 nm. After 2 ns equilibration, we have performed 10 ns simulations for the calculation of the potentials of mean force (PMF) between Na+ and Cl− ion using the Parrinello−Rahman barostat.76 The PMFs between the Na+ and Cl− ion pair are calculated by integrating the mean forces acting on them. The mean force at a separation r is the sum of the solute−solute direct force and the ensemble average of the solute−solvent forces. That is,

and Table S2 of the Supporting Information. The intermolecular interactions are taken to be pairwise additive and composed of the Lennard-Jones and Columbic terms. Uij(r ) =

Aij r12



Bij r6

+

qiqj r

(1)

Here, i and j denote a pair of interaction sites on different molecules, qi = charge located at site i and qj = charge located at site j, r = site−site separation. The terms Aij and Bij are determined from, Aij = 4 × (εij) × (σij)12

(2)

Bij = 4 × (εij) × (σij)6

(3)

whereas, εij and σij are calculated using the Lorentz− Berthelot mixing rules.67 εij = (εi × εj)1/2

(4)

⎛ σi + σj ⎞ σij = ⎜ ⎟ ⎝ 2 ⎠

(5)

The experimental critical temperatures of water−methanol mixtures as a function of methanol mole fraction have been reported by Marshall et al. but they have not reported the critical pressure data.68 Bulemela et al.69 have studied the volumetric behavior of water−methanol mixtures in the vicinity of the critical region and they have given an empirical formula for calculating critical pressure. In our simulations, we have used the experimental critical temperatures and densities and the critical pressures estimated by the method of Bulemela et al.69 These parameters are given in Table S3 in the Supporting Information. Our computed densities agree with experimental critical densities very well and thus we conclude that potential models reproduce critical conditions quite well. The densities obtained in our MD simulations and those reported experimentally are given in Table S4 in the Supporting Information.

F(r ) = Fd(r ) + ΔF(r )

(6)

where, ΔF(r) = ⟨F(r,t)⟩. The angular brackets denote the ensemble average. Since we are fixing the ion−ion distance, it reduces their entropy which in turn is equivalent to an entropic force of magnitude −2kBT/r. This is taken into account as77−79 C

DOI: 10.1021/acs.jpcb.5b05401 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B 2k T − d W (r ) = F (r ) − B dr r

pair (SShIP) and the third minimum (which is shallow) represents the solvent separated ion pair (SSIP).83 In each mixture, the position of the CIP is around 0.26 nm (Figure 1, parts a and b). Figure 1 also indicates that SShIP is not a stable structure at supercritical conditions. It is seen that with increase in the mole fraction of methanol, the contact ion-pair becomes more stable. This is because of the lower values of dielectric constants as we go on increasing the methanol mole fraction. The addition of alkali and alkaline earth metal halides in water and methanol also results in a decrease in the dielectric constant of the solvent.84 Low dielectric constants of solvents cannot screen the ion-pair effectively and this helps in the association of the ion-pair. As a result, CIPs become more stable. The CIPs in Figure 1a are more stable than the CIPs in Figure 1b. In case of supercritical conditions, the depths of the CIPs range from −151.7 kJ/mol (xMeOH = 1.00) to −54.3 kJ/ mol (xMeOH = 0.00) (Figure 1a). For the ambient case, the depths of the CIPs range from −31.2 kJ/mol (xMeOH = 1.00) to −11.0 kJ/mol (xMeOH = 0.00). The CIPs are about 5 times more stable under supercritical conditions than that under ambient conditions. This is mainly because of low dielectric screening by the solvent at higher temperatures which helps in ion-pair association. The position of CIP and SShIP in supercritical water agrees well with previous study of Yui et al.85 The depths of CIP as well as SShIP for xMeOH = 0.00 obtained in our case are higher than those obtained by Yui et al.85 The difference in depth may be a result of using different potential model for solvent molecules. The position and depth of CIP for xMeOH = 1.00 under supercritical conditions are in agreement with our previous study.62 Position of CIP and SShIP in pure methanol in ambient conditions is in good agreement with the results obtained by Dixit et al.86 The energy differences between CIPs and transition states (TSs) lie in the range of 25−40 kJ/mol (for supercritical conditions). The mixture compositions with xMeOH = 0.00 and 0.25 exhibit weak second minima at 0.51 nm which correspond to solvent shared ion-pairs (Figure 1a). The depths of SShIPs under supercritical condition are −26.6 kJ/mol (xMeOH = 0.00) and −47.6 kJ/mol (xMeOH = 0.25). There is a shoulder instead of an SShIP minimum when xMeOH = 0.50. Also, when xMeOH = 0.75 and 1.00, the PMFs indicate the presence of only the CIP species (Figure 1a). The minimum corresponding to SShIP, as observed in the case of solvent mixtures with lower values of xMeOH does not appear. From Figure 1b, it is observed that there is a presence of stable SShIPs in all compositions. The position of SShIP lies between 0.46 and 0.51 nm. The distances at which we observe the minima for Na+−Cl− ion pair in all compositions, and the corresponding values of PMFs are summarized in Table S7 (supercritical conditions) and Table S8 (ambient conditions) of the Supporting Information. The depths of SShIPs in ambient conditions range from −19.9 kJ/mol (xMeOH = 0.00) to −4.4 kJ/mol (xMeOH = 1.00). In each composition of water− methanol mixtures, the depth of the CIP minimum is more than the depth of the SShIP. This shows that the contact ion pair is more stable than the solvent shared ion pair in each composition of water−methanol mixtures in both ambient and supercritical conditions. In ambient conditions as seen from Figure 1b, we also observe a third minimum corresponding to solvent separated ion pair at around 0.72 nm. The depth of SSIPs range from −11.1 kJ/mol (xMeOH = 0.00) to −2.9 kJ/mol (xMeOH = 1.00). Depth of SSIPs also increase with increase in methanol mole fraction.

(7)

Integration of eq 7 yields the potential of mean force, W(r), W (r ) = W (r0) −

∫r

0

r

⎛r⎞ F(r ) dr + 2kBT ln⎜ ⎟ ⎝ r0 ⎠

(8)

The choice of W (r0) is required to be made in such a way that the calculated mean force potential matches the macroscopic Coulombic potential at long distances. qi × qj W (r0) = εr r0 (9) Here, εr is the dielectric constant of the solvent. It has been found that the ion - ion PMFs studied here are not sensitive to the choice of r0 and so r0 is chosen to be 1.0 nm.80 The error bars in PMFs for all cases are found to be around 0.2−0.5 kJ/ mol. We have used experimental densities of water−methanol binary mixtures.81,82 The details of the chosen solvent mixtures are listed in Table 1 and Table 2. The static dielectric permittivity (εr) was computed in the NPT ensemble. For the calculation of dielectric constant, we have performed MD simulations for 100 ns. Calculated dielectric constants of water−methanol binary mixtures are given in Table S5 (supercritical conditions) and Table S6 (ambient conditions) of Supporting Information. In the calculations of the ion pair distance residence times, we have performed 500 ps NPT trajectories for each of the 81configurations in the 0.20−1.0 nm region after releasing the constraint between the ions.

III. RESULTS AND DISCUSSION A. Potentials of Mean Force. The ion−ion PMFs for the five solvent compositions have been obtained by a direct integration of the total force according to eq 7. In parts a and b of Figure 1, we have shown the PMFs for the five solvent compositions.

Figure 1. Potentials of mean force of the Na+−Cl− ion-pair in water− methanol mixtures of different compositions at/near supercritical conditions [T = 623 K, P = 350 bar] (a) and ambient conditions [T = 298 K, P = 1 bar] (b).

Parts a and b of Figure 1 correspond to PMFs of Na+−Cl− ion-pair under supercritical conditions and ambient conditions, respectively. The PMFs for Na+−Cl− ion-pair under supercritical conditions are characterized by one deep minimum followed by a shallow minimum at low mole fractions of methanol (Figure 1a), whereas the PMFs for Na+−Cl− ion-pair in ambient conditions are characterized by three minima (Figure 1b). The first minimum represents the contact ion pair (CIP), the second minimum represents the solvent shared ion D

DOI: 10.1021/acs.jpcb.5b05401 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B A remarkable feature of the PMF plots of Na+−Cl− ion-pair in water−methanol mixtures under supercritical conditions is that, with increase in the mole fraction of methanol beyond 0.25, SShIPs vanish. SSIPs are not observed under supercritical conditions. Strong interaction between Na+ and C1− ion at low fluid density would prevent the ions from dissociating into the bulk, thus reducing the possibility of SShIP formation. A rationale for the absence of SShIP with increase in mole fraction of methanol molecules will be given after analyzing the hydrogen bonded networks that develop around the ion-pair in these solvent mixtures. Although the molarity of Na+−Cl− is ∼0.025 M in all the mixtures studied here, the molalities do differ. For the system containing one ion pair in 1411 water molecules (pure water, supercritical conditions), the molality is 0.04 m whereas for the system containing one ion pair in 462 methanol molecules (pure methanol supercritical conditions), the molality is 0.08 m. To check the effect of molality on the structural properties of the ion pair, we have performed simulations of 0.04 m methanolic solution (the system contains one ion pair in 800 methanol molecules, supercritical conditions). It is observed that there is marginal difference (less than about 1%) in the depth of CIP between the two molalities. We have also checked the effect of molality for xMeOH = 0.25, 0.50, and 0.75 and the difference observed is around 1 kJ/mol. The PMF curves as a function of molality for xMeOH = 0.25, 0.50, 0.75, and 1.00 are given in the Supporting Information in Figure S1a−d. The PMFs are confirmed by calculating ion pair distance residence times (IPDRTs) which is an approximate measure of the length of time the ion pair resides at a particular distance. The IPDRTs for ambient and supercritical conditions are given in parts a and b of Figure S2 of the Supporting Information. The formula used to calculate IPDRTs is described in detail elsewhere.62 Sharp peaks at 0.26 nm and small broad peaks at 0.51 nm of IPDRTs correspond to the CIPs and SShIPs at xMeOH = 0.00 and 0.25. In ambient conditions, we observe an intense peak at around 0.27 nm which is followed by two small broad peaks at 0.51 and 0.74 nm, respectively. These peaks confirms the presence of CIPs, SShIPs, and SSIPs, respectively. Peak heights of IPDRTs are also consistent with the stabilities of CIPs at both ambient as well as supercritical conditions. It is also noticed that IPDRTs have the highest value for xMeOH = 1.00. To understand the extent of dissociation of Na+−Cl− ion pair in ambient and supercritical conditions in water−methanol mixtures, we have calculated the ion-pair dissociation constants.62 The dissociation constants of Na+−Cl− ion pair in ambient and supercritical conditions are given in Tables S9 and S10 in the Supporting Information. According to TST, rate of a reaction slows down if the activation energy barrier is high. From the PMFs of Na+−Cl− ion pair in water−methanol mixtures, it is seen that the CIP is more stable in pure methanol. The rate of ion pair dissociation reaction for Na+− Cl− ion pair in water−methanol binary mixture also follows the same trend as seen from the Tables S9 and S10. The CIP is more stable in pure methanol, and the energy barrier between CIP and T.S. is very high which results in a lower rate of ionpair dissociation. B. Ion Solvent Radial Distribution Functions under Supercritical Conditions. In order to analyze the local solvation structures around Na+ and Cl− ions in five mixed solvents, we have computed the ion−solvent radial distribution functions, g(r)ion−solvent (RDFs) under supercritical conditions.

The calculations in ambient conditions have already been reported by Hawlicka et al.29 The ion−solvent radial distribution functions have been calculated for all the solvent mixtures. Parts a and b of Figure 2 represent the ion−solvent

Figure 2. Ion−solvent radial distribution functions for the Na+ ion in water−methanol mixtures under supercritical conditions at 623 K and 350 bar, (a) Na+−O (H2O) and (b) Na+−O (MeOH).

Figure 3. Ion−solvent radial distribution functions for the Cl− ion in water−methanol mixtures under supercritical conditions at 623 K and 350 bar, (a) Cl−−O (H2O) and (b) Cl−−O (MeOH).

RDFs for the Na+ ion and parts a and b of Figure 3 represent the ion−solvent RDFs for the Cl− ion. From parts a and b of Figure 2, we see that Na+−O (MeOH) and Na+−O (H2O) do not exhibit any change in the positions of peak maxima, but the intensities of the peaks change with change in solvent compositions. The peaks for Na+−O (H2O) and Na+−O (MeOH) are centered at around 0.24 nm. Position of the peak maxima are in good agreement with those reported by Hawlicka et al. in their study in ambient conditions.40,44 Thus, solvent composition does not alter peak position. From Figure 2a, we observe that when the mole fraction of water is the lowest, the peak heights of Na+−O (H2O) are the largest. The heights of the RDFs for these pairs decrease as the mole fraction of water increases. With increase in the water mole fraction, the bulk water density increases, but the local water density around the Na+ ion does not increase proportionately. This makes the peak height decrease with increasing water content in the mixtures. This will be rationalized further after noting the values of running coordination numbers. To study the effect of composition on the solvation structure around the ion pair, we have calculated the running coordination numbers around the ions. The running coordination number is defined as nαβ = 4πρβ E

∫0

R min

r 2gαβ (r ) dr

(10) DOI: 10.1021/acs.jpcb.5b05401 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B Table 3. Running Coordination Numbers (RCNs) of Na+ and Cl− Ions in Supercritical Water−Methanol Mixturesa Cl−

Na+

compositions xH2O

xMeOH

Rmin (MeOH)/nm

nMeOH

Rmin (H2O)/nm

nwater

Rmin (MeOH)/nm

nMeOH

Rmin (H2O)/nm

nwater

1.00 0.75 0.50 0.25 0.00

0.00 0.25 0.50 0.75 1.00

0.33 0.34 0.34 0.35 0.34

0.00 0.45 1.15 2.11 3.48

0.33 0.34 0.33 0.33 0.34

3.82 3.08 2.42 1.58 0.0

0.53 0.54 0.52 0.53 0.52

0.00 1.41 3.41 5.96 8.32

0.53 0.54 0.55 0.56 0.54

10.10 9.21 6.46 3.86 0.00

Key: Rmin (MeOH) = first minima in g(r) of Na+−O (MeOH) and Cl−−O(MeOH), Rmin (H2O) = first minima in g(r) of Na+−O (H2O) and Cl−− O(H2O), nMeOH = number of MeOH molecules, and nwater = number of H2O molecules. a

where nαβ represents the number of atoms of type β surrounding species α in a shell extending from zero to Rmin (the first minimum in the radial distribution function). The running co-ordination numbers of water and methanol molecules around Na+ and Cl− ion at CIPs are presented in Table 3. It is observed from Table 3 that when the mole fraction of methanol increases from 0.0 to 1.0, the number of methanol molecules increases around Na+ and Cl− ion in the CIP. Also, the number of water molecules decreases with increase in the mole fraction of methanol. From Table 3, we can also see that both Na+ and Cl− ion is preferentially solvated by water molecules for xMeOH = 0.25, 0.50, 0.75. The selective solvation of ions in binary mixtures is determined by an aggregate effect of two factors with opposite nature, namely, (i) attractive ion− dipole interactions between the ion and solvent molecules and (ii) the steric factor due to the size of solvent molecule. Methanol molecule being bulkier than water, the solvation shell is dominated by water molecules. A similar trend is observed for the cases of TS and SShIP configurations. The trends in the peak heights of Na+−O (H2O) RDFs can now be understood in terms of RCNs. For xH2O = 0.25, the ratio of Na/Ni at Rmin = 0.35 nm is 1.58/0.995 is nearly equal to 1.58. In the previous ratio, Na is the actual coordination number at the mole fraction under consideration and Ni is the ideal mole fraction that would result if the local mole fraction is the same as the bulk mole fraction. Since the RCN of pure H2O is 3.82, the value Ni is 3.82 × 0.25 = 0.995 for xH2O = 0.25. When xH2O increases to 0.50 and 0.75, the corresponding ratios are 2.42/1.91 (=1.26) and 3.08/2.86 (=1.07). As these ratios are less than 1.58 and are analogous to the ratios of local densities to bulk densities, the peak heights of Na+−O (H2O) and Na+− H (H2O) RDFs for xH2O = 0.50 and xH2O = 0.75 are smaller than the RDF peak height for xH2O = 0.25. Figure 2b presents the Na+−O (MeOH) RDFs. In the case of MeOH, the peak heights of Na+−O (MeOH) RDFs increase with an increase in the mole fractions of methanol as local densities increase with increasing bulk densities. The peak heights of Na+−O (H2O) RDFs are larger than the peak heights of Na+−O (MeOH) RDFs. A comparison of radial distribution functions for Na+ in aqueous and methanolic solutions shows that the average distance of cation to methanol’s oxygen is almost the same as the average distance of cation to water’s oxygen. This feature is in good agreement with the experimental results.87,88 The RDFs of Na+−H (H2O) and Na+−H (MeOH) are given in Supporting Information in Figure S3, parts a and b. The first peak of Na+−H (H2O) and Na+−H (MeOH) are broader and lower than that of Na+−O (H2O) and Na+−O (MeOH). The

peaks are centered at around 0.29 nm i.e. they are shifted to a longer distance of 0.05 nm. The trends in peak heights are similar to the case of Na+−O (H2O) and Na+−O (MeOH). In aqueous and methanolic solutions of Na+−Cl−, the first peak of Cl−−O (H2O) and Cl−−O (MeOH) is centered at 0.32 nm. This agrees well with the average distance deduced from Xray diffraction from the chloride ion to oxygen in aqueous89 and methanolic87,88 solutions. The peak position of Cl−−O (H2O) and Cl−−O (MeOH) is in excellent accord with the distance of 0.35 nm reported by Hawlicka et al.40,44 From Figure 3a, we observe that when the mole fraction of water is the lowest, the peak heights of Cl−−O (H2O) are the largest. The heights of the RDFs for these pairs decrease as the mole fraction of water increases. With increase in the water mole fraction, the bulk water density increases, but the local water density around the Cl− ion does not increase proportionately. This makes the peak height decrease with increasing water content in the mixtures. The trends in the peak heights of Cl−−O (H2O) RDFs can now be understood in terms of RCNs. For xH2O = 0.25, the ratio of Na/Ni at Rmin = 0.55 nm is 3.86/2.525 is nearly equal to 1.52. Na and Ni have been defined above. Since RCN of pure H2O is 10.10, the value Ni is 10.10 × 0.25 = 2.525 for xH2O = 0.25. When xH2O increases to 0.50 and 0.75, the corresponding ratios are 6.46/5.05 (= 1.27) and 9.21/7.575 (= 1.21). As these ratios are less than 1.52 and are analogous to the ratios of local densities to bulk densities, the peak heights of Cl−−O (H2O) and Cl−−H (H2O) RDFs for xH2O = 0.50 and xH2O = 0.75 are smaller than the RDF peak height for xH2O = 0.25. Figure 3(b represents the Cl−−O (MeOH) RDFs. The composition of the mixed solvent does not affect the peak positions, but it influences remarkably the peak heights. The trend observed is similar to that of Na+−O (MeOH) RDFs. The RDFs of Cl−−H (H2O) and Cl−−H (MeOH) are given in Supporting Information in Figure S4, parts a and b. There are two very close peaks in the RDFs of Cl−−H (H2O) and Cl−−H (MeOH). Such close peaks in the RDFs indicate the hydrogen bonding formed between Cl− ion and the hydrogen atom of water. The first peak here is centered at around 0.23 nm, which is about 0.09 nm shorter compared to the peaks in the RDFs of Cl−−O (H2O) and Cl−−O (MeOH). The peak positions agree well with those reported by Hawlicka et al. under ambient conditions.40,45 C. Preferential Solvation of the Na+−Cl− Ion-Pair under Supercritical Conditions. In order to examine the solvent preferences of the Na+−Cl− ion-pair in water− methanol binary mixtures, we have calculated the local mole fractions of water and methanol in the ion solvation shells as a function of the bulk solvent composition. The local mole F

DOI: 10.1021/acs.jpcb.5b05401 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B fraction of solvent molecules β in the first solvation shell around ions (α) is given by xαL(β ,

R min) =

ρβ ∫

R min

0

c

∑ j = 1 ρj ∫

0

fractions (in the bulk) are shown in Figure 4. The straight line on which xL is equal to xMeOH represents ideal behavior. We note that the departures from ideality are significant especially in the range of xMeOH from 0.25 to 0.75. The deviations from ideality are greater in ambient conditions in the above composition range (by about 10 to 20%). From Figure 4, it is observed that both the cation and the anion are preferentially solvated by water molecules for xMeOH = 0.25, 0.50 and 0.75. D. Excess Coordination Numbers (ECNs) under Supercritical Conditions. To investigate the preferential solvation of Na+ and Cl− ion pair in water−methanol mixtures, we calculate the excess coordination numbers (ECNs) of solvent around the ion pair. The excess coordination number is defined as90

gαβ (r )4πr 2 dr

R min

gαβ (r )4πr 2 dr

(11)

xLα(β,Rmin)

The plot of is shown in Figure 4. Here, c is the number of solvent components. From the figure, we find that

∫0

Nαβ(r ) = 4πρβ

R

r 2[gαβ (r ) − 1] dr

(12)

Here, 4πρβr gαβ(r) is the average number of β molecules around an α molecule in a spherical shell of width dr at radius r and 4πρβr2dr is the average number of β molecules in that spherical shell that are uninfluenced by α (as in a uniformly random distribution). Therefore, Nαβ(r) represents the excess number of α molecules around a β molecule measured up to a certain distance r = R. The ECNs of water and methanol around the ion pair is given in Figure S5a−f [for Na+] and Figure S6a−f [for Cl−] (Supporting Information). The Kirkwood−Buff theory gives relations between integrals of radial distribution functions and properties of solutions. The Kirkwood−Buff integral between species u and v is defined as90 2

Figure 4. Local mole fractions in water−methanol binary mixtures of different partial mole fractions under supercritical and ambient conditions.

both Na+ and Cl− ions are preferentially solvated by water molecules in the mixed solvents. When xH2O decreases from 1.00 to 0.75, the local mole fraction of water around Na+, xLwater (=(nwater/(nwater + nmethanol))) still remains high, around 0.87 (=(3.08/(3.08 + 0.45)), Table 3, row 2). As xH2O is further decreased to 0.50 and 0.25, the values of xLwater are 0.67 and 0.42 respectively. When xMeOH is increased from 0.00 to 0.25, the L local mole fraction of methanol around Na + , x methanol (=(nmethanol/(nmethanol + nwater))) decreases to a much lower value of 0.13 (=(0.45/(0.45 + 3.08)), Table 3, row −2). This is due to insertion of water molecules in the first solvation shell. For xMeOH = 0.50 and 0.75, the values of xLmethanol around Na+ are 0.32 and 0.56 respectively (Table 3, rows −3 and 4). Therefore, the above values indicate greater affinity of Na+ toward water molecules under supercritical conditions. When the values of xH2O are 0.75, 0.50 and 0.25, the corresponding xLwater values around Cl− ion are 0.88, 0.66, and 0.41 respectively (Table 3, rows -2, 3 and 4) respectively. On the other hand when xMeOH decreases from 1.00 to 0.75, the value of xLmethanol decreases to 0.63 (Table 3, row -4). With a further decrease in the values of xMeOH to 0.50 and 0.25, xLmethanol values drop to 0.34 and 0.13 (Table 3, row −2 and row −3) respectively. Thus, it is clear that Cl− ion has a greater affinity toward water molecules. Hawlicka et al. in their study have found selective solvation of anions (Cl−, I−) by methanol molecules and preferential hydration of Na+ was found in water deficit solutions.40,44 In their MD simulation study of the NaCl ion pair in water− methanol mixtures, they have observed that in water rich and in equimolar mixtures, Na+ ion is preferentially solvated by methanol and only in methanol-rich solvent is a weak preferential hydration of Na+ found.37 Day et al.45 have found preference for water for both positively (Na+) and negatively (Cl−) solutes in water−methanol mixtures. Thus, our results of preferential solvation of Na+ and Cl− ion are in good agreement with previous studies. The local mole fractions (xL) of the solvent molecules in the first solvation shell around the ions as a function of the mole

Guv(r ) = 4π

∫0

r

[guv(r′) − 1]r′2 dr′

(13)

where u and v represent solute and solvent molecules, guv(r) is the radial distribution function between species u and v. Chialvo and Cummings58 suggested that a more appropriate measure of the aggregation around a solute molecule might be to consider the difference between guv(r) and gvv(r), the solvent−solvent correlation function. Accordingly, they defined the excess number of solvent molecules around a solute molecule, Nex, to be, N ex(r ) = 4πρv

∫0



[guv(r ) − gvv(r )]r 2 dr

(14)

which measures the excess of solvent molecules over what would be obtained if the solute molecule was a solvent molecule. The Nex of solvent molecules around Na+ and Cl- ion are shown in Figure 5,parts a and b for xMeOH = 0.25, 0.50 and 0.75. In the case of water−methanol mixtures, for xMeOH = 0.25, 0.50, and 0.75 the Nex curve [Figure 5(a)] for Na+−O (H2O) shows positive values beyond 0.26 nm, indicating preferential solvation by water molecules. From Nex of methanol molecules around Na+ (Figure 5a), it is seen that the corresponding value of Nex is smaller than Nex of water molecules around Na+. It is noticed from Figure 5a that there is an excess of water molecules (a deficit of MeOH molecules) around Na+ ion in xMeOH = 0.25, 0.50, and 0.75. Therefore, Na+ is preferentially solvated by water in all methanol mole fractions. From Figure 5b, it is seen that for all methanol mole fractions, the high positive values of Nex beyond 0.34 nm indicate preferential solvation by water molecules. G

DOI: 10.1021/acs.jpcb.5b05401 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B

Ka value of 0.60 M−1, which is lower than what we have obtained (1.10 M−1). F. Hydrogen Bonding Analysis. Appropriate definitions of molecules forming hydrogen bond are usually based on either energetic or geometric criteria. Although there are several examples in the literature which suggests that geometrical criterion (based on only rO−H distance) is insufficient for the calculation of hydrogen bonding under supercritical conditions92−94 but the more general geometric criterion (based on distance and angle both) suggested by Ma et al.95 is found to be as efficient as an energetic criterion for calculation of hydrogen bond in the supercritical condition. In the present study, we have used the more general geometric definition of hydrogen bond. The following geometric criteria was used to detect the presence of hydrogen bonding in pure water,96 i) r(OWater···HWater) distance smaller than 2.74 Å ii) r(OWater···OWater) distance smaller than 3.30 Å iii) ∠H−O−O angle smaller than 30° Between methanol molecules, the geometric criteria were24 (i) r(OMeOH···HMeOH) distance smaller than 2.60 Å (ii) r(OMeOH···OMeOH) distance smaller than 3.50 Å (iii) ∠H−O−O angle smaller than 30° and between water and methanol molecules24 (i) r(OMeOH···HWater) or r(OWater···HMeOH) distance smaller than 3.10 Å (ii) r(OMeOH···OWater) distance smaller than 3.50 Å (iii) ∠H−O−O angle smaller than 30° The hydrogen bonding between methanol and water has also been calculated in both room temperature as well as under supercritical conditions. Hydrogen bonding between water and methanol molecules involves the configurations of the type shown in Figure 6.

Figure 5. Nex around Na+ (a) and Cl− (b) ion in water-MeOH mixtures with xMeOH = 0.25, 0.50, and 0.75 under supercritical conditions at 623 K and 350 bar.

E. Association Constants (Ka). To understand the extent of association of Na+−Cl− ion pair in water−methanol mixtures in ambient and supercritical conditions, we have calculated the association constant (Ka), by integrating the PMF to the first maxima which defines the outer limit of the solute−solute contact configuration. It is calculated using the formula given in our earlier work.62 The values of Ka for different systems (both under ambient and under supercritical conditions) simulated here are presented in Table S9. The uncertainties in the association constants are around 1%. As seen from the association constant values, we can conclude that at a particular temperature, as the methanol mole fraction increases, Ka increases. Association constant values are also consistent with the PMFs; stable CIPs correspond to higher Ka values. Association constants at room temperature range from 1.1 (xMeOH = 0.0) to 4.2 × 103 M−1 (xMeOH = 1.0). Their values under supercritical conditions are much larger, ranging from 7.02 × 102 (xMeOH = 0) to 8.25 × 1010 M−1 (xMeOH = 1.0). So we can say that ion pair association is strongly favored at high temperatures. Hawlicka et al.40,44 has also calculated association constants of ion pairs in their study. They have calculated association constants for concentrated solutions (containing 8 cations, 8 anions and 400 solvent molecules). For pure methanol, they have reported association constant of the order of 84.4 M−1, which is very low compared to those obtained by us. The large difference may be a result of different number of ion pairs as well as different potential model for solvent molecules used in the simulations. In our earlier studies, we have studied the thermodynamics of ion pair association for concentrated solutions.80 We have observed that, stability of contact ion pair decreases with increase in concentration of ion pairs which would result a decrease in the association constants. Hence, the association constants would be larger for a system containing one ion pair compared to those containing more than one pair. For Na+−Cl− ion pair in supercritical water, Chialvo et al.57 have calculated association constants with three different ion−water models (PRH, PR, CHJ). They have predicted the association constant of ion pair in supercritical water to be 4.36 ± 0.04 (PRH), 4.03 ± 0.04 (PR) and 3.71 ± 0.03 (CHJ) (association constants are in log KM a ). In terms of log KM a , we get a value of 2.84 in pure water under supercritical conditions. The difference in association constant may be a result of using different potential models for both solute and solvent as well as different conditions (temperature, pressure and density). Fennell and co-workers91 have studied the ion pairing of alkali halides in aqueous solutions using different potential models for water and ion pairs. For the Na+−Cl− ion pair in pure water using the TIP3P model, they have reported a

Figure 6. Hydrogen bonding between water and methanol molecules.

We have calculated the average number of hydrogen bonds per solvent molecule both in ambient as well as under supercritical conditions at CIPs. The results are given in the Table 4 and Table 5. From the above table, it is seen that under supercritical conditions, the average number of hydrogen bonds per water molecule when xMeOH = 0.00 is 2.32, whereas the same number at room temperature is 3.09. The average number of hydrogen bonds per methanol molecule under supercritical conditions when xMeOH = 1.00 is 0.73 and in room temperature it is 2.21. The average number of hydrogen bonds in pure water and pure methanol in ambient conditions are in very good agreement with the results reported for other models: the modified SPC water ⟨nHB⟩ = 3.5497 and H1 methanol ⟨nHB⟩ = 1.8723 respectively. In the methanol−water mixtures under both ambient and supercritical conditions, ⟨nHB⟩ decreases linearly with increasing xMeOH. Addition of electrolyte does not change this trend, although ⟨nHB⟩ numbers are smaller which indicates H

DOI: 10.1021/acs.jpcb.5b05401 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B

temperature. At 623 K, P = 30 MPa and a density of 0.365 g/ cm3, the average number of hydrogen bonds in the bulk was reported to be 0.50 for supercritical methanol62 which is close enough to what we have obtained for xMeOH = 1.00 in the bulk. We have also calculated the average number of hydrogen bonds per solvent molecule at rNa+Cl− = 0.50 nm (SShIP) and 0.74 nm (SSIP). Details are given in Tables S11−S14 of the Supporting Information. A possible reason for the absence of SShIP under supercritical conditions is due to the decrease in hydrogen bonding in the solvation shells of Na+ and Cl− as we increase xMeOH. There is a steep drop in the average number of hydrogen bonds per water molecule in the solvation shell of Na+ ion in the SShIP (from 1.38 in xMeOH = 0.00 to 0.39 in xMeOH = 0.25) and SSIP (from 1.41 in xMeOH = 0.00 to 0.41 in xMeOH = 0.25) configuration under supercritical conditions. So it is seen that addition of methanol molecules results in a decrease in stability of SShIP as well as SSIP configurations. In fact, SShIPs and SSIPs are absent in pure methanol. However, in ambient conditions, the corresponding decrease in average number of hydrogen bonds is from 1.70 in xMeOH = 0.00 to 1.43 in xMeOH = 0.25 (in SShIP) and from 1.76 in xMeOH = 0.00 to 1.48 in xMeOH = 0.25 (in SSIP). Thus, significant retention of hydrogen bonding in the solvation shells in ambient conditions aids in stabilizing SShIPs and SSIPs. G. Thermodynamic Properties under Supercritical Conditions. The entropy and enthalpy of association has been calculated from the temperature derivative of the PMFs via finite temperature difference method. 99 We have determined the enthalpy and entropy of the ion-pair at interionic distances from 0.2 to 1.0 nm with an interval of 0.01 nm. To calculate entropy and enthalpy of the ion-pair we have performed MD simulations of our system in NPT ensemble with varying temperature. There are small changes in molarities as we change temperature, but since the changes are of the order of 1%, they will not affect our results. Entropy is calculated from the finite temperature difference derivative of the PMF or ΔG(r) at each inter solute separation r,

Table 4. Average Number of Hydrogen Bonds in the Bulk and Solvation Shell Per Solvent Molecule under Supercritical Conditions at CIP average number of H-bonds in the bulk

average number of H-bonds in the solvation shell of Na+

average number of H-bonds in the solvation shell of Cl−

xMeOH

H2O

MeOH

H2O

MeOH

H2O

MeOH

0.00 0.25 0.50 0.75 1.00

2.32 0.89 0.50 0.32 _

_ 0.29 0.46 0.61 0.73

1.36 0.36 0.25 0.21 _

_ 0.18 0.35 0.41 0.52

1.56 0.50 0.39 0.28 _

_ 0.22 0.40 0.50 0.68

Table 5. Average Number of Hydrogen Bonds in the Bulk and Solvation Shell Per Solvent Molecule in Ambient Conditions at CIP average number of H-bonds in the bulk

average number of H-bonds in the solvation shell of Na+

average number of H-bonds in the solvation shell of Cl−

xMeOH

H2O

MeOH

H2O

MeOH

H2O

MeOH

0.00 0.25 0.50 0.75 1.00

3.09 2.13 2.02 1.79 _

_ 0.55 1.10 1.69 2.21

1.68 1.40 1.03 0.90 _

_ 0.16 1.02 1.40 1.56

1.78 1.57 1.23 1.01 _

_ 0.53 1.05 1.55 1.68

that in NaCl solutions some hydrogen bonds are broken. Addition of ionic solutes like NaCl creates a strong local structure around themselves and as a result H-bonding structure of both water and methanol is disrupted under both ambient and supercritical conditions. Decrease in the extent of average number of hydrogen bonds can be explained by weaker intermolecular interactions between solvent molecules resulting from a lower solvent density under supercritical conditions. We have also calculated the average number of hydrogen bonds in the solvation shells of both the cation and the anion. The solvation shell of the anion contains more number of hydrogen bonds than the solvation shell of the cation. In the solvation shell of both the cation and the anion, the average number of hydrogen bonds per water molecule decreases with increasing methanol mole fraction. The average number of hydrogen bonds per water molecule is more than the average number of hydrogen bonds per methanol molecule both in the bulk as well as in the solvation shell of the ions. Hawlicka et al.26 have studied the effect of NaCl on the hydrogen bond network of water−methanol mixtures in ambient conditions. According to their study the increase in methanol concentration causes a decrease in the average number of H-bonds. The average number of hydrogen bonds in the bulk as well as in the solvation shell of ions obtained by us is smaller than that reported by Hawlicka et al.26 The small differences found here can be attributed to different potentials (for solvent molecules) used in the simulations. Cochran et al.54 have found that the average number of hydrogen bonds per water molecule surrounding water molecule, at the dense supercritical state is about 1.0 and near supercritical state, the number is 0.8. They have reported that presence of ionic solutes does not affect the number of hydrogen bonds in the solvation shell of ions which is opposite to what we have obtained. Chalaris et al.98 studied the hydrogen bonding in supercritical methanol as a function of

ΔS(r ) = −

ΔG(r , T + ΔT ) − ΔG(r , T − ΔT ) ΔT

(15)

In the present calculation, values of T and ΔT are chosen to be 623 and 20 K respectively. The enthalpy contribution to the free energy, ΔH(r), can be obtained from entropy ΔS(r) and the PMF ΔG(r) at temperature T. ΔH(r ) = ΔG(r ) + T ΔS(r )

(16)

The temperature-dependent PMFs in all methanol mole fractions is given in Figure S7a−e (supercritical conditions) and Figure S8a−e (ambient conditions) in the Supporting Information. Generally, we see that as temperature increases, the depth of CIP increases and hence it becomes more stable. The increased depth of the CIP can be explained by the reduced dielectric screening at higher temperatures. Devlin et al.100 observed increased stability of contact ion-pair (CIP) with increasing temperature in the spectroscopic study of LiNO3 in DMSO. They suggested that the temperature effect is explained by “progressive partial desolvation” of the ions. An elevation of temperature promotes the desolvation of ions and hence the ion−ion interaction becomes stronger providing stability to CIP. The entropic (−TΔS(r)) and enthalpic (ΔH(r)) contributions to the PMF at 623 K for xMeOH = 0.50 is shown in Figure 7 along with the PMF. The entropic (−TΔS(r)) and enthalpic I

DOI: 10.1021/acs.jpcb.5b05401 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B

Figure 8. Density profiles of (a) MeOH and (b) H2O around the Na+−Cl − ion-pair at CIP in xMeOH = 0.50 (r = 0.26 nm) under supercritical conditions. Figure 7. Entropy and enthalpy contributions to PMFs of Na+−Cl− ion-pair in xMeOH = 0.50 under supercritical conditions.

(ΔH(r)) contributions to the PMF at supercritical conditions for all other mole fractions of methanol are shown in Figure S9a−d, and for ambient conditions these are shown in Figure S10a−e in the Supporting Information along with the PMFs. It is seen from Figure 7 that the CIP state is entirely stabilized by the favorable entropic contribution, enthalpic contribution being highly unfavorable. The entropic contribution at the CIP varies from −262.15 kJ/mol (xMeOH = 0.00) to −452.52 kJ/mol (xMeOH = 1.00). When the rNa+Cl− distance is shorter than 0.4 nm, a steady increase in ΔH(r) and a rapid decrease in −TΔS(r) are observed. As rNa+Cl− is decreased, solvent molecules form highly structured layers around the two solutes are released into the bulk, thereby increasing entropy (or decrease in − TΔS(r) as shown in Figure 7). On the other hand, as a result of expulsion of water molecules from the interionic region, attractive interactions between the solute and solvent in the confined region are lost. This results in a steady increase in enthalpy as the interionic distance is decreased. It is also noticed from Figure 7 that the formation of the CIP is driven by entropy and the process of ion-association is endothermic. The values of entropy, enthalpy and free energy at CIPs, TSs, SShIPs and SSIPs are given in Table S15 (supercritical conditions) and Table S16 (ambient conditions) of the Supporting Information. In every composition, this ionpair formation is endothermic in nature. It is also noticed that TΔS values are larger than ΔH values which suggests that ionpair formation in water−methanol mixtures under supercritical condition is driven by entropy rather than enthalpy. From the composition-dependent TΔS and ΔH values, we note that entropy of the system increases as we increase the methanol mole fraction. H. Density Profiles around the Ion-Pair. The details of solvent distributions around the ion pair can be studied by looking at density profiles of solvent molecules. The planes containing the interionic axis are considered for calculating the mean number densities of solvent sites. The densities of all the planes containing the interionic axis are projected on the x−y plane. The density profiles of water and MeOH around Na+− Cl− in xMeOH = 0.50 at the contact ion-pair (CIP) under both supercritical and ambient conditions are given in Figure 8, parts a and b) and Figure 9, parts a and b). The solvent densities shown in parts a and b of Figures 8 and 9 indicate that the both cation and anion are preferentially solvated by water. It is also observed from Figures 8 and 9 that, under supercritical conditions, both the cation and the anion are preferentially solvated by water and hardly solvated my methanol. In ambient conditions, although both the cation and anion are preferentially solvated by water, the density of

Figure 9. Density profiles of (a) MeOH and (b) H2O around the Na+−Cl − ion-pair at CIP in xMeOH = 0.50 (r = 0.26 nm) under ambient conditions.

methanol around Na+ is more than that under supercritical conditions. So, first solvation shell of Na+ contains both methanol and water molecules in ambient conditions. The density profiles of water and MeOH around Na+- Cl− in xMeOH = 0.50 at the SShIP and SSIP and for xMeOH = 0.25, 0.75 at CIP, SShIP, and SSIP under supercritical and ambient conditions is given in Supporting Information in Figure S11, parts a and b, to Figure S26, parts a and b. It is observed that, under ambient conditions, both the ions have well-defined solvent separated configurations. In the SSIPs, solvation shells are fully formed around both ions. H. Angular Distribution Functions. The angular orientations of solvent molecules around the ions can be characterized by the angle θ between (a) the vector connecting the ion to the oxygen site of solvents and (b) the dipole moment vector of solvent molecules. These vectors are shown in Figure 10. The angular distribution functions (ADFs), P(cos

Figure 10. Representation of ion−solvent vector and the dipole moment vectors of the solvents.

θ), as a function of cos θ are shown in Figure 11a−d for the Na+−Cl− ion pair under supercritical and ambient conditions for xmethanol = 0.25. The ADFs of solvent molecules present in the first coordination shell of Na+ (solid black) and Cl− (solid red) are shown at CIPs in the Figure 11a−d for xmethanol = 0.25. In the case of Cl−, the ADFs (for methanol and water) show a broad peak at cos θ = −0.66 which suggests that the solvent molecules are nearly in an antiparallel dipole orientation in the coordination shell of Cl−. The ADF peak of water and methanol for Cl− in ambient conditions are more compact J

DOI: 10.1021/acs.jpcb.5b05401 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B

IV. CONCLUSIONS We have used the constrained molecular dynamics technique to simulate solutions of the sodium chloride ion-pair in supercritical water−methanol mixtures of different compositions. The PMFs of a dilute (∼0.025 M) Na+−Cl− ion-pair under supercritical conditions exhibits some remarkable differences in comparison with that in ambient conditions. First, the interaction between two ions is much stronger under supercritical conditions than in ambient conditions. The derived potentials of mean force show the presence of both a CIP and an SShIP in xMeOH = 0.0 and 0.25 under supercritical conditions. The stability of CIP increases with increase in the mole fraction of methanol. The PMFs for ambient conditions have also been calculated and they show the presence of CIPs and SShIPs as well SSIPs. Under supercritical conditions, PMFs exhibit an inflection point near 0.47 nm in xMeOH = 0.50, 0.75, and 1.00. All the PMFs are confirmed by ion pair distance residence times. In all the mixtures considered in the present work, the local solvation shells around the ions are strong and do not show any changes in the RDF peak positions with changes in mixture compositions. It is shown that, under supercritical conditions, the peak heights become much larger due to the well-known local density inhomogeneity (local clustering) effect in SCFs. The height of the Na+−O (water) RDF peak is higher than that of the Na+−O (methanol) peak. The same trend is observed in the case of the Cl−−O peak. This indicates that the local density of solvent sites around ions changes very slowly in comparison with the bulk density of solvents. The preferential solvation analysis shows that both the cation and the anion preferentially interact with water molecules. The entropic and enthalpic contributions to the PMFs are calculated from the temperature dependence of the PMFs. The stabilization of the contact pair state is mainly due to an increase in entropy arising from the expulsion of the highly structured solvent molecules from the intersolute region. The unfavorable solvent induced contribution to the enthalpy change at the CIP state arises from the enhancements in solvent−solvent repulsions. It is observed that ion-pair association in water−methanol mixtures under supercritical conditions is endothermic and entropically driven. Even in ambient conditions, ion pair association is favored by entropy to a much smaller extent than under supercritical conditions.

Figure 11. ADFs of solvent molecules for Na+−Cl− ion pair under supercritical and ambient conditions for xMeOH = 0.25, (a) Na+− methanol and Cl−−methanol under supercritical conditions, (b) Na+− methanol and Cl−−methanol under ambient conditions, (c) Na+− water and Cl−−water under supercritical conditions and (d) Na+− water and Cl−−water under ambient conditions. The angle θ is defined in Figure 10.

compared to supercritical conditions. Therefore, solvent molecules are better oriented around Cl− in ambient conditions compared to supercritical conditions. As we go from CIP to SSIP, the ADF peaks of solvent molecules for Cl− become sharper. The peak height of water for Cl− (in all CIP, SShIP and SSIP configurations) is more than that of methanol. It is seen that water molecules are better oriented around Cl− compared to methanol molecules. Distribution functions of the angular orientation of solvent molecules in the primary shells of Na+ ion show the peak centered at cos θ = 1.0. This indicates that the dipole moment of solvent molecules (water and methanol) is in the same direction as the ion-oxygen vector. There is an enhancement in the orientations by about 20% at cos θ = 1.0 as we go from supercritical to ambient conditions. This peak height gets further enhanced as we go from CIP to SSIP. P(cos θ) increases proportionately with increase in methanol mole fraction for both the cation and the anion in all CIP, SShIP and SSIP configurations (both in ambient and supercritical conditions). The ADFs for CIP configurations for xmethanol = 0.50 and 0.75 are given in Supporting Information in Figures S27a−d and S28a−d, respectively. The ADFs for SShIP and SSIP configurations for xmethanol = 0.25, 0.50, and 0.75 are given in Supporting Information in Figures S29a−d to S34a−d. Angular distribution functions of Na+, Cl− and I− ions in water−methanol mixtures have been studied by Hawlicka et al.40,44 in ambient conditions. They observed a sharp peak in the ADF at cos θ = 1.0 for the Na+ ion. In their study, they observed a narrow peak for methanol molecules around Na+ [P(cos θ) = 0.14] compared to water [P(cos θ) = 0.045]. But in our study, we observe the peak for water molecules around Na+ [P(cos θ) = 0.15] to be narrower than that for methanol [P(cos θ) = 0.021]. For the Cl− ion, in both aqueous and methanolic solutions, they observed a broad peak at cos θ = −0.65 which is similar to what we have obtained.



ASSOCIATED CONTENT

* Supporting Information S

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.jpcb.5b05401. Force field and geometric parameters of solvent and solute molecules (Tables S1 and S2), critical parameters of water−methanol mixtures (Table S3), comparison of densities of water−methanol mixtures obtained from MD simulations with experimental critical densities (Table S4), dielectric constants at supercritical and ambient conditions (Tables S5 and S6), characteristics of PMFs at supercritical and ambient conditions (Tables S7 and S8), association constants (Table S9), dissociation constants (Table S10), hydrogen bonding at SShIP and SSIP at supercritical and ambient conditions (Tables S11−S14), characteristics of thermodynamic properties at supercritical and ambient conditions (Table S15 and S16), PMFs as a function of molality (Figure S1), ion pair distance residence times (Figure S2), radial K

DOI: 10.1021/acs.jpcb.5b05401 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B



(12) Hawlicka, E. Self-Diffusion in Water − Alcohol Mixtures. Part 1. Water − Methanol Solution of NaI. Ber. Bunsenges. Phys. Chem. 1983, 87, 425−428. (13) Versmold, H. NMR Studies of the Reorientational Motion in Methanol and Methanol − Water Mixtures. Ber. Bunsenges. Phys. Chem. 1980, 84, 168−173. (14) Derlacki, Z. J.; Easteal, A. J.; Edge, A. V. J.; Woolf, L. A.; Roksandic, Z. Diffusion Coefficients of Methanol and Water and the Mutual Diffusion Coefficient in Methanol - Water Solutions at 278 and 298 K. J. Phys. Chem. 1985, 89, 5318−5322. (15) Ferrario, M.; Haughney, M.; McDonald, I. R.; Klein, M. L. Molecular - Dynamics Simulation of Aqueous Mixtures: Methanol, Acetone, and Ammonia. J. Chem. Phys. 1990, 93, 5156−5166. (16) Tanaka, H.; Gubbins, K. E. Structure and Thermodynamic Properties of Water−Methanol Mixtures: Role of the Water − Water Interaction. J. Chem. Phys. 1992, 97, 2626−2634. (17) Wensink, E. J. W.; Hoffmann, A. C.; van Maaren, P. J.; van der Spoel, D. Dynamic Properties of Water/Alcohol Mixtures Studied by Computer Simulation. J. Chem. Phys. 2003, 119, 7308−7317. (18) Ono, T.; Horikawa, K.; Ota, M.; Sato, Y.; Inomata, H. Insight into the Local Composition of the Wilson Equation at High Temperatures and Pressures through Molecular Simulations of Methanol−Water Mixtures. J. Chem. Eng. Data 2014, 59, 1024−1030. (19) Rahman, A.; Stillinger, F. H. Molecular Dynamics Study of Liquid Water. J. Chem. Phys. 1971, 55, 3336−3359. (20) Zielkiewicz, J. Structural Properties of Water: Comparison of the SPC, SPCE, TIP4P, and TIP5P Models of Water. J. Chem. Phys. 2005, 123, 104501. (21) Jorgensen, W. L. Structure and Properties of Liquid Methanol. J. Am. Chem. Soc. 1980, 102, 543−549. (22) Jorgensen, W. L.; Ibrahim, M. Pressure Dependence of Hydrogen Bonding in Liquid Methanol. J. Am. Chem. Soc. 1982, 104, 373−378. (23) Haughney, M.; Ferrario, M.; McDonald, I. R. MolecularDynamics Simulation of Liquid Methanol. J. Phys. Chem. 1987, 91, 4934−4940. (24) Pagliai, M.; Cardini, G.; Righini, R.; Schettino, V. Hydrogen Bond Dynamics in Liquid Methanol. J. Chem. Phys. 2003, 119, 6655− 6662. (25) Batista da Silva, J. A.; Moreira, F. G. B.; Leite dos Santos, V. M.; Longo, R. L. On the Hydrogen Bond Networks in the Water− Methanol Mixtures: Topology, Percolation and Small-World. Phys. Chem. Chem. Phys. 2011, 13, 6452−6461. (26) Hawlicka, E.; Swiatla-Wojcik, D. Molecular Dynamics Studies of NaCl Solutions in Methanol−Water Mixtures. An Effect of NaCl on Hydrogen Bonded Network. Chem. Phys. 1998, 232, 361−369. (27) Duffie, J. A.; Beckman, W. A. Solar Engineering of Thermal Processes; Wiley: New York, 1992. (28) Aliev, M. M.; Magee, J. W.; Abdulagatov, I. M. PVTx and Isochoric Heat Capacity Measurements for Aqueous Methanol Solutions. Int. J. Thermophys. 2003, 24, 1551−1579. (29) Antal, M. J., Jr.; Allen, S. G.; Schulman, D.; Xu, X.; Divilio, R. J. Biomass Gasification in Supercritical Water. Ind. Eng. Chem. Res. 2000, 39, 4040−4053. (30) Saka, S.; Kusdiana, D. Biodiesel Fuel from Rapeseed Oil as prepared in Supercritical Methanol. Fuel 2001, 80, 225−231. (31) Genta, M.; Iwaya, T.; Sasaki, M.; Goto, M. Supercritical Methanol for Polyethylene Terephthalate Depolymerization: Observation using Simulator. Waste Manage. 2007, 27, 1167−1177. (32) Dixit, S.; Crain, J.; Poon, W. C. K.; Finney, J. L.; Soper, A. K. Molecular Segregation observed in Concentrated Alcohol − Water Solution. Nature 2002, 416, 829−832. (33) Ikushima, Y.; Hatakeda, K.; Saito, N.; Arai, M. An in situ Raman Spectroscopy Study of Subcritical and Supercritical water: The Peculiarity of Hydrogen Bonding near the Critical Point. J. Chem. Phys. 1998, 108, 5855−5860. (34) Owczarek, E.; Rybicki, M.; Hawlicka, E. Influence of CaCl2 on Micro Heterogeneity of Methanol−Water Mixtures. Chem. Phys. 2009, 363, 78−87.

distribution functions of cation and anion with hydrogen site of water and methanol (Figures S3 and Figure S4), excess coordination numbers (Figures S5 and S6), PMFs as a function of temperature at supercritical and ambient conditions (Figures S7 and S8), decomposition curve of PMFs at supercritical conditions (Figure S9) and room temperature (Figure S10), density profiles of solvent molecules around solutes in both supercritical and ambient conditions (Figures S11−S26), and angular distribution functions (Figures S27−S34) (PDF)

AUTHOR INFORMATION

Corresponding Author

*Telephone: +91-22-2576-4199. Fax: +91-22-2576-7152 Email: [email protected]. Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS We would like to thank high computing facility of Indian Institute of Technology Bombay and Chemistry Department of IIT Bombay. We would also like to thank Mayank Kumar Dixit for his generous help in calculating density profiles and angular distribution functions.



REFERENCES

(1) Eckert, C. A.; Knutson, B. L.; Debenedetti, P. G. Supercritical Fluids as Solvents for Chemical and Materials Processing. Nature 1996, 383, 313−318. (2) Akgerman, A. In Supercritical fluids: Extraction and Pollution Prevention; Abraham, M. A., Sunol, A. K., Eds.; ACS Symposium Series; American Chemical Society: Washington, DC, 1997; Chapter 15, pp 208−231. (3) Goldman, S.; Gray, C. G.; Li, W.; Tomberli, B.; Joslin, C. G. Predicting Solubilities in Supercritical Fluids. J. Phys. Chem. 1996, 100, 7246−7249. (4) Lira, C. T. In Supercritical Fluid Extraction and Chromatography: Techniques and Applications; Charpentier, B. A., Sevenants, M. R., Eds.; ACS Symposium Series; American Chemical Society: Washington, DC, 1988; Chapter 1, pp 1−25. (5) Chialvo, A. A.; Cummings, P. T. Solute Induced Effects on the Structure and Thermodynamics of Infinitely Dilute Mixtures. AIChE J. 1994, 40, 1558−1573. (6) Chialvo, A. A.; Kusalik, P. G.; Cummings, P. T.; Simonson, J. M.; Mesmer, R. E. Molecular Approach to High Temperature Solvation. Formal, Integral Equation and Experimental Results. J. Phys.: Condens. Matter 2000, 12, 3585−3593. (7) Chialvo, A. A. In Fluctuation Theory of Solutions: Applications in Chemistry, Chemical Engineering and Biophysics; Matteoli, E., O’Connell, J. P., Smith, P. E., Eds.; CRC Press: Boca Raton, FL, 2013, Chapter 8, pp 191−224. (8) Niessen, H. G.; Woelk, K. Investigations in Supercritical Fluids. Top Curr. Chem. 2007, 276, 69−110. (9) Hu, J.; Duan, Z.; Zhu, C.; Chou, I-M. PVTx Properties of the CO2−H2O and CO2−H2O−NaCl Systems below 647 K: Assessment of Experimental Data and Thermodynamic Models. Chem. Geol. 2007, 238, 249−267. (10) Woolf, L. A. Insights into Solute - Solute - Solvent Interactions from Transport Property Measurements with Particular Reference to Methanol - Water Mixtures and their Constituents. Pure Appl. Chem. 1985, 57, 1083−1090. (11) Roux, A. H.; Desnoyers, J. E. Association Models for Alcohol − Water Mixtures. Proc. Indian Acad. Sci. (Chem. Sci.) 1987, 98, 435− 451. L

DOI: 10.1021/acs.jpcb.5b05401 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B (35) Hawlicka, E.; Swiatla-Wojcik, D. Molecular Dynamics Studies on the Structure of Methanol-Water Solutions of NaCl. Chem. Phys. 1995, 195, 221−233. (36) Hawlicka, E.; Swiatla-Wojcik, D. Molecular Dynamics Simulation of NaCl Solutions in Methanol-Water Mixtures. Intramolecular Vibrations of the Solvent Components. Chem. Phys. 1997, 218, 49−55. (37) Hawlicka, E.; Swiatla-Wojcik, D. Solvation of Ions in Binary Solvents − Experimental and MD Simulation Studies. J. Mol. Liq. 1998, 78, 7−18. (38) Hawlicka, E.; Swiatla-Wojcik, D. Dynamic Properties of NaCl − Methanol − Water Systems − MD Simulation Studies. Phys. Chem. Chem. Phys. 2000, 2, 3175−3180. (39) Hawlicka, E.; Swiatla-Wojcik, D. Effect of Sodium Halides on the Structure of Methanol − Water mixtures − MD Simulation Studies. J. Mol. Liq. 2002, 98-99, 357−367. (40) Hawlicka, E.; Swiatla-Wojcik, D. MD Simulation Studies of Selective Solvation in Methanol-Water Mixtures: An Effect of the Charge Density of a Solute. J. Phys. Chem. A 2002, 106, 1336−1345. (41) Bopp, P.; Jancso, G.; Heinzinger, K. An Improved Potential for Non − Rigid Water Molecules in the Liquid Phase. Chem. Phys. Lett. 1983, 98, 129−133. (42) Jancso, G.; Bopp, P.; Heinzinger, K. Molecular Dynamics Study of High Density Liquid Water Using a Modified Central-Force Potential. Chem. Phys. 1984, 85, 377−387. (43) Pálinkás, G.; Hawlicka, E.; Heinzinger, K. A Molecular Dynamics Study of Liquid Methanol with a Flexible Three-Site Model. J. Phys. Chem. 1987, 91, 4334−4341. (44) Hawlicka, E.; Swiatla-Wojcik, D. Aggregation of Ions in Methanol−Water Solutions of Sodium Halides. J. Chem. Phys. 2003, 119, 2206−2213. (45) Day, T. J. F.; Patey, G. N. Ion Solvation Dynamics in Water− Methanol and Water− Dimethyl Sulfoxide Mixtures. J. Chem. Phys. 1999, 110, 10937−10944. (46) Feakins, D.; Watson, P. Studies in Ion Solvation in NonAqueous Solvents and their Aqueous Mixtures. Part I. Properties of Ion Constituents. J. Chem. Soc. 1963, 905, 4734−4741. (47) Bujnicka, K.; Hawlicka, E. Solvation of Ca2+ in Aqueous Methanol-MD Simulation Studies. J. Mol. Liq. 2006, 125, 151−157. (48) Rybicki, M.; Hawlicka, E. Solvation of Mg2+ Ions in Methanol− Water Mixtures: Molecular Dynamics Simulation. Chem. Phys. 2012, 400, 79−85. (49) Pálinkás, G.; Hawlicka, E.; Heinzinger, K. Molecular Dynamics Simulations of Water − Methanol Mixtures. Chem. Phys. 1991, 158, 65−76. (50) Covington, A. K.; Dunn, M. Nuclear Magnetic Resonance Study of Preferential Solvation. J. Chem. Soc., Faraday Trans. 1 1989, 85, 2827−2834. (51) Holz, M.; Weingärtner, H.; Hertz, H. G. Nuclear Magnetic Relaxation of Alkali Halide Nuclei and Preferential Solvation in Methanol + Water Mixtures. J. Chem. Soc., Faraday Trans. 1 1977, 73, 71−83. (52) Hawlicka, E. Self-Diffusion of Sodium, Chloride and Iodide Ions in Methanol-Water Mixture. Z. Naturforsch., A: Phys. Sci. 1986, 41a, 939−943. (53) Chialvo, A. A.; Cummings, P. T.; Simonson, J. M.; Mesmer, R. E. Thermodynamics and Kinetics of Ion Speciation in Supercritical Aqueous Solutions: A Molecular-based Study. Fluid Phase Equilib. 1998, 150−151, 107−115. (54) Cochran, H. D.; Cummings, P. T.; Karaborni, S. Solvation in Supercritical Water. Fluid Phase Equilib. 1992, 71, 1−16. (55) Cummings, P. T.; Cochran, H. D.; Simonson, J. M.; Mesmer, R. E.; Karaborni, S. Simulation of Supercritical Water and of Supercritical Aqueous Solutions. J. Chem. Phys. 1991, 94, 5606−5621. (56) Chialvo, A. A.; Cummings, P. T.; Simonson, J. M.; Mesmer, R. E. Molecular Simulation Study of Speciation in Supercritical Aqueous NaCl Solutions. J. Mol. Liq. 1997, 73-74, 361−372.

(57) Chialvo, A. A.; Cummings, P. T.; Cochran, H. D.; Simonson, J. M.; Mesmer, R. E. Na+−Cl− Ion Pair Association in Supercritical Water. J. Chem. Phys. 1995, 103, 9379−9387. (58) Cummings, P. T.; Chialvo, A. A.; Cochran, H. D. Molecular Simulation Study of Solvation Structure in Supercritical Aqueous Solutions. Chem. Eng. Sci. 1994, 49, 2735−2748. (59) Lee, S. H.; Cummings, P. T.; Simonson, J. M.; Mesmer, R. E. Molecular Dynamics Simulation of the Limiting Conductance of NaCl in Supercritical Water. Chem. Phys. Lett. 1998, 293, 289−294. (60) Lee, S. H.; Cummings, P. T. Molecular Dynamics Simulation of Limiting Conductances for LiCl, NaBr, and CsBr in Supercritical Water. J. Chem. Phys. 2000, 112, 864−869. (61) Chialvo, A. A.; Cummings, P. T. Hydrogen Bonding in Supercritical Water. J. Chem. Phys. 1994, 101, 4466−4469. (62) Keshri, S.; Sarkar, A.; Tembe, B. L. Molecular Dynamics Simulation of Na+−Cl− Ion Pair in Supercritical Methanol. J. Supercrit. Fluids 2015, 103, 61−69. (63) Hess, B.; Kutzner, C.; van der Spoel, D.; Lindahl, E. GROMACS 4: Algorithms for Highly Efficient, Load - Balanced, and Scalable Molecular Simulation. J. Chem. Theory Comput. 2008, 4, 435−447. (64) Jorgensen, W. L.; Chandrasekhar, J.; Madura, J. D.; Impey, R. W.; Klein, M. L. Comparison of Simple Potential Functions for Simulating Liquid Water. J. Chem. Phys. 1983, 79, 926−935. (65) Jorgensen, W. L.; Maxwell, D. S.; Tirado-Rives, J. Development and Testing of the OPLS All-Atom Force Field on Conformational Energetics and Properties of Organic Liquids. J. Am. Chem. Soc. 1996, 118, 11225−11236. (66) Martínez, L.; Andrade, R.; Birgin, E. G.; Martínez, J. M. Packmol: A Package for Building Initial Configurations for Molecular Dynamics Simulations. J. Comput. Chem. 2009, 30, 2157−2164. (67) Allen, M. P.; Tildesley, D. J. Computer Simulations of Liquids; Claredon: Oxford, U.K., 1987. (68) Marshall, W. L.; Jones, E. V. Liquid Vapor Critical Temperatures of Several Aqueous Organic and Organic-Organic Solutions Systems. J. Inorg. Nucl. Chem. 1974, 36, 2319−2323. (69) Bulemela, E.; Tremaine, P.; Ikawa, S. Volumetric Behavior of Water - Methanol Mixtures in the Vicinity of Critical Region. Fluid Phase Equilib. 2006, 245, 125−133. (70) Darden, T.; York, D.; Pedersen, L. Particle mesh Ewald: An N· log(N) Method for Ewald Sums in Large Systems. J. Chem. Phys. 1993, 98, 10089−10092. (71) van Gunsteren, W. F.; Berendsen, H. J. C. Computer Simulation of Molecular Dynamics: Methodology, Applications, and Perspectives in Chemistry. Angew. Chem., Int. Ed. Engl. 1990, 29, 992−1023. (72) Ryckaert, J. P.; Ciccotti, G.; Berendsen, H. J. C. Numerical Integration of the Cartesian Equations of Motion of a System with Constraints: Molecular Dynamics of n-Alkanes. J. Comput. Phys. 1977, 23, 327−341. (73) Bussi, G.; Donadio, D.; Parrinello, M. Canonical Sampling through Velocity Rescaling. J. Chem. Phys. 2007, 126, 014101−014107. (74) Berendsen, H. J. C.; Postma, J. P. M.; DiNola, A.; Haak, J. R.; van Gunsteren, W. F. Molecular Dynamics with Coupling to an External Bath. J. Chem. Phys. 1984, 81, 3684−3690. (75) van Gunsteren, W. F.; Berendsen, H. J. C. A Leap-frog Algorithm for Stochastic Dynamics. Mol. Simul. 1988, 1, 173−185. (76) Parrinello, M.; Rahman, A. Polymorphic Transitions in Single Crystals: A New Molecular Dynamics Method. J. Appl. Phys. 1981, 52, 7182−7190. (77) Trzesniak, D.; Kunz, A. P. E.; van Gunsteren, W. A Comparison of Methods to Compute the Potential of Mean Force. ChemPhysChem 2007, 8, 162−169. (78) Hess, B.; Holm, C.; van der Vegt, N. Osmotic Coefficients for Atomistic NaCl(aq) Force Fields. J. Chem. Phys. 2006, 124, 164509. (79) Neumann, R. M. Entropic Approach to Brownian movement. Am. J. Phys. 1980, 48 (5), 354−357. (80) Patil, U. N.; Keshri, S.; Tembe, B. L. Thermodynamics of Na+ Cl− Ion − Pair Association in Acetonitrile−Dimethyl Sulfoxide Mixtures. Chem. Phys. Lett. 2015, 620, 134−138. M

DOI: 10.1021/acs.jpcb.5b05401 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B (81) Ono, T.; Kobayashi, S.; Hoshina, T.; Sato, Y.; Inomata, H. Volumetric Behavior and Solution Microstructure of Methanol−Water Mixture in Sub- and Supercritical State via Density Measurement and MD Simulation. Fluid Phase Equilib. 2011, 302, 55−59. (82) Arce, A.; Blanco, A.; Soto, A.; Vidal, I. Densities, Refractive Indices, and Excess Molar Volumes of the Ternary Systems Water + Methanol + 1-Octanol and Water + Ethanol + 1-Octanol and their Binary Mixtures at 298.15 K. J. Chem. Eng. Data 1993, 38, 336−340. (83) Marcus, Y.; Hefter, G. Ion Pairing. Chem. Rev. 2006, 106, 4585− 4621. (84) Wang, P.; Anderko, A. Computation of Dielectric Constants of Solvent Mixtures and Electrolyte Solutions. Fluid Phase Equilib. 2001, 186, 103−122. (85) Yui, K.; Sakuma, M.; Funazukuri, T. Molecular Dynamics Simulations on Ion Pair Association of NaCl from Ambient to Supercritical Water. Fluid Phase Equilib. 2010, 297, 227−235. (86) Dixit, M. K.; Tembe, B. L. Potentials of Mean Force of Sodium Chloride Ion Pair in Dimethyl Sulfoxide − Methanol Mixtures. J. Mol. Liq. 2013, 178, 78−83. (87) Megyes, T.; Grósz, T.; Radnai, T.; Bakó, I.; Pálinkás, G. Solvation of Calcium Ion in Polar Solvents: An X-ray Diffraction and Ab Initio Study. J. Phys. Chem. A 2004, 108, 7261−7271. (88) Neilson, G. W.; Enderby, J. E. Neutron and X-Ray Diffraction Studies of Concentrated Aqueous Electrolyte Solutions. Annu. Rep. Prog. Chem., Sect. C: Phys. Chem. 1979, 76, 185−220. (89) Yu, H.; Whitfield, T. W.; Harder, E.; Lamoureux, G.; Vorobyov, I.; Anisimov, V. M.; MacKerell, A. D.; Roux, B. Simulating Monovalent and Divalent Ions in Aqueous Solution using a Drude Polarizable Force Field. J. Chem. Theory Comput. 2010, 6, 774−786. (90) Kirkwood, J. G.; Buff, F. P. The Statistical Mechanical Theory of Solutions. 1. J. Chem. Phys. 1951, 19, 774−777. (91) Fennell, C. J.; Bizjak, A.; Vlachy, V.; Dill, K. A. Ion Pairing in Molecular Simulations of Aqueous Alkali Halide Solutions. J. Phys. Chem. B 2009, 113, 6782−6791. (92) Swiatla-Wojcik, D. Evaluation of the Criteria of Hydrogen Bonding in Highly Associated Liquids. Chem. Phys. 2007, 342, 260− 266. (93) Kalinichev, A. G.; Bass, J. D. Hydrogen Bonding in Supercritical Water: A Monte Carlo Simulation. Chem. Phys. Lett. 1994, 231, 301− 307. (94) Kalinichev, A. G.; Bass, J. D. Hydrogen Bonding in Supercritical Water. 2. Computer Simulations. J. Phys. Chem. A 1997, 101, 9720− 9727. (95) Ma, H.; Ma, J. Density Dependence of Hydrogen Bonding and Translational-Orientational Structure Order in Supercritical Water: A Molecular Dynamics Study. J. Chem. Phys. 2011, 135, 054504. (96) Mezei, M.; Beveridge, D. L. Theoretical Studies of Hydrogen Bonding in Liquid Water and in Dilute Aqueous Solutions. J. Chem. Phys. 1981, 74, 622−632. (97) Martí, J.; Padro, J. A.; Guárdia, E. Molecular Dynamics Simulation of Liquid Water along the Coexistence Curve: Hydrogen Bonds and Vibrational Spectra. J. Chem. Phys. 1996, 105, 639−649. (98) Chalaris, M.; Samios, J. Hydrogen Bonding in Supercritical Methanol: A Molecular Dynamics Investigation. J. Phys. Chem. B 1999, 103, 1161−1166. (99) Pettitt, B. M.; Rossky, P. J. Alkali Halides in Water: Ion−Solvent Correlations and Ion−Ion Potentials of Mean Force at Infinite Dilution. J. Chem. Phys. 1986, 84, 5836−5844. (100) Wooldridge, R.; Fisher, M.; Ritzhaupt, G.; Devlin, J. P. Dependence of Ion Pairing and Solvation on Solution Temperature: Evidence from DMSO Solution and Matrix-Solvation Spectra. J. Chem. Phys. 1987, 86, 4391−4395.

N

DOI: 10.1021/acs.jpcb.5b05401 J. Phys. Chem. B XXXX, XXX, XXX−XXX