15586
J. Phys. Chem. C 2007, 111, 15586-15595
On the Microheterogeneity in Neat and Aqueous Amides: A Molecular Dynamics Study† Larisa Zoranic´ ,‡,| Redha Mazighi,‡ Franjo Sokolic´ ,§,| and Aure´ lien Perera*,‡ Laboratoire de Physique The´ orique de la Matie` re Condense´ e (UMR CNRS 7600), UniVersite´ Pierre et Marie Curie, 4 Place Jussieu, F75252, Paris cedex 05, France, Laboratoire de Spectrochimie Infrarouge et Raman (UMR CNRS 8516), Centre d’Etudes et de Recherches Lasers et Applications, UniVersite´ des Sciences et Technologies de Lille, F59655 VilleneuVe d’Ascq Cedex, France, and Department of Physics, Faculty of Natural Sciences, Mathematics and Kinesiology, UniVersity of Split, Nikole Tesle 12, 21000, Split, Croatia ReceiVed: May 14, 2007; In Final Form: August 17, 2007
A molecular dynamics study of three amides, namely, formamide, N-methyl formamide, and dimethyl formamide, and their mixtures with water is conducted to sort out the nature of the local microstructure of both the pure amides and the aqueous mixtures. For pure amides, both the cluster distribution and pair structure factors are investigated. The former reveals no specific cluster formation for all three amides, while the latter shows that formamide has a weak prepeak at about k ≈ 1.32 Å-1 reminiscent of that observed for water, thus enforcing the puzzle of the nature of the hydrogen-bonding self-clustering in this type of liquids. The analysis of aqueous mixtures shows that all three mixtures are strongly microheterogeneous at the nanometer/nanosecond scales. This is confirmed in two ways: the direct visual observation and through the analysis of the pair distribution functions. The observation of microheterogeneity through usually measured quantities, namely the Kirkwood-Buff integrals, is discussed in light of the present analysis, and these latter quantities are shown to be weak indicators of microheterogeneity mainly because all three aqueous amides have nearly ideal such integrals in contrast to displaying strong microheterogeneity. As a direct consequence of this fact two features are reported. First, the Kirkwood-Buff integrals of the amide-amide correlations are obtained in relatively good agreement with the experimental ones, and second, those involving correlations with water are considerably affected by the microheterogeneity beyond the 15 Å range and need prohibitively larger statistics despite the large system sizes used in this study.
1. Introduction The development of a molecular theory of complex liquids has not reached yet the status of that of simple liquids, despite excellent textbooks have been made available in recent years.1,2 One of the main reasons is that the local symmetry breaking, induced by the nonspherical symmetry of the pair interactions in molecular liquids, has far-reaching consequences that are not made apparent by the statistical theory of such liquids. This remark does not concern obvious cases, such as liquid crystalline molecules, for which the relationship between microscopical and macroscopical symmetry breaking has been decrypted outside the realms of the theory of liquids within the framework of symmetry-breaking phase transitions. Rather, it concerns molecular liquids, principally associated liquids such as water, alcohols, and formamide for example, for which the form of the anistropic pair interaction does not induce a global symmetry-breaking transition within the liquid phase; hence, these liquids and their mixtures must be treated by the statistical theory of the liquid state. The major problem is the following: associated liquids present a local microstructure that is a direct consequence of the strong directionality of the pair interactions, namely the hydrogen bonding. This microstructure is to be distinguished from the density or concentration fluctuations: the former is a †
Part of the “Keith E. Gubbins Festschrift”. * To whom correspondence should be addressed. ‡ Universite ´ Pierre et Marie Curie. § Universite ´ des Sciences et Technologies de Lille. | University of Split.
ground state property due to the direct interactions between molecules, while the latter is due to thermal excitations. The statistical theory of liquids allows the latter to be related to response functions such as the compressibility and more generally to the small-k behavior of the structure factors.2 However, we lack a theory for the manifestation of the former. The main reason is the inherent disorder of liquid matter at macroscopic level. Yet, the manifestations of H-bond are apparent, as for example in water, such as the density anomaly, or the hypothetical existence of (possibly metastable) a liquidliquid-phase coexistence at low temperatures.3,4 Because liquids are disordered, we must spot manifestations of local order in the same response functions, namely through the pair distribution function. The situation becomes more complicated with mixtures of associated and/or simple liquids, because the microstructure turns into microheterogeneity (MH), which is the structure that emerges from the mixture of local ordering tendency (microstructure of neat liquids) and local disorder common to all liquids. Perhaps the most convincing evidence of this situation has been brought by the recent work in the group of A. K. Soper by showing that even methanol-water mixtures showed considerable microheterogeneity.5,6 In recent works,7,8 we have analyzed the microstructure of neat alcohols, namely methanol and tert-butanol, and shown that these H-bonding liquids had specific clustering tendencies, measurable through both an effective one-body function and the pair-level response function. The latter showed a prepeak while the former had a distinct peak at around cluster size 5 in the probability distribution of cluster sizes. Surprisingly enough,
10.1021/jp0736894 CCC: $37.00 © 2007 American Chemical Society Published on Web 10/10/2007
MD Study of Microheterogeneity in Neat and Aqueous Amides a similar analysis conducted for water did not reveal specific clustering in a convincing way;7,8 there was no peak in cluster distribution probabilities and only a weak peak in the structure factor, the latter which is equally in agreement with numerous experimental findings from various scattering experiments. It was speculated that the absence (or weak evidence) of clustering does not mean that this feature is absent, but only a different form of clustering, perhaps more “diffuse”, such as in a polymeric form, or rather a “polymorphic” clustering.9 In the present work, we would like to apply the same analysis to amides and mainly to explore the MH in water-amides mixtures. Amides are often used as template to model peptide bonds and their interaction with water is of capital importance to understand the solvation of peptides. Simple amides are generically constituted as R1OCNR2R3 with Ri being radicals and most commonly methyl groups or hydrogen atoms. The specificity of such amides is to have the same skeletal backbone OCN and geometry (that of a “bonelike” shape), so that one can focus only on the influence of different type of radicals, hydrophobic methyl groups, or hydrophilic ones. We investigate formamide (FA) (HOCNH2), N-methyl formamide (NMF) (HOCNHCH3), and dimethyl formamide (DMF) (HOCN(CH3)2) in which all three have the same R1 ) H radical and differ through the R2 and R3 radicals. All three amides are fully miscible with water under ambient conditions. Because the presence of increasing methyl groups induces higher solvophobicity, one expects all three amides to microsegregate differently in water. Moreover, FA and NMF are protic solvents, while DMF is an aprotic polar solvent, meaning FA and NMF have self-H-bonding ability, twice as much for FA than for NMF. In addition, FA has some similarity with water as it has some builtin hydrophobicity mechanism that allows this liquid to form microemulsions with appropriate surfactants, although to a lesser extent than water.10 We expect the difference in both the geometry and the chemical properties of the molecules to show up in the structural and thermodynamical properties of the neat liquid and the aqueous mixtures. Amides have been studied through computer simulations by many authors. Models for amides come in different flavors, such as the OPLS,11 modified OPLS,12 GROMOS,13 the models of Cordeiro et al.,14,15 the Kirkwood-Buff integrals (KBI) compliant model of Kang and Smith,16 and other refinements such as polarizability17 and cis/trans group asymmetry.18,19 The focus of the recent studies were on the specificity of H-bonding in neat and aqueous amides, whereas our present focus will be rather local structures induced by the pair interaction specificities such as H-bonds. One of the issues in building force field models concerns the case of aqueous mixtures. Recent computer simulations results by many authors20,16,21,22,27 have revealed that these are extremely sensitive to force field parameters, notably that of the solute, and many choices lead to more extensive selfclustering and segregation than experimentally observed, leading often to demixing-like behavior. Recent progress has put forward the necessity of fitting the solute force field parameters on the experimental KBI.16 This latter method necessitates that the radial distribution functions are calculable through computer simulations. It turns out that this calculation is directly related to the problems encountered in equilibrating and sampling these mixtures. The three different amides we have selected should then allow us to directly test the influence of hydrophilic versus hydrophobic groups with keeping a similar overall geometry of the molecule, let aside, of course, the differences in chemical properties induced by the different radicals.
J. Phys. Chem. C, Vol. 111, No. 43, 2007 15587 The remainder of this paper is then as follows. In the next section, we review the difficulties posed by force fields and simulations in describing aqueous mixtures. In the results section, we describe the clustering observed in our simulations both for pure and aqueous amides. Finally, the last section contains an overall discussion and our conclusions. 2. Theoretical and Computer Simulation Details 2.1. KBI and Microheterogeneity. KBI are the integrals of the radial distribution functions (RDF) gij(r) for a multicomponent mixture (indices i, j being related to a pair of components)which are themselves the angle averages of the total pair
Gij ) 4π
∫0∞ dr r2(gij(r) - 1)
(1)
distribution function, and the averages performed over all possible rotation of the pair of molecules. Kirkwood and Buff have shown24 that experimentally available quantities such as the compressibility, the partial molar volumes, and chemical potential density derivatives could be obtained in terms of such integrals, in the grand canonical level of description. Ben-Naim has proposed25 that, conversely, some information on the structure, albeit in an integrated form, can be obtained by reversing the procedure and extracting the Gij from the experimental data. It turns out that this procedure is not always easy, mostly because large errors are introduced by the measurement of the chemical potentials themselves or the quantities that give access to them and their density derivatives.26-28 A supplementary problem arises in computer simulations, because most of the contribution to the integral arises from the medium to long-range part of gij(r), namely that which corresponds precisely to the nanoscale where the MH develops. This is clearly related to the necessity to sample the RDFs not only in the vicinity of first neighbors but also in the range of the MH, hence subjecting it to larger statistical errors, influence of periodic boundary conditions, and other such constraints that arise naturally in finite size simulations. Hence, the need of large systems. Even then, it is not at all obvious that this procedure is entirely proper because of inherent statistical errors. A final difficulty comes from the fact that the analysis of Kirkwood and Buff pertains to concentration fluctuations, as attested by the fact that the Gij are quite simply the k ) 0 limit of the partial structure factors29
Sij(k) )
r ij(r) - 1) ∫ dbr exp(ikB‚b)(g
(2)
which can be obtained directly from various scattering experiments, such as neutron or X-ray scatterings. This fact implies that the KBIs sample in fact thermal effects of the disorder29 and not directly the MH itself. Indeed, one expects that this MH should show up at some finite nonzero k-vector in the Sij(k). It is well confirmed that in the case of microemulsions,30 a prepeak appears that attests the presence of well-formed structures. These microemulsions are formed namely with an aqueous mixture of large alcohols with a number of methyl groups above 10.31 There is a major difference between these structures and the aqueous mixtures we consider: the latter are disordered, while the former are ordered structures. Indeed, microemulsions have a nontrivial one-body density function as order parameter, whereas in our case, this is simply the number density. This fact is directly related to the puzzle of the MH: how can such systems be macroscopically disordered while being microsegregated, hence ordered at a nanoscale level? This is still an open problem of the aqueous mixtures.
15588 J. Phys. Chem. C, Vol. 111, No. 43, 2007
Zoranic´ et al.
Finally, we note that the site-site radial distributions functions (SSRDF) are equally interesting because they allow for an indirect description of the angular dependence of the total pair correlation function. We use the latter quantities in the result section. 2.2. Clustering, Microstructure, and Microheterogeneity. To palliate for the lack of true order parameter in our systems, namely for the one component systems, we have investigated an effective one-body function, namely the cluster distribution. This is defined here analogously to that in our recent works.7,8 For each configuration, we simply count how many particles appear connected in a cluster of a given size n, which is a nice way to account for local heterogeneity that will not average to zero when computed over many realizations. The probability of finding a cluster of size n is defined as
p(n) )
∑ks(n,k) ∑m∑ks(m,k)
(3)
where s(n,k) represents the number of clusters of size n in the configuration k. This definition involves an arbitrariness through the definition as to when we consider that two particles are mutually bonded within a cluster; in the case of the Stillinger cluster32,33 this is simply a distance criteria lC. In practice, this distance is chosen to correspond to the pair contact probability as described by the corresponding radial distribution function. Many such distances can be defined, as for example the one corresponding to the first peak, or anywhere around until the first minima. It has been empirically found34 that the distance corresponding to the first minima gave the largest temporal stability for clusters. For neat amides, we have checked several types of contact pair: between the center of the masses like atom pairs as well as contact pairs such as N-H·‚‚O and C-H‚‚‚O, which are directly related to H-bonding. In doing so, we have extended the usual algorithm for cluster identification to unlike sites. 2.3. Computer Simulation Details. In this work, we have systematically investigated mixtures with N ) 2048 number of particles. We have conducted constant N,P,T ensemble simulations with a Berendsen thermostat and barostat (with relaxation times of 0.1 and 0.5 ps, respectively) to ensure ambient conditions with T ) 300 K and a pressure of 1 atm. The time step was taken as 2 fs. We have taken the simple point chargeextended (SPC/E) model for water35 and the Cordeiro et al. models for amides.14 We have equally tested the OPLS models,11 which gave very similar overall results, but found that these models were more suited for transferable intermolecular potential (TIP)-like water models,36 as was equally confirmed by recent simulations of aqueous amides.12 Our focus here is not to test the model for accuracy on various dynamical quantities but rather for overall description of MH. The force field parameters are given in Table 1. For neat amides, equilibration steps were performed over few 104 steps and production runs over 20 000 steps. For mixtures, longer equilibrations were needed, and we have performed several batches of 25 000 steps (50 ps) to trace the evolution of the range of the RDF that lies in the MH part of the sampling. The cluster distribution was calculated by using at least 103 configurations that were sampled with a time step of 0.1 ps. Finally, all our results were obtained by using the DL-POLY2 package.37 In that respect, we note that the long-range part of each of the like-like sites SSRDF that are the outcome of this package has been here corrected for an unjustified (N - 1)/N (where N is the total number of sites) factor that is apparently
TABLE 1: Force Field Parameters for the Amides Simulated in This Work FA mol-1
(kJ σ (nm) q (e)
)
H
O
C
N
Hc
Ht
0.1591 0.275 0.12
0.8792 0.296 -0.46
0.4396 0.375 0.34
0.7117 0.325 -0.83
0.000 0.000 0.415
0.000 0.000 0.415
H
O
C
N
Ht
M
0.8792 0.296 -0.46
0.4396 0.375 0.34
0.7117 0.325 -0.70
0.000 0.000 0.415
0.543 0.38 0.285
NMF (kJ mol-1) 0.1591 σ (nm) 0.275 q (e) 0.12
DMF H -1
(kJ mol ) 0.1591 σ (nm) 0.275 q (e) 0.12
O
C
N
Mc
Mt
0.8792 0.296 -0.46
0.4396 0.375 0.34
0.7117 0.325 -0.57
0.543 0.38 0.285
0.543 0.38 0.285
O (kJ mol-1) σ(nm) q (e)
0.650 0.3165 -0.8476
SPC/E H 0.000 0.000 0.4238
H 0.000 0.000 0.4238
inherently implemented in this package. This correction is essential to recover the correct asymptotic behavior of the RDF.38 3. Results As hinted in the Introduction, the purpose of the present calculations is not to (re)model amides to reproduce (better) their known physical properties, as is often done in the simulation literature. This type of modeling is of course necessary to achieve the necessary reliability; however, it has been recently admitted that aqueous mixtures need a supplementary modeling ingredient, which is the reproduction of the correct microheterogeneity, namely through the reproduction in the correct integrals of the RDFs. It turns out that this step is rather a crux because it couples inextricably the microscopic pair interaction to the nanometer scale structure and leads to both a considerable slowing down of the dynamics and the need for large scale simulations. Needless to say, any data related to short-range feature of the structure will be nearly insensitive to the nanoscale structure. This concerns particularly the internal energy.39 Indeed, the statistical formulation of the energy involves the integral of the product of the pair correlation function with the pair interaction; for polar molecules, this product is always short ranged over few molecular diameters; hence, the integrand is slightly sensitive to the structure of the MH beyond this range. Table 2 shows the results of our simulations for neat and equimolar aqueous amide mixtures, together with the experimental values. The agreement is relatively good, which confirms that the models for amides are suited to our investigations. 3.1. The Structure of Neat Amides. 3.1.1. The Pair Correlations in Neat Amides. Figure 1 shows selected sitesite RDF for FA (top panel), NMF (middle panel), and for DMF (lower panel). Because the vertical scale is the same in all plots, we can compare the SSRDFs directly. The center-of-mass RDF is equally displayed in all plots and serves as a measure of the packing effects of the liquid state. It also serves as a measure of which sites govern the local packing effects. Indeed, for FA it is the carbon-carbon site correlations that lead the packing effects, while it is the nitrogen-nitrogen and oxygen-oxygen
MD Study of Microheterogeneity in Neat and Aqueous Amides
Figure 1. Selected site-site radial distribution functions for neat formamide (top panel), NMF (middle panel), and DMF (bottom panel). Carbon-carbon RDFs are in magenta, nitrogen-nitrogen are in green, oxygen-oxygen are in cyan, and oxygen-nitrogen (related to Hbonding) are in red. The center of mass RDF is shown in black.
TABLE 2: Thermodynamic Properties of Neat and Equimolar Aqueous Amidesa FA mol-1
∆H (kJ ) Vm (cm3 mol-1) F (g cm-3)
FA-W
Expt
Sim
Expt
Sim
-61.546 39.89 1.13
-61.897 41.32 1.09
29.82 1.09
-54.31 29.76 1.05
NMF -1
∆H (kJ mol ) Vm (cm3 mol-1) F (g cm-3)
NMF-W
Expt
Sim
-57.65 57.56 1.02
-57.23 61.77 0.96
DMF mol-1
∆H (kJ ) Vm (cm3 mol-1) F (g cm-3) a
Expt
Sim -54.51 39.43 0.98
DMF-W
Expt
Sim
Expt
Sim
-46.47 77.66 0.94
-46.27 79.43 0.92
45.75 0.976
-48.64 47.94 0.949
Experimental values are from refs 14, 41.
for the NMF and the nitrogen-nitrogen for the DMF. These differences are related to the local preferred relative orientations of the molecules, which in turn are related to the H-bonding differences. Indeed, the C-C correlations vary from monopeaked to bipeaked configuration when going from FA to DMF through NMF. Similarly, the O-O correlations are monopeaked for NMF but bipeaked for both FA and DMF. Because double peaks indicate the existence of two preferred orientations of near neighbors, we see that for the case of FA that it is the dual H-bonding possibility between O and N atoms that creates the double-peaked structure of these two peaks, while for DMF it is the dual orientations (C or N end parts) that create the double peak. Both peaks are asymmetric, but the H-bonding for FA is more asymmetric than the neutral orientation-based double peak
J. Phys. Chem. C, Vol. 111, No. 43, 2007 15589
Figure 2. Selected site-site structure factor for neat formamide (top panel), NMF (middle panel), and DMF (bottom panel). The color conventions are as in Figure 1.
in the case of DMF. This asymmetry for FA has been equally noticed and commented on in terms of pair energy probability distribution.12 As far as concerns of the H-bonding correlations displayed through the O-N RDFs, it is quite apparent that FA and NMF amides tend to H-bond through nitrogen-oxygen, that is in the N-H·‚‚O configuration, as expected. The hydrogen in the hydoxyl H-CdO end group of all three amides is a rather weak hydrogen bond. The status of this H atom is different in the OPLS models,11 where the charge is entirely neglected and the Cordeiro models14 where it is made explicit. Lei et al.40 have investigated the OPLS model and the status of the weak H-bonding, which is an issue equally addressed by Cordeiro.14 The partial structure factors corresponding to all the correlations of Figure 1 are displayed in Figure 2. In accordance to what is seen for the corresponding RDFs, we see that the C-C, N-N, and to some extent the O-O structure factors indicate more specifically the packing structure of the molecules, because these structure factors do look like those of very dense LennardJones liquids. This confirms the disordered structure of the overall molecular distribution through the inert sites. The first peak of the structure factor SCC(k) indicates the packing modulation at wave-vector km ) 2π/σm, where σm is the spherical diameter of each molecule. Conversely, the O-N structure factor (in red) shows a specific H-bonding structures, except of course for DMF. It can be seen that there is a notable tendency to exhibit a prepeak that is at wave vector smaller than km for respective amides, which is about k ) 1.32 Å-1 for both FA and NMF. DMF shows a very weak shoulder. This prepeak corresponds to some local structure that is driven by the H-bond. We have shown recently7,8 that in alcohols such a structure had specific morphology related to the topology of H-bonds (chain, ring,...). In the present case, the prepeaks are rather weak and very much reminiscent of that observed for water,8 which is a strong H-bonding liquid. It is important to realize that these prepeaks correspond to structures in the RDFs that are buried
15590 J. Phys. Chem. C, Vol. 111, No. 43, 2007
Zoranic´ et al.
Figure 3. (a) Cluster probabilities for neat formamide. Left panel is for like sites (lC ) 4.0 Å). Right panel is for unlike sites (main panel lC ) 2.9-3.2 Å and inset lC ) 3.1 Å). (b) Cluster probabilities for neat NMF (right panel) and DMF (left panel). Main panels display like clusters (lC ) 4.5 Å for NMF and lC ) 5.0 Å for DMF). Insets: (lC ) 3.0 Å for NMF and lC ) 3.0-4.0 Å for DMF).
after the first peaks, such that the calculation of the structure factors is an important way to address the question of these local structures. 3.1.2. Cluster Distribution in Neat Amides. Because we have shown and wondered about the absence of specific cluster morphology in water,7,8 it is interesting to see if this is the case for FA and NMF. For that, we have calculated the clusters directly through eq 3. The results are shown in Figure 3a for formamide and Figure 3b for NMF (left panel) and DMF (right panel). First of all, let us recall what the cluster analysis reveals in the case of a simple Lennard-Jones liquid. One finds that the probability for monomer is nearly unity, and the cluster distribution decreases rapidly (almost exponentially) with higher cluster sizes, as one expects for any unassociated liquid. The same analysis for water revealed a distribution similar to that of a simple liquid in contrast to alcohols, which showed a small peak around cluster size 5.7,8 A rapid look at Figure 3a,b shows that we have a cluster distribution analogous to that of unassociated Lennard-Jones type liquid. In other words, the cluster distribution of a highly associated liquid such as FA is found to be analogous to that of DMF that is closer to a simple polar fluid because it is a weakly hydrogen-bonded liquid. This is more surprising than the results we have previously found
for water7,8 and reveals that the clustering in water and formamide is very complex and goes undetected by a straightforward cluster analysis. Let us now look into some detail the clustering differences between different sites. In all cases, the lC parameter in eq 3 was taken to be between the first peak and the first minimum of the corresponding site-site RDF. For formamide, the left panel of Figure 3a shows the self-site-site clustering and the clustering related to H-bonding. The main panel of left side shows that the carbon-carbon cluster distribution is identical to that of the center of masses. The inset shows the self-cluster distribution of oxygen, nitrogen, and that of the hydrogen linked with the carbon atom, and this over the entire range of cluster sizes (with N ) 256 molecules). It reveals again that there is no small preferred cluster size, while there is some preferred clustering at sizes close to that of the simulated system, which is always of very small probability (note the log scale in the inset). In the right panel, the O-N clustering, which should show a peak as formamide having some preferred cluster size, is shown in the main panel for several lC parameter values, varied between distances before the first peak and after the first minimum. It is very clear that there is no preferred clustering, because monomer probability is always the highest, and this
MD Study of Microheterogeneity in Neat and Aqueous Amides finding is robust with respect to cluster definition. The inset displays cluster distributions between the oxygen atom and the cis-, trans-, and carbon-attached hydrogen sites. Again, no specific clustering is observed in the small cluster region. In Figure 3b, the left (NMF) and right (DMF) panels again confirm the absence of preferred clustering. For NMF, the main panel shows that there is some difference between the methylmethyl and center of masses cluster distribution, while the former is nearly identical to the N-N and O-O distributions. The inset shows the large differences between the clusters of the two hydrogen. For DMF, the right panel shows the relative isotropy of the molecule with respect to cluster formation because all site-site probabilities are very similar, indicating that this molecule is a polar entity with very weak H-bonding. This confirms the absence of specific clustering for DMF but raises the question of the type of clustering for FA and NMF. As for the case of water, we conclude that amides have a cluster morphology that has no simple topology and is rather branchedlike so that it escapes a direct determination. It is notable that we could not isolate any specific cluster shapes, as for water. It is important to note that some previous studies14 have insisted on dimer formations in chain type clustering for amides as a sign of specific clustering. Our analysis shows that dimer formation is not anything specific because no peak is seen at cluster size 2, but rather that the entire distribution is very much like that of a simple liquid. Similarly, if chain formation is occurring in amides, it has no specific size, similar to what we found for methanol.8 3.2. The Structure of Aqueous Amides. Because both water and amides have nontrivial clustering, it is interesting to see how well they can mix. Indeed, because FA and NMF can H-bond to water and are rather weakly polar molecules, one expects good mixing conditions for these solutes. A look at the experimental KBIs is a good place to start to get an idea of the mixing features of aqueous amides. 3.2.1. The Experimental Kirkwood-Buff Integrals. The KBI of aqueous amides have been obtained experimentally by Marcus42 and are displayed in Figure 4 as symbols. These data deserve some remarks. First of all, the expressions of the KBIs for a binary mixture, in terms of the molar volume V, partial molar volumes V h i, isothermal compressibility χT, and density derivatives of the chemical potentials µi (i ) 1,2 being the component indices) can be cast as follows:
(
Gij ) - G12δij + (1 - δij) kBTχT -
) (
)
V h 1V δij V h2 hi + - V (4) VD xi D
where D ) (ni/nj)(∂µi/∂ni)T,P,nj, ni, nj are the number of moles of species i and j, and δij is the Kronecker symbol. In practice, this expression can be simplified by two considerations. First, that the compressibility of dense liquids is nearly zero (and this term is always of orders of magnitude smaller than all others), and second that the volume is a linear function of x the solute (or symmetrically the solvent) mole fraction V(x) ) (1 - x)V(0) + xV(1). In the latter case, one has simply V h i ) V(xi ) 0). This latter approximation is generally quite good as soon as the volumes of the two components differ by more than a factor two. Table 2 shows that this is true for amides because the molar volume of water is about 18 cm3/mol. There is a third approximation that is to assume that the mixture is ideal, meaning that the excess chemical potentials are zero, which means D ) 1 in the expressions above. This last approximation is much more difficult to justify, because
J. Phys. Chem. C, Vol. 111, No. 43, 2007 15591
Figure 4. Experimental KBIs for amides: top panel for FA, middle for NMF, and bottom for DMF. The symbols are from ref 42: in each panel Gww is shown in open dots, Gws in crosses, and Gss in squares. The lines are the ideal mixing KBIs discussed in the text: dashes for Gww, dots for Gws, and full line for Gss.
even if there is no direct interaction between the two components they still interact through the excluded volume effects. It turns out that in practice this is often a good approximation when two components mix well, which is to say they do not show any sorts of microheterogeneity. We realize that this is a rather stringent approximation because we have hinted all along that we have MH in aqueous amides, as indeed is the case. If we nevertheless admit these hypotheses, then the KBIs become very simple functions of the solute mole fraction x. These functions have been drawn as lines in Figure 4, along with the experimental data. A striking fact is then observed: for DMF and NMF, the simplified Gss and Gws (where s stands for solute and w for water) are nearly superposed to the experimental data. This fact calls for several inevitable conclusions. First of all, for DMF and NMF, because D is the same for all three KBIs, the near good agreement mentioned above means that the same agreement should equally be observed for Gww and therefore the experimental Gww are erroneous. The exact reasons for possible errors can be deduced from the form of the KBI in eq 4. Indeed, it can seen that Gww is apparently singular at x ) 1 (the same way Gss is singular at x ) 0). These singularities are only apparent, as can be seen from the simplified expressions, because the form of V avoids these. Hence, the experimental data that is used to obtain V (or the fitting expressions if any) do not have the appropriate limiting behavior near x ) 0 and x ) 1 endpoints. This type of error is often achieved when inaccurate calculation of partial molar volumes are conducted numerically from experimental molar volumes.28 Figure 4 shows for DMF a comparison between the KBIs obtained from constant partial molar volumes as assumed above
15592 J. Phys. Chem. C, Vol. 111, No. 43, 2007
Figure 5. Selected site-site radial distribution functions for equimolar aqueous FA, NMF, and DMF. Top panel, gww(r); middle panel, gwC(r); bottom panel, gCC(r). Line conventions: full for FA, dotted for NMF, and dashed for DMF. Each data contain two sets corresponding to two successive runs (see text).
and from accurate experimental data from ref 41. The difference is hardly seen from bottom panel of Figure 4, which proves indirectly that accurate KBIs can be obtained with linear behavior for V.28 With this explanation in mind, it is instructive to now examine the experimental KBIs for formamide in the top panel of Figure 4. It seems quite clear that both Gww and Gss have very similar errors at their respective pseudo-singular endpoints. The singularity in Gws is harder to understand and must be attributed to erroneous manipulation of data. Can we guess the correct KBIs for formamide? Knowing that DMF and NMF are more hydrophobic than FA and that these two components are quite accurately ideal, then it seems quite reasonable to suppose that FA should be even more ideal. This leads to accept that the approximate KBI drawn as lines should be quite close to the real ones. Once these points are settled, it still remains to explain why the KBIs for DMF (and to a lesser extent that of NMF and FA) are so nearly ideal, which is to say that D ) 1 for all x. 3.2.2. The Pair Correlation Functions of the Aqueous Amides. Figure 5 shows selected SSRDF for aqueous mixtures at equimolar mixture. Top panel shows the water-water RDFs, the middle panel the water-carbon site RDF, and the lower panel shows the carbon-carbon RDF. For each of the amides, RDFs obtained from two successive runs separated by 50 ps are shown. Because these are almost indistinguishable, one may be tempted to assume that the configurations are well converged. The most notable features between all three RDFs is that gww(r) has the highest peaks, while gwC(r) and gCC(r) have more pronounced oscillatory structure at next neighbors. Let us first analyze water-water correlations. The main peak is higher for DMF and in decreasing order when going from DMF to FA. Knowing that the height of the first peak is about 3 for pure SPC/E water, we conclude that the height of the peak corre-
Zoranic´ et al. sponds to the increase in the strength of the water-water near neighbor correlations, which are higher when the body is more hydrophobic. This means that water mixes better with FA than with DMF. On the other hand, it is quite remarkable that the next neighbor and beyond structure is nearly the same for all three amides with however a small second peak increase for DMF. This is certainly indicative of the particularity of the water structure in the presence of a solute. Next, we analyze the solute-solute correlations from the lower panel. Using a similar argument about differences in first peak heights, we can say that higher pairing probability will translate into a higher peak and thus confirm what is already known through direct H-bonding analysis, namely that formamide H-bonds more than NMF, which in turn H-bonds more than DMF. The double peak structure observed for DMF is similar to that of the neat DMF (see Figure 1c). Perhaps the most important feature is that all three C-C correlations are very similar to that of the respective neat amides. This is a surprising result that indicates clearly that the amides are segregated between them such that their first neighbor correlations are quite similar than those in the neat liquids. The above conclusion is correlated by the water-carbon site correlation in the middle panel of Figure 5. We see indeed that the first peaks are surprisingly closer in height and width for all three amides, at least comparatively to the differences observed in the two other cases. This indicates that the amide water contact is governed by the interface between the segregated domains and thus is nearly insensitive to the interaction details between the different amides. The other striking feature is that the water-NMF peak is higher than water-FA. Indeed, from the arguments developed in the paragraphs above, we should expect peak heights to diminish from water-FA to water-DMF. The fact that NMF has more affinity toward water is rather surprising. However, a look at the correlations at larger distances reveals important differences. First, NMF has almost no secondary structure, which indicates that only high first neighbor correlations are important, and supports the idea of a strongly microsegregated situation with sharp interface. On the contrary, water-DMF has a more pronounced secondary structure, which indicates alternate layering and enforcing the idea of a better mixing behavior. Formamide is somewhat intermediary. The fact that NMF appears more segregated with water has an experimental support, namely that it is known that NMF does not mix with water as well as FA and DMF.14,43 To confirm the global picture that emerges from our analysis above, it is necessary to confront this to the direct calculation of the KBIs from the simulation SSRDFs. 3.2.3. Kirkwood-Buff Integrals from Simulations. Figure 6 shows the running KBIs (RKBIs) corresponding to the RDFs shown in Figure 5. The running KBIs are defined as Gij(r) ) 4π ∫r0 dt t2(gij(t) - 1). From eq 1, one has Gij(r f ∞) ) Gij; hence, the asymptotic behavior of the RKBIs should provide an indication about the calculated KBIs. For this purpose, the ideal KBIs obtained in Figure 4 are displayed as a horizontal line. It is important to note that the KBIs obtained through any pair of sites on respective components are identical to that of the corresponding component pair.2 As in Figure 5, the RKBIs corresponding to successive calculations have been displayed and can be clearly observable in Figure 6. It is seen that there is not much variation between the two runs within the intermolecular distance range smaller than some range that is smaller for water than for amides. For the latter, this distance is about 1 nm. For water, this range varies with amides: it is about 2.7 Å for FA, which is in fact the first peak for gww(r),
MD Study of Microheterogeneity in Neat and Aqueous Amides
Figure 6. KBI corresponding to the SSRDF shown in Figure 5. The horizontal lines correspond to the ideal mixing KBI shown in Figure 4 for all three mixtures. Line color conventions as in Figure 5.
and gets larger as amides have less H-bond ability because it is about 5 Å for DMF. We see that because the RDFs in Figure 4 did not show any variation does not mean that there is not some evolution between two runs, which is a direct consequence of the evolution of the MH at a larger scale than tens of ps. These two observations indicate the nature of the respective microheterogeneity of water and amides in these mixtures. We come back to this point below through the analysis of the snapshots. The next most apparent feature is the differences in the asymptotic behavior of the RKBIs. First of all, it is quite apparent that the amide-amide KBIs (shown in magenta) seem to converge to the corresponding experimental values, in contrast with the water-water ones (and consequently the water-amides ones, because these two are obviously coupled). This convergence is relatively satisfactory for DMF and NMF, knowing
J. Phys. Chem. C, Vol. 111, No. 43, 2007 15593 that the Cordeiro force fields have not been adjusted to reproduce the experimental KBI. In passing, this fact indicates that this latter fitting procedure is not really required, and that the true problems lie elsewhere. This convergence is not so good for FA, which could mean several things: either the mixture is not ideal or the force field parameters are not good, or even larger simulation statistics are needed. Turning to the RKBIs for water, one sees that there is clearly no sign that the waterwater correlations are converged. On the other hand, two successive runs are not that far apart to indicate that the system is not stabilized yet. Our conclusion is that water is more strongly microheterogeneous than amides, and therefore one needs to perform statistics over the slow evolution of the microheterogeneity itself. This procedure could take unreasonably long times such as of the order of several nanoseconds.20,16 Because the inherent time scale of the MH is not known a priori, one may wonder if in addition the distance scale used here (box size about 50 Å) is enough to properly stabilize the asymptotic part of the KBI. It is clear that much larger system sizes would take even longer times to sample properly, perhaps tens or hundreds of nanoseconds. With that respect, the asymmetry of the MH of both components is quite remarkable. To confirm this difference in MH, we have displayed in Figure 7a snapshot of the equimolar water and formamide configurations and similarly in Figure 8 equimolar water and DMF. In each shot, the adverse species is shown in skeletal representation, which allows us to appreciate the relative spatial distributions. Although snapshots are only partial representations of one realization, it is quite clearly seen that water is more heterogeneously distributed and looks nearly identical in both pictures, despite the fact that one should naively expect it to be more homogeneously distributed in the case of mixing with FA because of the more possibilities to H-bond. However, it is clear that this is not the case, and we will accept here and for future issues that water is invariably microheterogeneous in mixing conditions and whatever the solute. We find that water looks very spongelike in structure, and we have observed this in many other mixing conditions. This is consistent with what is found and claimed in ref 5 in the case of aqueous methanol. In contrast to water, both amides look more homogeneous, and an inspection of the snapshots indicates that amides occupy both the voids left by water, as well as inside water microdomains. This is where both amides differ, because FA tends to occupy water
Figure 7. Snapshots of aqueous formamide equimolar mixture. Left and right shots are for water and formamide molecules, respectively, and in each of which the alter species is shown in skeletal representation. The color codes for atoms are O (red), H (white), N (blue), and C (cyan).
15594 J. Phys. Chem. C, Vol. 111, No. 43, 2007
Zoranic´ et al.
Figure 8. Snapshots of aqueous DMF equimolar mixture. Left and right shots are for water and DMF molecules, respectively. The representations and color code conventions are as in Figure 7, and the methyl groups are in cyan.
domains more than DMF does. With that respect, DMF is slightly more microheterogeneous than FA. The observation of the microstructure of aqueous amide has two far-reaching consequences. First, it seems almost impossible to claim a proper statistical sampling in the range of the MH, which is to say from 15 Å and beyond. The main reason for this is as follows. The MH is quite consistent in shape, Swiss cheese or spongelike, and the sampling through all possible variations of such shape seems to require both very large systems and long times, the latter of which is definitely required to sample all possible variation of shape within a given range (e.g., in the 15-25 Å range). It seems quite improbable to achieve this within a few nanoseconds, as claimed in some cases.20,16 Probably, a proper sampling would take few hundreds of nanoseconds if not few microseconds. In such case, the problem of propagation of numerical errors seems unavoidable. The second consequence is related to the puzzle of the near ideal KBIs while these systems are clearly microsegregated. We are led to conclude that the KBIs are not a good indicator of the microheterogeneity. This conclusion is, however, not entirely satisfying, because Matteoli and Leopori26 have shown that solute molecules possessing hydrophobic groups, such as acetone and tetrahydrofuran for example, tend to exhibit high KBIs, which is interpreted as the sign of MH. One possible way to resolve this puzzling issue is to realize that ideal mixing does not necessarily imply homogeneity, hence D ) 1 in the expression eq 4 for the KBIs, despite the fact that the mixture is heterogeneous. In other words, there may be compensating effects in the chemical potentials that make their density derivatives nearly zero. Physically, this interpretation means that the system restructures itself locally to minimize the response to density fluctuations, thus introducing the microheterogeneity. Despite this argument, the relationship between the KBIs and the MH remains to be settled. 4. Conclusion In this study, we have focused on the microstructure of neat amides and the microheterogeneous structure of aqueous amides. Three amides were studied: FA, which is analogous to water in that it is both donor and acceptor of H-bond and can form microemulsions, NMF, and DMF. For neat amides, our study has indicated that no specific morphology of cluster arises for FA and NMF, hinting some analogy with similar behavior in water. As for water, we expect the clustering to be polymeri-
like and possibly branched in such a way that the resulting clusters have always a probability lower than the monomer. This fact also implies a strong dynamics/kinetics of the clustering, because monomers are found to have the largest probability, and yet the liquid is very associated. For aqueous amides, we have found that water displayed a very similar spongelike structure for all three aqueous amide mixtures, apparently insensitive to the differences between the force fields of the three amides. This is equally observed for many other aqueous mixtures with very different polar solutes.44 Conversely, the microheterogeneous structure of the three amides is found to be somewhat different. The description of the MH through the calculation of the KBI reveals that the proper calculation of the range of the RDF in the microheterogeneous part is not properly converged. This lack of convergence is herein attributed to the inherent difficulty of sampling across the larger level of disorder of the MH, and we have argued that it was more reasonable to consider a time scale of tens of nanoseconds, with the inherent numerical errors due to computer statistics on such scale of time. This argument poses the limits of computer simulations of aqueous mixtures, at least the obtention of proper structural information beyond the nanometer scale. Because amides are quite simple solutes that do not have large KBIs and hence are supposed to be quite homogeneously soluble in water, it comes as a surprise that even this type of solute should have microsegregation problems that are expected with solutes with larger hydrophobic groups. The results displayed in this work tend to point toward water as the universal MH maker, rather than the solutes themselves. In contrast, the KBIs are an indicator of global concentration fluctuations. This study opens several questions about the way concentration fluctuations could be coupled to MH, which deserves clear answers for aqueous solutions to be fully understood. Supporting Information Available: Centre de Calcul de l’Universite´ Pierre et Marie Curie is acknowledged for generous allocation of computer time. L.Z. gratefully acknowledges the financial support from The National Foundation for Science, Higher Education and Technological Development of Republic of Croatia Project Number 03.01/05 . References and Notes (1) Gray, C.; Gubbins, K. Theory of Molecular Fluids; Oxford University Press: London, 1984.
MD Study of Microheterogeneity in Neat and Aqueous Amides (2) Hansen, J. P.; McDonald, I. R. Theory of Simple Liquids; Academic: London, 1986. (3) Stanley, H. E.; Cruz, L.; Harrington, S. T.; Poole, P. H.; Sastry, S.; Sciortino, F. F.; Starr, W.; Zhang, R. Physica A 1997, 236, 19. (4) Errington, J. R.; Debenedetti, P. Nature 2001, 409, 318. (5) Dixit, S.; Crain, J.; Poon, W. C.; Finney, J. L.; Soper, A. K. Nature 2002, 416, 829. (6) Dougan, L.; Bates, S. P.; Hargreaves, R.; Fox, J. P.; Crain, J.; Finney, J. L.; Reat, V.; Soper, A. K. J. Chem. Phys. 2004, 121, 6456. (7) Perera, A.; Sokolic, F.; Zoranic, L. Phys. ReV. E 2007, 75, 060502(R). (8) Zoranic, L.; Sokolic, F.; Perera, A. J. Chem. Phys. 2007, 127, 024502. (9) Koga, Y. J. Phys. Chem. 1996, 100, 5172. (10) Schubert, K. V.; Buse, G.; Strey, R.; Kahlweit, M. J. Phys. Chem. 1993, 97, 248. (11) Jorgensen, W. L.; Swenson, C. J. J. Am. Chem. Soc. 1985, 107, 569. (12) Elola, M. D.; Ladanyi, B. M. J. Chem. Phys. 2006, 125, 184506. (13) Zielkiewicz, J. Phys. Chem. Chem. Phys. 2000, 2, 2925. (14) Cordeiro, J. M. Int. J. Quantum Chem. 1997, 65, 709. (15) Cordeiro, M. A. M.; Santana, W. P.; Cusinato, R.; Cordeiro, J. M. M. THEOCHEM 2006, 759, 159. (16) Kang, M.; Smith, P. E. J. Comput. Chem. 2006, 27 (13), 1477. (17) Gao, J.; Pavelites, J. J.; Habibollazadeh, D. J. Phys. Chem. 1996, 100, 2689. (18) Krienke, H.; Fisher, R.; Barthel, J. J. Mol. Liq. 2002, 98-99, 329. (19) Skarmoutsos, I.; Samios, J. Chem. Phys. Lett. 2004, 384, 108. (20) Chitra, R.; Smith, P. E. J. Chem. Phys. 2001, 114, 426. (21) Lee, M. E.; van der Vegt, N. F. J. Chem. Phys. 2005, 122, 114509. (22) Sokolic´, F.; Idrissi, A.; Perera, A. J. Chem. Phys. 2002, 116, 1636. (23) Perera, A.; Sokolic´, F. J. Chem. Phys. 2004, 121, 11272. (24) Kirkwood, J. G.; Buff, F. P. J. Chem. Phys. 1951, 19, 664.
J. Phys. Chem. C, Vol. 111, No. 43, 2007 15595 (25) Ben-Naim, A. J. Chem. Phys. 1977, 67, 4884. (26) Matteoli, E.; Lepori, L. J. Chem. Phys. 1984, 80, 2856. (27) Perera, A.; Sokolic, F.; Alma`sy, L.; Westh, P.; Koga, Y. J. Chem. Phys. 2005, 123, 024503. (28) Perera, A.; Sokolic, F.; Alma`sy, L.; Koga, Y. J. Chem. Phys. 2006, 124, 124515. (29) Bathia, A. B.; Thornton, D. E. Phys. ReV. B 1970, 2, 3004. (30) Teubner, M.; Strey, R. J. Chem. Phys. 1987, 87, 3195. (31) D’Arrigo, G.; Giordano, R.; Teixeira, J. Eur. Phys. J. E 2003, 10, 135. (32) Stillinger, F. H. J. Chem. Phys. 1963, 38, 1486. (33) Pugnaloni, L. A.; Vericat, F. J. Chem. Phys. 2002, 116, 3. (34) Senger, B.; Schaaf, P.; Corti, D. S.; Bowles, R.; Pointu, D.; Voegel, J.-C.; Reiss, H. J. Chem. Phys. 1999, 110, 6438. (35) Berendsen, J. C.; Postma, J. P. M.; Von Gusteren, W. F.; Hermans, J. In Intermolecular Forces; Pullman, B., Ed.; Reidel: Dortrecht, 1981. (36) Jorgensen, W. L.; Chandrasekhar, J.; Madura, J. D.; Impey, R. W.; Klein, M. L. J. Chem. Phys. 1983, 79, 926. (37) DL_POLYpackagefromCCLRCDaresburyLaboratory.www.ccp5.ac.uk/ DL_POLY/. (38) Perera, A.; Zoranic´, L.; Sokolic´, F. Unpublished, available through http://arxiv.org/abs/0708.2346. (39) Zoranic, L.; Perera, A.; Sokolic, F. J. Mol. Liquids, submitted for publication, 2007. (40) Lei, Y.; Li, H.; Pan, H.; Han, S. J. Phys. Chem. A 2003, 107, 1574. (41) Torres, R. B.; Marchiore, A. C. M.; Volpe, P. L. O. J. Chem. Thermodyn. 2006, 38, 526. (42) Marcus, Y. J. Chem. Soc., Faraday Trans. 1990, 86 (12), 2215. (43) Lide, D. R. CRC Handbook of Physics and Chemistry, 73rd edition; CRC Press: Cleveland, Ohio, 1992. (44) Perera, A.; Zoranic´, L.; Sokolic´, F.; Mazighi, R. To be submitted for publication.