- Email: [email protected]
Free energy of solubility of non-polar particles in water: The role of pressure
Chemical Physics 434 (2014) 15–19
Contents lists available at ScienceDirect
Chemical Physics journal homepage: www.elsevier.com/locate/chemphys
Free energy of solubility of non-polar particles in water: The role of pressure C. Gastón Ferrara a,b, J. Raul Grigera c,⇑ a
Instituto de Investigaciones Fisicoquímicas Teóricas y Aplicadas, CONICET La Plata – UNLP, La Plata, Argentina Universidad Nacional Arturo Jauretche, Florencio Varela, Argentina c CEQUINOR (CONICET La Plata – UNLP) and Departamento de Cs. Biológicas, Facultad de Ciencias Exactas, Universidad Nacional de La Plata, Argentina b
a r t i c l e
i n f o
Article history: Received 21 October 2013 In final form 30 January 2014 Available online 28 February 2014 Keywords: Free energy Solubility of non-polar compounds in water Non-polar aggregates Hydrophobic interaction High pressure
a b s t r a c t The behavior of aqueous solutions of non-polar substances at different temperatures has been studied extensively. The interpretation of the behavior of these systems under pressure is somewhat controversial when the attempt is made to connect the hydrocarbon model with the folding–unfolding of proteins under pressure. We consider the free energy of solubility of non-polar substances in water by molecular dynamics simulation to compute changes in the free energy of transference of a Lennard-Jones particle to an existent cluster of such particles in water from 1 to1500 bar and 300 ± 4 K. The critical cluster size needed to form an aggregate increase with pressure. Using the values of the free energy of transference, an analytical expression is proposed to describe its behavior. The formation of clusters requires a tetrahedral structure to allow hydrophobic interaction, then, decrease of cluster size is due to the change of water properties under pressure. Ó 2014 Elsevier B.V. All rights reserved.
1. Introduction The hydrophobic behavior of hydrocarbons in water and its dependence on temperature show similar characteristics to those of the stabilization of many structures of biological origin, at least with respect to some thermodynamic aspects. However, the attempt to link the model of hydrocarbons in water with protein folding and its stability when facing temperature and pressure changes is not yet clear [1] and is mainly based on the role of pressure. The solubility changes of methane and other hydrocarbons in water when temperature changes have been studied experimentally [2–4] and by computer simulation [5–10]. The influence of pressure [11–14], however, has deserved less attention. Studies done by means of simulations of methane in water show that the aggregates tend to stabilize at atmospheric pressure, while at higher pressures, within a range of 7–8 kbar, the predominant components are monomers, dimmers and trimmers, and the entropy change becomes insignificant at 7 kbar [11,13]. This behavior has been interpreted as a result of the changes in water structure due to pressure [15–17]. In order to extend the analysis of the solubility of non-polar substances in water at different pressures, we computed the free energy of solubility of simple Lennard-Jones (L-J) particles in water, ⇑ Corresponding author. Tel.: +54 (0) 221 425 27 67; fax: +54 (0) 221 425 9485. E-mail address: [email protected] (J.R. Grigera). http://dx.doi.org/10.1016/j.chemphys.2014.01.018 0301-0104/Ó 2014 Elsevier B.V. All rights reserved.
using a perturbation protocol for a molecular dynamics under different conditions. As L-J particles are pure non-polar, they are an excellent model for the studies of hydrophobic interaction, and hence of solubility. 2. Materials and methods 2.1. MD simulations The molecular dynamics simulation was carried out using the GROMACS 4.0.5 package [18] in which the equations of motion are solved using a leap-frog integration step. Water molecules were constrained using the SETTLE algorithm [19]. Long-range interactions were computed using the Particle Mesh Ewald method [20,18]. GNU/Linux was used for all the simulation runs and MS Window particles or GNU/Linux for all plots and graphics. The simple points charge/extended (SPC/E) [21] water model was used for all systems. Lennard-Jones particles interactions were of 0.42 nm diameter and minimum energy of 0.65017 kJ/mol; the cut-off radius was setting at 0.9 nm. All simulations were run under periodic boundary conditions and a different number seed was used to assign the initial velocities. The systems were weakly coupled to a Berendsen’s thermal and hydrostatic bath in order to work in the isothermal–isobaric ensemble [22], and equilibrated at each pressure and temperature condition for 3000 ps. In all cases, the bath temperature was
16
C.G. Ferrara, J.R. Grigera / Chemical Physics 434 (2014) 15–19
maintained at 300 ± 4 K and the bath pressure at 1, 235 and 1500 bar. The integration step used was 2 fs. Initially, all the systems consisted of a simulation cubic box with different concentrations of Lennard-Jones particles, thereby obtaining simulation boxes with different initial volumes. 2.2. Determination of cluster size in the non-polar solute After equilibration, the simulation was continued for 2000 ps more, in which data were collected (4000 points) for the analysis of the cluster sizes at each pressure and concentration. Clusters were defined using a criterion of instantaneous direct connectivity, according to which a connecting time is not required and connectivity is checked only for each simulation frame. We used the direct connectivity criterion of Stillinger [23], which is based solely upon the distance between particles. We considered that two solute particles were directly connected if the distance between them was less than 1.5 times their radius. A radius of 0.42 nm was selected for solute particles. 2.3. Determination of free energy Free energy computation was carried out by adding a particle into an already formed cluster, using the soft-core interactions [18,14]. The procedure is based on the knowledge of the cluster configuration. Having fixed the solute particles to a reference position from a known configuration, we were able to choose one of the particles that compose the aggregate and apply a slow-growth mechanism [18,24], keeping the other solute particles fixed in the reference position and allowing water to move freely. The duration of the simulations used for calculating the free energy for this mechanism was 2000 ps. The integration step used was 0.2 fs, and the order of the perturbation, delta lambda, was 0.000001. All the simulations were performed maintaining the same temperature and pressure used during stabilization. This process was performed for three different configurations for each size of aggregate and for each pressure condition. 3. Results and discussion 3.1. Kinetic behavior of the system and effect of pressure In Fig. 1 we can observe the evolution of the formation or cluster size for a system of 60 L-J particles in 2000 SPC/E water
molecules along 5000 ps of simulation time. The simulation system starts with a homogenous distribution of the particles in the box. The aggregation proceeds until reaching the minimal size; then the mean cluster size fluctuates around a value determined by the ratio between solute and solvent areas [25]. The existence of a minimum critical size in the formation of non-polar solutes aggregates in water has been observed previously for water–methane solutions at low pressure and temperature [8,26], and at low temperature and high pressure [11]. The critical size of the aggregates obtained changes with the methane concentration. The experimental results for the solubility of hydrocarbons in water [27–29], in this range of temperature and pressure, show a phase separation in the solution, consistent with the formation of aggregates for which the size depends on concentration. We did not observe a significant change in the size of large aggregates at low pressures throughout the simulations, as shown in Fig. 2. When the pressure increase became significant, we observed important changes in the cluster size; monomers, dimmers or trimmers were the predominant species. This result is consistent with the classic works of Wallqvist [11]. These results clearly show that, as pressure increases, there are a loss of the cooperative effect, and the critical size of the cluster decreases. The analysis of the cluster size shows that, regardless of the pressure and concentration condition, the species that were observed more frequently during the analysis stage are single-particle ones (over 60–70% of the total). Observing this, we studied the free energy necessary to incorporate a particle into an already formed aggregate, choosing those with the more frequent size (see Tables 1 and 2). The formation of clusters depends not only on the concentration, but also on the pressure [30]. To analyze this dependence, we measured the free energy of cluster formation at three different pressures: 1, 235 and 1500. The choice of the pressure level for the simulations was based on our previous results [30], where we defined a coexistence curve between the aggregation state and the non-aggregation state for different conditions of pressure and temperature. The 22 systems simulated in the process of stabilization are shown in Tables 1–3. Each table shows the composition of each system as well as the range of more frequent cluster sizes after stabilization. The data of Table 1 (1 bar) corresponds to a system which clearly shows two different phases, one of which is composed by
Fig. 1. Mean cluster size versus time along the formation process for a system of Lennard-Jones particle in water at 300 K and 1 bar. Each curve corresponds to the same system with different initial conditions. The system consists of 60 Lennard-Jones particles in 2000 SPC/E water molecules.
17
C.G. Ferrara, J.R. Grigera / Chemical Physics 434 (2014) 15–19
Fig. 2. Normalized mean cluster size for a system of 60 Lennard-Jones particles in 2000 water molecules (SPC/E) at 300 K vs. pressure.
single solute particles. Table 2 (235 bar) corresponds to a condition of pressure and temperature in which the system is unstable and does not show a definite phase. When a small pressure variation occurs, the system fluctuates between one or two phases. For the last case (Table 3), with pressure of 1500 bar, the system cannot form large aggregates unless we drastically increase the concentration of the system.
3.2. Pressure-dependence of the free energy of solubility We computed the free energy on 113 systems in different conditions of pressure and concentration. The calculation was repeated in each case for three different initial configurations in which the particle that was selected to be incorporated to the aggregate was different in each case (a total of 339 systems).
Table 1 Systems stabilized at 1 bar. System Pressure A1 A2 A3 A4 A5 A6 A7 A8 A9
Number of LJ
Number of SPCE
1 bar and 300 K 80 2584 75 2606 67 2665 60 2000 60 2665 40 999 40 1099 30 1000 25 900
Fraction of particles solute
Cluster size occur more frequently
0.030 0.027 0.024 0.029 0.022 0.038 0.035 0.029 0.027
75–79 68–70 60–64 49–53 1–5 36–40 26–29 18–22 15–21
Fraction of particles solute
Cluster size occur more frequently
0.038 0.047 0.038 0.024 0.032 0.029 0.05
52–55 and 57, 59 34–40 32–38 1–6 50–53 47–50 and 54 68–73
Table 2 Systems stabilized at 235 bar. System Pressure B1 B2 B3 B4 B5 B6 B7
Number of LJ
Number of SPCE
235 bar and 300 K 60 1500 40 800 40 1000 40 1600 60 1800 60 2000 80 1500
Table 3 Systems stabilized at 1500 bar. System Pressure C1 C2 C3 C4 C5 C6
Number of LJ
Number of SPCE
1500 bar and 300 K 60 1200 60 2665 60 600 40 800 40 900 80 1500
Fraction of particles solute
Cluster size occur more frequently
0.047 0.022 0.09 0.047 0.042 0.05
44–46 1–6 58–60 27–28 and 32–35 7–9 and 11 71–73 and 77–78
Fig. 3 shows the mean energy necessary to incorporate a solute into an aggregate with respect to the final cluster size. The energy values correspond to the average between the three simulated configurations for each cluster size and pressure. Analyzing the solubility value at 1 bar (first point, Fig. 3-A) and the thermal energy (kT) we can see that the ratio solubility/kT is around 3.1, a value that is closer to the values obtained by Widom et al. [5]. The results shown in Figs. 3 and 2 clearly show the influence of pressure on aggregate formation. This is seen in aggregates in the number of particles of the minimum critical size for which the change in free energy becomes negative. We define the minimum critical size as the size corresponding to the smallest aggregate in which the free energy of insertion becomes favorable. These results clearly show that the formation of large aggregates is favored and that this effect is lost when the pressure is increased significantly. The study of the hydrophobic interaction dependence on the pressure from the formation of dimmers or trimmers does not appear to be fully adequate; therefore, this description appears to be suitable only in a few cases. To describe the free energy change in terms of the final cluster size and its pressure dependence, we propose an exponential expression such as:
DG ¼ DG0 þ AeN=b
ð1Þ
where DG is the free energy necessary to incorporate a particle into an existing aggregate; DG0 is the free energy of solubility, and N the final number of particles composing the aggregate. If we assume that, for the model, the critical minimal size of the aggregate is
18
C.G. Ferrara, J.R. Grigera / Chemical Physics 434 (2014) 15–19
Fig. 3. Free energy change needed of insertion of a particle into an aggregate at the final cluster size. (A) 1 bar; (B) 235 bar; (C) 1500 bar.
obtained when DG is equal to zero, we see that the critical value depends on N and b. The model correctly shows the dependence of the minimum critical size (see Table 4) and the effect of pressure as a disaggregation effect on formed clusters. Comparing the exponential fit to the model presented by Raschke and Levitt [26] for methane and benzene simulations in water at 1 bar, the results are similar, with a very small improvement at low pressure. This small difference between the models vanishes at higher pressures (see Fig. 3C). Eq. (1) describes correctly the free energy of cluster formation at different pressures, indicating clearly the point at which the aggregates become unstable. The driving force of cluster formation for the present system (LJ particles in water) is the hydrophobic
Table 4 Minimum critical size parameters obtained with the adjustment made with the exponential model at different pressures. Pressure/bar
Minimal critical size
1 235 1500
41 43 132
interaction, which depends on the pre-existent hydrogen bond lattice (HBL); hydrophobic contacts produce a partial disassembling of the HBL, increasing the entropy. The hydrogen bond lattice is altered when pressure increases, moving from the regular open
C.G. Ferrara, J.R. Grigera / Chemical Physics 434 (2014) 15–19
tetrahedral structure of low density to a more compact hexagonal one, losing the hydrogen bonds gradually. The transition from the low density liquid (LDL) to the high density liquid (HDL) has been characterized, allowing to detect the change of one state to the other [31]. The loss of the hydrogen bonds due to pressure inhibits the hydrophobic interaction, causing the aggregates to be unstable. The present results are compatibles with previous works, related with protein denaturation. Using theoretical [32], computer simulation [33–35], and experimental approach [36], it was shown that increasing pressure produce the disturbance on the tetrahedral hydrogen bond lattice which, consequently weakened the hydrophobic interaction. The simple system studied here reproduces the effect showing clearly that at high pressure the hydrophobic interaction is not viable. 4. Conclusions The results obtained with the present method at 1 bar are consistent with previous reports [26,4] as regards the change of free energy versus the size of the aggregates, the solubility of a single non polar particle in water, and the fitting of the change of energy needed to add a particle to already formed aggregates. We therefore consider this method as reliable at other pressures. We have defined the critical cluster size as the minimum size at which the free energy became negative. The criterion used to fix the size and formation of aggregates corresponds to what is usually in the literature known as conditions of minimum contact (CM). Our results show that at low pressures (1 and 235 bar), changes in the energy required to incorporate a particle to an already formed cluster and the critical minimum sizes do not differ significantly. However, when pressure is substantially greater (in our case 1500 bar), the energy required to incorporate a solute particle to an aggregate changes dramatically with respect to the conditions of low pressure. This change is due to the disruption of water–water H-bonds avoiding the formation of tetrahedral water structure and, consequently, dramatically weakening the hydrophobic Due to the simplicity of the proposed model it may be of interest in technological problems dealing with water/non-polar systems.
19
Acknowledgements We wish to thank Dr. Osvaldo Chara for his useful comments. The work was funded by the National Scientific Council of Argentina (CONICET) and the University of La Plata. C.G. Ferrara was supported also by the Bunge y Born Foundation. JRG is member of the Research Carrere of CONICET. References [1] [2] [3] [4] [5] [6] [7] [8] [9] [10] [11] [12] [13] [14] [15] [16] [17] [18] [19] [20] [21] [22] [23] [24] [25] [26] [27] [28] [29] [30] [31] [32] [33] [34] [35] [36]
W. Kauzmann, Nature 325 (1987) 763. K. Lekvan, P.R. Bishno, Fluid Phase Equilib. 131 (1997) 297. K.S. Pedersen, J. Milter, C.P. Rasmusen, Fluid Phase Equilib. 189 (2001) 85. P. Servio, P. Englezos, J. Chem. Eng. Data 47 (2002) 87. B. Widom, P. Bhimalapuram, K. Koga, Phys. Chem. Chem. Phys. 5 (2003) 3085. R. Zangi, B.J. Berne, J. Phys. Chem. B 112 (2008) 8634. S.L. Lee, P.G. Debenedetti, J.R. Errington, J. Chem. Phys. 122 (2005) 204511. C. Oostembrink, W.F. van Gunsteren, Phys. Chem. Chem. Phys. 7 (2004) 53. S.W. Rick, J. Phys. Chem. B 107 (2003) 9853. N.T. Southall, K.A. Dill, A.D.J. Haymet, J. Phys. Chem. 106 (2002) 521. A. Wallqvist, J. Phys. Chem. 95 (1991) 8921. K. Koga, Chem. Phys. 121 (2004) 7304. S.W. Rick, J. Phys. Chem. 104 (2000) 6488. T. Ghosh, A.E. Garcia, S. Garde, J. Am. Chem. Soc. 123 (2001) 10997. A.Y. Wu, E. Whalley, G. Dolling, Mol. Phys. 47 (1982) 603. A.V. Okhulkov, Y.N. Demianets, Y.E. Gorbaty, J. Chem. Phys. 100 (1994) 1578. M. Bellisent-Funel, L.J. Bosio, J. Chem. Phys. 102 (1995) 3727. D. Van Der Spoel, E. Lindahl, B. Hess, G. Groenhof, A.E. Mark, H.J.C. Berendsen, J. Comput. Chem. 26 (2005) 1701. S. Miyamoto, P.A. Kollman, J. Comp. Chem. 13 (1992) 952. T. Darden, D.M. York, L.G. Pedersen, J. Chem. Phys. 98 (10089–10092) (1993) 21. H.J.C. Berendsen, J.R. Grigera, T.J. Straatsma, J. Phys. Chem. 91 (1987) 6269. H.J.C. Berendsen, J.P.M. Postma, A. DiNola, J.R. Haak, J. Phys. Chem. 81 (1984) 3684. F.H. Stillinger, J. Chem. Phys. 38 (1963) 1486. P. Kollman, Chem. Rev. 93 (1993) 2395. G. Ferrara, A.N. McCarthy, J.R. Grigera, J. Chem. Phys. 127 (2007) 104502. T.M. Raschke, J. Tsai, M. Levitt, Proc. Natl. Acad. Sci. U.S.A 98 (2001) 5965. A. Dhima, J.C. Hemptinne, G. Moracchini, Fluid Phase Equilib. 145 (1998) 129. A. Dhima, J.C. Hemptinne, J. Jose, Ind. Eng. Chem. Res. 38 (1999) 3144. T.R. Rettlch, P. Hand, R. Battino, E. Wllhelm, J. Chem. Phys. 85 (1981) 3230. C.G. Ferrara, O. Chara, J.R. Grigera, J. Chem. Phys. 137 (2012) 135104. O. Chara, A.N. McCarthy, J.R. Grigera, Phys. Lett. A 375 (2011) 572. G. Hummer, S. Garde, A.E. Garcia, M.E. Paulaitis, L.R. Pratt, Proc. Natl. Acad. Sci. U.S.A 95 (1998) 1552. T. Ghosh, A.E. Garcia, S. Garde, J. Chem. Phys. 116 (2002) 2480. S. Sarupria, T. Ghosh, A.E. Garcia, S. Garde, Proteins 78 (2009) 1641. J.R. Grigera, A.N. McCarthy, Biophys. J. 98 (1626) (2010) 1631. J.B. Rouget, T. Aksel, J. Roche, J.L. Saldana, A.E. García, D. Barrick, C. Royer, J. Am. Chem. Soc. 133 (2011) 6020.
Contents lists available at ScienceDirect
Chemical Physics journal homepage: www.elsevier.com/locate/chemphys
Free energy of solubility of non-polar particles in water: The role of pressure C. Gastón Ferrara a,b, J. Raul Grigera c,⇑ a
Instituto de Investigaciones Fisicoquímicas Teóricas y Aplicadas, CONICET La Plata – UNLP, La Plata, Argentina Universidad Nacional Arturo Jauretche, Florencio Varela, Argentina c CEQUINOR (CONICET La Plata – UNLP) and Departamento de Cs. Biológicas, Facultad de Ciencias Exactas, Universidad Nacional de La Plata, Argentina b
a r t i c l e
i n f o
Article history: Received 21 October 2013 In final form 30 January 2014 Available online 28 February 2014 Keywords: Free energy Solubility of non-polar compounds in water Non-polar aggregates Hydrophobic interaction High pressure
a b s t r a c t The behavior of aqueous solutions of non-polar substances at different temperatures has been studied extensively. The interpretation of the behavior of these systems under pressure is somewhat controversial when the attempt is made to connect the hydrocarbon model with the folding–unfolding of proteins under pressure. We consider the free energy of solubility of non-polar substances in water by molecular dynamics simulation to compute changes in the free energy of transference of a Lennard-Jones particle to an existent cluster of such particles in water from 1 to1500 bar and 300 ± 4 K. The critical cluster size needed to form an aggregate increase with pressure. Using the values of the free energy of transference, an analytical expression is proposed to describe its behavior. The formation of clusters requires a tetrahedral structure to allow hydrophobic interaction, then, decrease of cluster size is due to the change of water properties under pressure. Ó 2014 Elsevier B.V. All rights reserved.
1. Introduction The hydrophobic behavior of hydrocarbons in water and its dependence on temperature show similar characteristics to those of the stabilization of many structures of biological origin, at least with respect to some thermodynamic aspects. However, the attempt to link the model of hydrocarbons in water with protein folding and its stability when facing temperature and pressure changes is not yet clear [1] and is mainly based on the role of pressure. The solubility changes of methane and other hydrocarbons in water when temperature changes have been studied experimentally [2–4] and by computer simulation [5–10]. The influence of pressure [11–14], however, has deserved less attention. Studies done by means of simulations of methane in water show that the aggregates tend to stabilize at atmospheric pressure, while at higher pressures, within a range of 7–8 kbar, the predominant components are monomers, dimmers and trimmers, and the entropy change becomes insignificant at 7 kbar [11,13]. This behavior has been interpreted as a result of the changes in water structure due to pressure [15–17]. In order to extend the analysis of the solubility of non-polar substances in water at different pressures, we computed the free energy of solubility of simple Lennard-Jones (L-J) particles in water, ⇑ Corresponding author. Tel.: +54 (0) 221 425 27 67; fax: +54 (0) 221 425 9485. E-mail address: [email protected] (J.R. Grigera). http://dx.doi.org/10.1016/j.chemphys.2014.01.018 0301-0104/Ó 2014 Elsevier B.V. All rights reserved.
using a perturbation protocol for a molecular dynamics under different conditions. As L-J particles are pure non-polar, they are an excellent model for the studies of hydrophobic interaction, and hence of solubility. 2. Materials and methods 2.1. MD simulations The molecular dynamics simulation was carried out using the GROMACS 4.0.5 package [18] in which the equations of motion are solved using a leap-frog integration step. Water molecules were constrained using the SETTLE algorithm [19]. Long-range interactions were computed using the Particle Mesh Ewald method [20,18]. GNU/Linux was used for all the simulation runs and MS Window particles or GNU/Linux for all plots and graphics. The simple points charge/extended (SPC/E) [21] water model was used for all systems. Lennard-Jones particles interactions were of 0.42 nm diameter and minimum energy of 0.65017 kJ/mol; the cut-off radius was setting at 0.9 nm. All simulations were run under periodic boundary conditions and a different number seed was used to assign the initial velocities. The systems were weakly coupled to a Berendsen’s thermal and hydrostatic bath in order to work in the isothermal–isobaric ensemble [22], and equilibrated at each pressure and temperature condition for 3000 ps. In all cases, the bath temperature was
16
C.G. Ferrara, J.R. Grigera / Chemical Physics 434 (2014) 15–19
maintained at 300 ± 4 K and the bath pressure at 1, 235 and 1500 bar. The integration step used was 2 fs. Initially, all the systems consisted of a simulation cubic box with different concentrations of Lennard-Jones particles, thereby obtaining simulation boxes with different initial volumes. 2.2. Determination of cluster size in the non-polar solute After equilibration, the simulation was continued for 2000 ps more, in which data were collected (4000 points) for the analysis of the cluster sizes at each pressure and concentration. Clusters were defined using a criterion of instantaneous direct connectivity, according to which a connecting time is not required and connectivity is checked only for each simulation frame. We used the direct connectivity criterion of Stillinger [23], which is based solely upon the distance between particles. We considered that two solute particles were directly connected if the distance between them was less than 1.5 times their radius. A radius of 0.42 nm was selected for solute particles. 2.3. Determination of free energy Free energy computation was carried out by adding a particle into an already formed cluster, using the soft-core interactions [18,14]. The procedure is based on the knowledge of the cluster configuration. Having fixed the solute particles to a reference position from a known configuration, we were able to choose one of the particles that compose the aggregate and apply a slow-growth mechanism [18,24], keeping the other solute particles fixed in the reference position and allowing water to move freely. The duration of the simulations used for calculating the free energy for this mechanism was 2000 ps. The integration step used was 0.2 fs, and the order of the perturbation, delta lambda, was 0.000001. All the simulations were performed maintaining the same temperature and pressure used during stabilization. This process was performed for three different configurations for each size of aggregate and for each pressure condition. 3. Results and discussion 3.1. Kinetic behavior of the system and effect of pressure In Fig. 1 we can observe the evolution of the formation or cluster size for a system of 60 L-J particles in 2000 SPC/E water
molecules along 5000 ps of simulation time. The simulation system starts with a homogenous distribution of the particles in the box. The aggregation proceeds until reaching the minimal size; then the mean cluster size fluctuates around a value determined by the ratio between solute and solvent areas [25]. The existence of a minimum critical size in the formation of non-polar solutes aggregates in water has been observed previously for water–methane solutions at low pressure and temperature [8,26], and at low temperature and high pressure [11]. The critical size of the aggregates obtained changes with the methane concentration. The experimental results for the solubility of hydrocarbons in water [27–29], in this range of temperature and pressure, show a phase separation in the solution, consistent with the formation of aggregates for which the size depends on concentration. We did not observe a significant change in the size of large aggregates at low pressures throughout the simulations, as shown in Fig. 2. When the pressure increase became significant, we observed important changes in the cluster size; monomers, dimmers or trimmers were the predominant species. This result is consistent with the classic works of Wallqvist [11]. These results clearly show that, as pressure increases, there are a loss of the cooperative effect, and the critical size of the cluster decreases. The analysis of the cluster size shows that, regardless of the pressure and concentration condition, the species that were observed more frequently during the analysis stage are single-particle ones (over 60–70% of the total). Observing this, we studied the free energy necessary to incorporate a particle into an already formed aggregate, choosing those with the more frequent size (see Tables 1 and 2). The formation of clusters depends not only on the concentration, but also on the pressure [30]. To analyze this dependence, we measured the free energy of cluster formation at three different pressures: 1, 235 and 1500. The choice of the pressure level for the simulations was based on our previous results [30], where we defined a coexistence curve between the aggregation state and the non-aggregation state for different conditions of pressure and temperature. The 22 systems simulated in the process of stabilization are shown in Tables 1–3. Each table shows the composition of each system as well as the range of more frequent cluster sizes after stabilization. The data of Table 1 (1 bar) corresponds to a system which clearly shows two different phases, one of which is composed by
Fig. 1. Mean cluster size versus time along the formation process for a system of Lennard-Jones particle in water at 300 K and 1 bar. Each curve corresponds to the same system with different initial conditions. The system consists of 60 Lennard-Jones particles in 2000 SPC/E water molecules.
17
C.G. Ferrara, J.R. Grigera / Chemical Physics 434 (2014) 15–19
Fig. 2. Normalized mean cluster size for a system of 60 Lennard-Jones particles in 2000 water molecules (SPC/E) at 300 K vs. pressure.
single solute particles. Table 2 (235 bar) corresponds to a condition of pressure and temperature in which the system is unstable and does not show a definite phase. When a small pressure variation occurs, the system fluctuates between one or two phases. For the last case (Table 3), with pressure of 1500 bar, the system cannot form large aggregates unless we drastically increase the concentration of the system.
3.2. Pressure-dependence of the free energy of solubility We computed the free energy on 113 systems in different conditions of pressure and concentration. The calculation was repeated in each case for three different initial configurations in which the particle that was selected to be incorporated to the aggregate was different in each case (a total of 339 systems).
Table 1 Systems stabilized at 1 bar. System Pressure A1 A2 A3 A4 A5 A6 A7 A8 A9
Number of LJ
Number of SPCE
1 bar and 300 K 80 2584 75 2606 67 2665 60 2000 60 2665 40 999 40 1099 30 1000 25 900
Fraction of particles solute
Cluster size occur more frequently
0.030 0.027 0.024 0.029 0.022 0.038 0.035 0.029 0.027
75–79 68–70 60–64 49–53 1–5 36–40 26–29 18–22 15–21
Fraction of particles solute
Cluster size occur more frequently
0.038 0.047 0.038 0.024 0.032 0.029 0.05
52–55 and 57, 59 34–40 32–38 1–6 50–53 47–50 and 54 68–73
Table 2 Systems stabilized at 235 bar. System Pressure B1 B2 B3 B4 B5 B6 B7
Number of LJ
Number of SPCE
235 bar and 300 K 60 1500 40 800 40 1000 40 1600 60 1800 60 2000 80 1500
Table 3 Systems stabilized at 1500 bar. System Pressure C1 C2 C3 C4 C5 C6
Number of LJ
Number of SPCE
1500 bar and 300 K 60 1200 60 2665 60 600 40 800 40 900 80 1500
Fraction of particles solute
Cluster size occur more frequently
0.047 0.022 0.09 0.047 0.042 0.05
44–46 1–6 58–60 27–28 and 32–35 7–9 and 11 71–73 and 77–78
Fig. 3 shows the mean energy necessary to incorporate a solute into an aggregate with respect to the final cluster size. The energy values correspond to the average between the three simulated configurations for each cluster size and pressure. Analyzing the solubility value at 1 bar (first point, Fig. 3-A) and the thermal energy (kT) we can see that the ratio solubility/kT is around 3.1, a value that is closer to the values obtained by Widom et al. [5]. The results shown in Figs. 3 and 2 clearly show the influence of pressure on aggregate formation. This is seen in aggregates in the number of particles of the minimum critical size for which the change in free energy becomes negative. We define the minimum critical size as the size corresponding to the smallest aggregate in which the free energy of insertion becomes favorable. These results clearly show that the formation of large aggregates is favored and that this effect is lost when the pressure is increased significantly. The study of the hydrophobic interaction dependence on the pressure from the formation of dimmers or trimmers does not appear to be fully adequate; therefore, this description appears to be suitable only in a few cases. To describe the free energy change in terms of the final cluster size and its pressure dependence, we propose an exponential expression such as:
DG ¼ DG0 þ AeN=b
ð1Þ
where DG is the free energy necessary to incorporate a particle into an existing aggregate; DG0 is the free energy of solubility, and N the final number of particles composing the aggregate. If we assume that, for the model, the critical minimal size of the aggregate is
18
C.G. Ferrara, J.R. Grigera / Chemical Physics 434 (2014) 15–19
Fig. 3. Free energy change needed of insertion of a particle into an aggregate at the final cluster size. (A) 1 bar; (B) 235 bar; (C) 1500 bar.
obtained when DG is equal to zero, we see that the critical value depends on N and b. The model correctly shows the dependence of the minimum critical size (see Table 4) and the effect of pressure as a disaggregation effect on formed clusters. Comparing the exponential fit to the model presented by Raschke and Levitt [26] for methane and benzene simulations in water at 1 bar, the results are similar, with a very small improvement at low pressure. This small difference between the models vanishes at higher pressures (see Fig. 3C). Eq. (1) describes correctly the free energy of cluster formation at different pressures, indicating clearly the point at which the aggregates become unstable. The driving force of cluster formation for the present system (LJ particles in water) is the hydrophobic
Table 4 Minimum critical size parameters obtained with the adjustment made with the exponential model at different pressures. Pressure/bar
Minimal critical size
1 235 1500
41 43 132
interaction, which depends on the pre-existent hydrogen bond lattice (HBL); hydrophobic contacts produce a partial disassembling of the HBL, increasing the entropy. The hydrogen bond lattice is altered when pressure increases, moving from the regular open
C.G. Ferrara, J.R. Grigera / Chemical Physics 434 (2014) 15–19
tetrahedral structure of low density to a more compact hexagonal one, losing the hydrogen bonds gradually. The transition from the low density liquid (LDL) to the high density liquid (HDL) has been characterized, allowing to detect the change of one state to the other [31]. The loss of the hydrogen bonds due to pressure inhibits the hydrophobic interaction, causing the aggregates to be unstable. The present results are compatibles with previous works, related with protein denaturation. Using theoretical [32], computer simulation [33–35], and experimental approach [36], it was shown that increasing pressure produce the disturbance on the tetrahedral hydrogen bond lattice which, consequently weakened the hydrophobic interaction. The simple system studied here reproduces the effect showing clearly that at high pressure the hydrophobic interaction is not viable. 4. Conclusions The results obtained with the present method at 1 bar are consistent with previous reports [26,4] as regards the change of free energy versus the size of the aggregates, the solubility of a single non polar particle in water, and the fitting of the change of energy needed to add a particle to already formed aggregates. We therefore consider this method as reliable at other pressures. We have defined the critical cluster size as the minimum size at which the free energy became negative. The criterion used to fix the size and formation of aggregates corresponds to what is usually in the literature known as conditions of minimum contact (CM). Our results show that at low pressures (1 and 235 bar), changes in the energy required to incorporate a particle to an already formed cluster and the critical minimum sizes do not differ significantly. However, when pressure is substantially greater (in our case 1500 bar), the energy required to incorporate a solute particle to an aggregate changes dramatically with respect to the conditions of low pressure. This change is due to the disruption of water–water H-bonds avoiding the formation of tetrahedral water structure and, consequently, dramatically weakening the hydrophobic Due to the simplicity of the proposed model it may be of interest in technological problems dealing with water/non-polar systems.
19
Acknowledgements We wish to thank Dr. Osvaldo Chara for his useful comments. The work was funded by the National Scientific Council of Argentina (CONICET) and the University of La Plata. C.G. Ferrara was supported also by the Bunge y Born Foundation. JRG is member of the Research Carrere of CONICET. References [1] [2] [3] [4] [5] [6] [7] [8] [9] [10] [11] [12] [13] [14] [15] [16] [17] [18] [19] [20] [21] [22] [23] [24] [25] [26] [27] [28] [29] [30] [31] [32] [33] [34] [35] [36]
W. Kauzmann, Nature 325 (1987) 763. K. Lekvan, P.R. Bishno, Fluid Phase Equilib. 131 (1997) 297. K.S. Pedersen, J. Milter, C.P. Rasmusen, Fluid Phase Equilib. 189 (2001) 85. P. Servio, P. Englezos, J. Chem. Eng. Data 47 (2002) 87. B. Widom, P. Bhimalapuram, K. Koga, Phys. Chem. Chem. Phys. 5 (2003) 3085. R. Zangi, B.J. Berne, J. Phys. Chem. B 112 (2008) 8634. S.L. Lee, P.G. Debenedetti, J.R. Errington, J. Chem. Phys. 122 (2005) 204511. C. Oostembrink, W.F. van Gunsteren, Phys. Chem. Chem. Phys. 7 (2004) 53. S.W. Rick, J. Phys. Chem. B 107 (2003) 9853. N.T. Southall, K.A. Dill, A.D.J. Haymet, J. Phys. Chem. 106 (2002) 521. A. Wallqvist, J. Phys. Chem. 95 (1991) 8921. K. Koga, Chem. Phys. 121 (2004) 7304. S.W. Rick, J. Phys. Chem. 104 (2000) 6488. T. Ghosh, A.E. Garcia, S. Garde, J. Am. Chem. Soc. 123 (2001) 10997. A.Y. Wu, E. Whalley, G. Dolling, Mol. Phys. 47 (1982) 603. A.V. Okhulkov, Y.N. Demianets, Y.E. Gorbaty, J. Chem. Phys. 100 (1994) 1578. M. Bellisent-Funel, L.J. Bosio, J. Chem. Phys. 102 (1995) 3727. D. Van Der Spoel, E. Lindahl, B. Hess, G. Groenhof, A.E. Mark, H.J.C. Berendsen, J. Comput. Chem. 26 (2005) 1701. S. Miyamoto, P.A. Kollman, J. Comp. Chem. 13 (1992) 952. T. Darden, D.M. York, L.G. Pedersen, J. Chem. Phys. 98 (10089–10092) (1993) 21. H.J.C. Berendsen, J.R. Grigera, T.J. Straatsma, J. Phys. Chem. 91 (1987) 6269. H.J.C. Berendsen, J.P.M. Postma, A. DiNola, J.R. Haak, J. Phys. Chem. 81 (1984) 3684. F.H. Stillinger, J. Chem. Phys. 38 (1963) 1486. P. Kollman, Chem. Rev. 93 (1993) 2395. G. Ferrara, A.N. McCarthy, J.R. Grigera, J. Chem. Phys. 127 (2007) 104502. T.M. Raschke, J. Tsai, M. Levitt, Proc. Natl. Acad. Sci. U.S.A 98 (2001) 5965. A. Dhima, J.C. Hemptinne, G. Moracchini, Fluid Phase Equilib. 145 (1998) 129. A. Dhima, J.C. Hemptinne, J. Jose, Ind. Eng. Chem. Res. 38 (1999) 3144. T.R. Rettlch, P. Hand, R. Battino, E. Wllhelm, J. Chem. Phys. 85 (1981) 3230. C.G. Ferrara, O. Chara, J.R. Grigera, J. Chem. Phys. 137 (2012) 135104. O. Chara, A.N. McCarthy, J.R. Grigera, Phys. Lett. A 375 (2011) 572. G. Hummer, S. Garde, A.E. Garcia, M.E. Paulaitis, L.R. Pratt, Proc. Natl. Acad. Sci. U.S.A 95 (1998) 1552. T. Ghosh, A.E. Garcia, S. Garde, J. Chem. Phys. 116 (2002) 2480. S. Sarupria, T. Ghosh, A.E. Garcia, S. Garde, Proteins 78 (2009) 1641. J.R. Grigera, A.N. McCarthy, Biophys. J. 98 (1626) (2010) 1631. J.B. Rouget, T. Aksel, J. Roche, J.L. Saldana, A.E. García, D. Barrick, C. Royer, J. Am. Chem. Soc. 133 (2011) 6020.