Solvation of halogen ions in aqueous solutions at 500 K–600 K under 100 atm
Shen Haoa),b), Hao Tinga),b), Zhang Feng-Shou†a),b),c)
The Key Laboratory of Beam Technology and Material Modification of the Ministry of Education, College of Nuclear Science and Technology, Beijing Normal University, Beijing 100875, China
Beijing Radiation Center, Beijing 100875, China
Center of Theoretical Nuclear Physics, National Laboratory of the Heavy Ion Accelerator of Lanzhou, Lanzhou 730000, China

Corresponding author. E-mail: fszhang@bnu.edu.cn

*Project supported by the National Natural Science Foundation of China (Grand Nos. 11025524 and 11161130520), the National Basic Research Program of China (Grant No. 2010CB832903), and the European Commission’s 7th Framework Programme (FP7-PEOPLE-2010-IRSES) (Grant No. 269131).

Abstract

Structural properties of the pure water and halogen solutions at high temperatures and pressures are studied by using the molecular dynamics simulations and quantum molecular simulations. The related characters are calculated as functions of temperature and pressure. The results show that the hydrogen bonded networks become looser as temperature increases, with the collapse of the traditional tetrahedral structure. It is similar to the concentration-dependent collapse in the NaCl solutions. However, adding other halogen elements has no further effects on the already weakly bonded water molecules. At the phase changing points, the process of hydration is evident for the bigger ions, so that the bigger the ion is, the smaller a cluster is formed.

PACS: 36.20.Ey; 36.40.Sx; 87.10.Tf
Keyword: water solution; molecular dynamics simulation; quantum chemical simulation; high temperature and pressure
1. Introduction

The hydration processes of biological molecules, the formation of ion clusters, and the chemical reactions in aqueous media, especially for the solutions at high temperatures and pressures, are all intimately associated with water.[17] It is generally considered that the key roles in both biophysical and biochemical processes of water, either the solvation of polar and ionic substances or the repulsion of hydrophobic materials, are all determined by the hydrogen bonding ability. Therefore, the studies on the dynamical characters and microstructures of ionic solvation are very important and intriguing. Many experiments and theoretical studies have been reported on the physical and chemical natures of the hydration processes.[810] A series of water models and the corresponding functions of phase transitions are also proposed.[11, 12] However, the mechanisms behind their behaviors and functions are still poorly investigated.

As a solvent, its own unique nature and characteristics affect the solvation of ions and biological molecules.[13, 14] Corresponding researches on the hydration of DNA, proteins have also been revealed in detail.[1519] Compared to other simple organic liquids, water has more complex phase diagrams, in that a number of triple points and one, probably two, critical points can be found.[2023] Critical water, [2426] near the critical point of the phase diagram, is such a special fluid that many strange properties, such as the decrease of dielectric constant, the special density change, and the coexisting phases, are all observed.[27] A slight change of the temperature or pressure may shift the phase of the system from one state to another. Therefore, the relationship between ions and water molecules at extreme conditions is of theoretical and experimental importance. Since the halogen elements are one of the most important singly charged anions, most of the experiments and theoretical studies of their solutions are focused on the vapor and gas– liquid mixture phase. A thorough investigation of the changes of water structure and hydration processes of halogen elements at the phase changing points is desperately needed.

Molecular dynamics (MD) simulation, one of the classical computation methods, [28] allows one to relate macroscopic to microscopic phenomena. In our present work, the hydration processes and cluster formation have been extensively studied by using the MD simulation, [29, 30] especially for its temperature-dependent properties.[31] The anomalous properties of water are also partly illustrated.[32, 33]

2. Model and method

In our MD simulations the SPC water model is used.[34] 256 water molecules and 8 pairs of ions, NaF, NaCl, NaBr or NaI, are contained in the system, respectively. The corresponding concentration is 1.74  M, which is smaller than the corresponding saturation concentrations. The ensemble for the systems is NPT, revealed by using the Noosé – Hoover method. The time step for the simulations is 0.1  fs, with the simulation time of 5  ns. The pressure for the simulations is 100  atm (1  atm = 1.01325× 105  Pa), which is the same order of magnitude with the critical point: 220  atm.

Table  1 lists the parameters for the SPC water model and ions. The parameters for ions are determined by Smith and Dang.[35] For different particles, the parameters were treated by σ mn = (σ m + σ n)/2 and

Table 1. Parameters of the SPC water model and ions in the MD simulations.

The intermolecular energies in the MD simulations for all particles are shown as,

in which

is the electrostatic interaction, where q is the charge of the atom, ε 0 the vacuum permittivity, ε the relative permittivity, and the rmn the distance between atom m and n. The other one,

is the Lennard-Jones interaction, where ε and σ determine the energy and distance scale, respectively.

In our MD simulations, the parameters used are empirical and whether it is possible to correctly describe the hydration processes at extreme conditions is still under consideration. Since the development of Density Function Theory (DFT) and its application in current computer simulations, the quantum molecular (QM) method is developed as an instrument to explore the properties of liquids.[36, 37] Therefore, comparative analysis between the MD simulations and QM simulations would judge the accuracy of the MD simulations and, further, shed light on the improvement of force field parameters for better describing the bulk water solutions.[38] In the QM simulations, the CP2K code is employed.[39] The self-consistent-charge density-functional tight-binding method (SCCDFTB) is used for the systems.[40, 41] The exchange– correlation energy function is expanded to the second order via using the SCCDFTB method. The reference density is calculated by . The energy in the QM calculations is shown as:

3. Results and discussion
3.1. Validation of the classical parameters

Comparisons of the results from the QM and MD simulations are shown in Figs.  1(a) and 1(b). The inset of Fig.  1(a) is the density against the simulation time, averaged over the adjacent steps in the QM simulations (at 298  K under 1  atm). After 8.2  ps the system reaches equilibrium and the density of the system is basically around 0.96  g/cm3, which is closer to that of the real water, 1.00  g/cm3. This result is better than that of the density function theory together with BLYP or PBE (the density of the system is 0.75  g/cm3 for BLYP and 0.88  g/cm3 for PBE).[42]

Fig.  1. RDFs of O– O, O– H, H– H in the QM and MD simulations at room temperature and pressure (a) and at critical point (b). The density in the QM simulations is shown in the inset of panel  (a).

The Radial Distribution Functions (RDFs) of the QM simulations are derived from the trajectories between 8.2  ps and 12  ps. As shown in Fig.  1(a), the position of the RDF peak for O– H is about 0.1  Å larger than that in the MD simulation, in accordance with the larger O– O distance (O– H· · · O) and the lower density of around 0.96  g/cm3. This is because that parts of the long-term weak interactions are hard to be included in the QM method.[39, 40] However, for the hydrogen bond, the second peak of g(r)OH (H– O· · · H), the results in both the QM and MD simulations are nearly the same.

For a system at high temperature and pressure: 220  atm and 680  K, the RDFs via using the two methods are almost coincidence as shown in Fig.  1(b). The density of the systems is changed from the initial 1.00  g/cm3 to just 0.26  g/cm3. The RDF peaks are also lower and broader compared to that at the room temperature and pressure. Furthermore, the second peaks of g(r)OO almost disappear, indicating that the hydrogen bonds between water molecules become looser.

Actually, the dispersion part of long-term weak interaction is not significant at high temperatures and pressures. This is the main reason for the special properties such as a low dielectric constant, and large ions clusters in these specific conditions[43, 44] and so do almost overlapped curves. Through the comparative analysis, the corresponding distributions at high temperature and pressure are descried consistently and the two methods are both suitable for mimicking the bulk water at high temperatures and pressures though the cost is expensive for the former one.

3.2. Effects of added salts and increased temperature

RDFs of O– O in pure water and NaCl solutions are calculated and shown in Fig.  2, respectively. The value of the first peak in Fig.  2(a) decreases from 3.1 to 2.5 monotonously when the temperature increases from 298  K to 400  K, with the decreased coordinated molecules from 3.25 to 3.08. Consequently, the decrease of the surrounding coordinated water molecules results in the diminish of the number of the hydrogen bonds. This is in agreement with the decreased second RDF peaks of g(r)OO in Fig.  2(a) and g(r)OH in Fig.  1(a), respectively.

Fig.  2. RDFs of O– O in pure water (a) and NaCl solutions (b).

Figure  2(b) is the g(r)OO of water molecules in sodium chloride solutions under different concentrations at 298  K under 1  atm. For the concentration of 1.74  M, the RDF peak at 4.5  Å almost disappears. Regarding the concentration of 2.38  M or 2.82  M, the reduction of the second peak is more significant. The position of the second peak in pure water solution is approximately 1.7 times to that of the first peak, indicating the existence of the typical tetrahedral liquid of water molecules.[13] However, in our present work, the second peak is extended from 4.5  Å to 6  Å as concentration increases, especially for the concentration of 2.82  M that the position of the second peak is nearly 2 times that of the first peak. It is already the typical Lennard-Jones liquid with the spherical symmetry structure. Surprisingly, the first peak of g(r)OO in each case is hardly affected. Consequently, adding NaCl or increasing the temperature can both largely affect the bulk structure of water solution, however, the first hydration shell frame in NaCl solutions is largely reserved. It seems that what the ions strongly influence is the bulk H-bonds networks rather than the water molecules in the hydration shell. This is different from the structure changes when decreasing the polarity of water molecules.[19]

3.3. Hydration processes and cluster formation in different halogen solutions

In order to investigate the effects of halogen elements at high temperatures and pressures, the hydration processes and cluster formation at 500  K under 100  atm are explored. As shown in Fig.  3(a), the g(r)OO and g(r)OH in these simulations almost coincide, in that replacing the chloridion by other different halogen elements has no further effects to the already weakly bonded water molecules. However, the hydration processes of Na+ are different in each case that the RDF peak in the NaF solution is the lowest and getting bigger when the size of anion is bigger. For the hydration of halogen elements in Fig.  3(c), both the positions of the first and second shells increase gradually with the atomic number. In short, there are no further effects to the bulk water structure by adding different halogen elements and the hydration of the cations and anions are regularly distributed.

Fig.  3. RDFs of O– O (a), O– Na (b), H– anions (c), Na– anions (d), Na– Na (e), anions-anions (f) in different solutions: NaF, NaCl, NaBr, NaI, under 100  atm at 500  K.

Figure  3(d) is the RDFs between the Na+ and the halogen elements in different solutions. The values of the first and second peaks gradually decrease when the atomic number increases. Compared to the oxygen of water molecular, the Na+ is more preferable to be interacted with the halogen, especially for F. However, when the atomic number increases, the force for hydration (Fig.  3(b), Fig.  3(c)) is comparable with the force for cluster formation (Fig.  3(d)). The halogen elements are therefore ‘ robbed’ by the water molecules. The structured cluster is also significant in the NaF solution that the values of the peaks are biggest for both g(r)NaNa and g(r)anionsanions in Figs.  3(e) and 3(f). On the contrary, the RDF peaks in other solutions are all less than 1, reflecting a lower probability for forming clusters.

For a better understanding, RDFs of Na– F, F– F, I– I in NaF and NaI solutions under 100  atm at 500  K, 550  K, and 600  K are shown in Fig.  4. The solutions are all controlled at liquid state by increasing or decreasing the temperature abruptly near by the phase changing point. Regarding Fig.  4(a), the value of g(r)NaF is extremely high at 600  K, above the phase changing point of water. The corresponding coordination number of Na– F is near 1.0, so that almost all the solvents are bonded. The big clusters in these conditions are shown in Fig.  5 that parts of the clusters are already precipitated.

Fig.  4. RDFs of Na– F (a), F– F (b) in NaF solutions, and I– I (c) in NaI solutions under 100  atm at 500  K, 550  K, and 600  K, respectively.

Fig.  5. Clusters of Na+ and F in NaF solution at 600  K under 100  atm. The water molecules are ignored for clearity.

However, the situations are different for the bigger solvents. When the g(r)FF in Fig.  4(b) significantly increases from 2.3 to 4.8 as the temperature increases to 600  K, the cluster formation is hardly seen in the NaI solution. The same results can be found in the NaBr solution that no increase of its RDF peak is found as the temperature increases. It seems that the smaller ions are more likely to cluster together at the phase changing point, suggesting the existence of the solvent precipitation in these solutions.

4. Conclusions

In summary, we have presented the MD and QM simulations of the pure water under high pressures and temperatures. The QM simulation at ambient condition is consistent with the MD simulation, convinced that the MD simulation method is suitable for a water system at high pressures and temperatures. The corresponding characters of pure water and NaCl solutions are investigated by using radial distribution functions and coordination numbers.

Both adding NaCl and increasing temperature can affect the bulk water structure. However, no further effects to the already weakly bonded water molecules are found by adding different halogen ions. The hydration of the ions are regularly distributed, suggesting the non-dominated role to the change of the bulk water structure.

At the phase changing point, the hydration force of sodium is comparable with that for cluster formation in NaI solution. Therefore, the sodium ions in NaBr or NaI solutions are more easy to be ‘ robbed’ by the water molecules. Meanwhile, the sodium ions in NaF solution are more inclined to be clustered together, suggesting the existence of the solutes precipitation in such conditions. As to the impact of salt to the structure of water at high temperature and pressure, whether there is a synergistic effect is still not concluded, and the detailed mechanisms should be more complicated. Nonetheless, the impact of salt should be limited for the already weakly bonded water molecules at high temperature and pressures. In addition, the force field and the related parameters should be further improved for these extreme conditions in the future.

Reference
1 Skinner J L 2010 Science 328 985 DOI:10.1126/science.1190093 [Cited within:1]
2 Kropmen M F and Bakker H J 2001 Science 291 2118 DOI:10.1126/science.1058190 [Cited within:1]
3 Pinilla C, Pópolo M G D, Lynden-Bell R M and Kohanoff J 2005 J. Phys. Chem. B 109 17922 DOI:10.1021/jp052999o [Cited within:1]
4 Gu B, Zhang F S, Wang Z P and Zhou H Y 2008 J. Chem. Phys. 129 184505 DOI:10.1063/1.3002485 [Cited within:1]
5 Schwegler E, Galli G and Gygi F 2000 Phys. Rev. Lett. 84 2429 DOI:10.1103/PhysRevLett.84.2429 [Cited within:1]
6 Stiopkin I V, Weeraman C, Pieniazek P A, Shalhout F Y and Skinner J L 2011 Nature 474 192 DOI:10.1038/nature10173 [Cited within:1]
7 Vaida V 2011 J. Chem. Phys. 135 020901 DOI:10.1063/1.3608919 [Cited within:1]
8 Ji M B, Odelius M and Gaffney K J 2010 Science 328 1003 DOI:10.1126/science.1187707 [Cited within:1]
9 Errington J R and Debenedetti P G 2001 Nature 409 318 DOI:10.1038/35053024 [Cited within:1]
10 McBride C, Noya E G, Aragones J L and Vega C 2012 Phys. Chem. Chem. Phys. 14 10140 DOI:10.1039/c2cp40962c [Cited within:1]
11 Choukrouna M and Grasset O 2007 J. Chem. Phys. 127 124506 DOI:10.1063/1.2768957 [Cited within:1]
12 Dunaeva A N, Antsyshkin D V and Kuskov O L 2010 Solar System Res. 44 202 DOI:10.1134/S0038094610030044 [Cited within:1]
13 Tielrooij K J, Garcia-Araez N, Bonn M and Bakker H J 2010 Science 328 1006 DOI:10.1126/science.1183512 [Cited within:2]
14 Goncharov A F, Goldman N, Fried L E, Crowhurst J C, Kuo I F, Mundy C J and Zaug J M 2005 Phys. Rev. Lett. 94 125508 DOI:10.1103/PhysRevLett.94.125508 [Cited within:1]
15 Gu B, Zhang F S, Wang Z P and Zhou H Y 2008 Phys. Rev. Lett. 100 088104 DOI:10.1103/PhysRevLett.100.088104 [Cited within:1]
16 Zhang F S and Lyneden-bell R M 2003 Mol. Phys. 101 1641 DOI:10.1080/0026897031000078685 [Cited within:1]
17 Chen B, Ivanor I, Klein M L and Parrinello M 2003 Phys. Rev. Lett. 91 215503 DOI:10.1103/PhysRevLett.91.215503 [Cited within:1]
18 Shen H, Cheng W and Zhang F S 2014 Phys. Rev. E 89 022717 DOI:10.1103/PhysRevE.89.022717 [Cited within:1]
19 Shen X, Gu B, Che S A and Zhang F S 2011 J. Chem. Phys. 135 034509 DOI:10.1063/1.3610549 [Cited within:2]
20 Klug D D, Tulk C A and Svensson E C 1999 Phys. Rev. Lett. 83 2584 DOI:10.1103/PhysRevLett.83.2584 [Cited within:1]
21 Torre R, Bartlini P and Righini R 2004 Nature 428 296 DOI:10.1038/nature02409 [Cited within:1]
22 Gaiduk V I and Crothers D S F 2006 J. Mol. Struct. 798 75 DOI:10.1016/j.molstruc.2006.03.052 [Cited within:1]
23 Bernabei M, Botti A, Bruni F, Ricci M A and Soper A K 2008 Phys. Rev. E 78 021505 DOI:10.1103/PhysRevE.78.021505 [Cited within:1]
24 Ohtaki H 2003 J. Mol. Liq. 103 3 [Cited within:1]
25 Wernet Ph, Testemale D, Hazemann J L, Argoud R, Glatzel P, Pettersson L G M, Nilsson A and Bergmann U 2005 J. Chem. Phys. 123 154503 DOI:10.1063/1.2064867 [Cited within:1]
26 Lu J, Brown J S, Liotta C L and Eckert C A 2001 Chem. Commun. 2001 665 [Cited within:1]
27 Koneshan S and Rasaiah J C 2000 J. Chem. Phys. 113 8125 DOI:10.1063/1.1314341 [Cited within:1]
28 Noya E G, Sesé L M, Ramírez R, McBride C, Conde M M and Vega C 2011 Mol. Phys. 109 149 DOI:10.1080/00268976.2010.528202 [Cited within:1]
29 Perera A 2011 Mol. Phys. 109 2433 DOI:10.1080/00268976.2011.617712 [Cited within:1]
30 Lynden-Bell R M and Debenedetti P G 2005 J. Phys. Chem. B 109 6527 [Cited within:1]
31 Mirzoev A and Lybartsev A P 2011 Phys. Chem. Chem. Phys. 13 5722 DOI:10.1039/c0cp02397c [Cited within:1]
32 Koneshan S and Rasaiah J C 2000 J. Chem. Phys. 113 8125 DOI:10.1063/1.1314341 [Cited within:1]
33 Noworyta J P, Koneshan S and Jayendran C R 2000 J. Am. Chem. Soc. 122 11194 DOI:10.1021/ja001979r [Cited within:1]
34 Jorgensen W L and Tirado-Rivers J 2005 Proc. Natl. Acad. Sci. USA 102 6665 DOI:10.1073/pnas.0408037102 [Cited within:1]
35 Smith D E and Dang L X 1994 J. Chem. Phys. 100 3757 DOI:10.1063/1.466363 [Cited within:1]
36 Liu J, Wang Z P, Zhang M, Cukier R I and Bu X Y 2013 Phys. Rev. Lett. 110 107602 DOI:10.1103/PhysRevLett.110.107602 [Cited within:1]
37 Leung K and Rempe S B 2006 Phys. Chem. Chem. Phys. 8 2153 DOI:10.1039/b515126k [Cited within:1]
38 Lkeda T, Katayama Y, Saitoh H and Aoki K 2010 J. Chem. Phys. 132 121102 DOI:10.1063/1.3374812 [Cited within:1]
39 The cp2k developers group 2000-2009 http://cp2k.berlios.de/ [Cited within:2]
40 Elstner M, Porezag D, Elsner J, Hargk M and Frauenheim Th 1998 Phys. Rev. B 58 7260 DOI:10.1103/PhysRevB.58.7260 [Cited within:2]
41 Perdew J P, Burke K and Ernzerhof M 1996 Phys. Rev. Lett. 77 3865 DOI:10.1103/PhysRevLett.77.3865 [Cited within:1]
42 Schmidt J, Vondele J V, Kuo W, Sebastiani D and Mundy C J 2009 J. Phys. Chem. B 113 11959 DOI:10.1021/jp901990u [Cited within:1]
43 Wang S, Lu H J, Tu Y S, Wang C L and Fang H P 2009 Chin. Phys. Lett. 26 68702 DOI:10.1088/0256-307X/26/6/068702 [Cited within:1]
44 Li J Y, Yang Z X, Fang H P, Zhou R H and Tang X W 2007 Chin. Phys. Lett. 24 2710 [Cited within:1]