This is an open access article published under an ACS AuthorChoice License, which permits copying and redistribution of the article or any adaptations for non-commercial purposes.
Article Cite This: ACS Omega 2018, 3, 1060−1068
Generic Mechanism for Pattern Formation in the Solvation Shells of Buckminsterfullerene Sanwardhini Pantawane,†,∥ Dibyendu Bandyopadhyay,‡ and Niharendu Choudhury*,†,§ †
Theoretical Chemistry Section and ‡Heavy water Division, Bhabha Atomic Research Centre, Mumbai 400085, India § Homi Bhabha National Institute, Anushaktinagar, Mumbai 400094, India S Supporting Information *
ABSTRACT: Accurate description of solvation structure of a hydrophobic nanomaterial is of immense importance to understand protein folding, molecular recognition, drug binding, and many related phenomena. Moreover, spontaneous pattern formation through self-organization of solvent molecules around a nanoscopic solute is fascinating and useful in making template-directed nanostructures of desired morphologies. Recently, it has been shown using polarizable atomistic models that the hydration shell of a buckminsterfullerene can have atomically resolved ordered structure, in which C60 atomic arrangement is imprinted. In analyzing any peculiar behavior of water, traditionally, emphasis has been placed on the long-ranged and orientation-dependent interactions in it. Here, we show through molecular dynamics simulation that the patterned solvation layer with the imprints of the hydrophobic surface atoms of the buckminsterfullerene can be obtained from a completely different mechanism arising from a spherically symmetric, shortranged interaction having two characteristic lengthscales. The nature of the pattern can be modified by adjusting solvent density or pressure. Although solute−solvent dispersion interaction is the key to such pattern formation adjacent to the solute surface, the ordering at longer lengthscale is a consequence of mutual influence of short-range correlations among successive layers. The present study thus demonstrates that the formation of such patterned solvation shells around the buckminsterfullerene is not restricted to water, but encompasses a large class of anomalous fluids represented by two-lengthscale potential.
1. INTRODUCTION Understanding many important phenomena like protein folding/misfolding, drug binding, molecular recognition, and micelle and bilayer formation requires accurate description of water structure around a nanoscopic solute. Recently, carbon nanotubes (CNT)1−3 and buckminsterfullerene (C60)4−9 have attracted huge attention in the field of nanomaterials and nanomedicines because of their unique structures and nonreactive nature. The use of fullerenes as biosensors8,9 and drugdelivery agents4,5 requires knowledge on the behavior of fullerenes in their aqueous solution. Moreover, as proteins have nanoscopic, hydrophobic patches in contact with water, the nanoscopic fullerene−water or CNT−water interface is generally used as a model system in the investigations of hydrophobic hydration and dewetting-induced aggregation.10−20 In depicting the solvation structure of a fullerene or a CNT, the azimuthally averaged radial density of the solvent particles around the solute is generally evaluated and thus the effect of detailed atomic arrangement of the surface atoms of the solute on the solvation structure is neglected. In a recent study, using molecular dynamics (MD) simulation, Levitt and co-workers,21 however, have observed azimuthally resolved C60 hydration shell with remarkable patterns arising due to selforganization of water molecules around the fullerene molecule. Although pattern formation22−27 through self-organization is well known in Langmuir monolayers,24 lipid monolayers,25 and polymer films,26,27 the same at the nanoscale can have many © 2018 American Chemical Society
technological applications in nanolithography and other related techniques.28 The formation of a self-organized structure with specific patterns is generally a consequence of a competition between two different interactions.29 Many a time, this competition arises between short-range attractive interactions and long-ranged Coulombic repulsions, as in the case of magnetic materials.30 When a solute with a specific surface structure is immersed in a solvent, arrangement of the solvent molecules surrounding the solute is governed by a competition between solute−solvent and solvent−solvent interactions. The patterned hydration layers of the fullerene molecule have been obtained21 using highly sophisticated polarizable atomistic models for water and C60. It has been shown that the pattern formation in the immediate vicinity of the fullerene molecules is due to solute−water attractive dispersion interaction and that the same at longer lengthscale is due to long-range interaction arising from polarization.21 As the water molecules are polar and have orientation-dependent interactions, it is natural to associate31,32 such pattern formation with the specific properties of water. It is therefore challenging and intriguing to examine whether a spherically symmetric model with two characteristic lengthscales can yield such a patterned solvation structure. It is now well researched that a spherically symmetric Received: November 27, 2017 Accepted: January 16, 2018 Published: January 26, 2018 1060
DOI: 10.1021/acsomega.7b01858 ACS Omega 2018, 3, 1060−1068
Article
ACS Omega
Figure 1. Short-range ADFs for local density of waterlike solvent in front of the first solvation shell of C60. The 60 carbon atoms of C60 are shown as black dots. Solvent densities are calculated by considering waterlike solvent particles in a spherical layer of width 0.1 Å around the C60. The middle of the spherical layer is at 6.63 Å (r* = 2.09) (corresponding to the first vertical dotted line in Figure 2C) from the C60 center, and the local waterlike solvent density patterns shown in (A)−(F) are for different bulk densities of the CG waterlike solvent. The pattern formation due to preferential solvation of the C60 surface is evident. (A) For bulk density ρ0 = 0.07, local waterlike solvent density is only on the 20 hexagonal faces, whereas all of the 12 pentagonal faces of C60 remain empty, giving rise to an icosahedral pattern of peaks.21 (B) At ρ0 = 0.16, all of the hexagonal and pentagonal faces are covered with CG waterlike solvent particles, but the density on the hexagonal face is higher than that on the pentagonal face. (C) At ρ0 = 0.24, both the hexagonal and pentagonal faces are covered with almost equal-density waterlike solvent leading to 32 vertex object known as a pentakis dodecahedron.21 (D−F) At higher densities, the pattern is different from that in (A)−(C). The ADF presented in the present study is scaled by the highest density.
properties of waterlike fluids reasonably well, we expect it to be successful in reproducing patterned solvation layers around a fullerene (C60) molecule immersed in a CS fluid. In essence, our main objective here is to verify whether a simple, spherically symmetric, two-lengthscale model for waterlike anomalous fluids can reproduce hydration patterns around the fullerene molecule as obtained from polarizable atomistic water model. Moreover, the long-range order observed by Levitt and co-workers21 was due to polarization, and it is therefore of interest whether the present nonpolar model can reproduce the long-range order. If such pattern formation at a longer lengthscale is observed for the present nonpolar model, then the mechanism of forming such long-range order will also be of interest.
core-softened (CS) pair potential with two lengthscales corresponding to a hard core and a soft corona can reproduce many anomalous properties of water. 33−42 Initially, a discontinuous water model with a hard core and a repulsive ramp, parameterized by Jagla,33,34 was shown to exhibit33−40 waterlike structural, thermodynamic, and dynamic anomalies.43 The continuous version of the CS potential, originally designed by Hemmer and Stell,41 and subsequently parameterized by many groups, has been shown to exhibit41,42,44−46 a waterlike cascade of structural, thermodynamic, and dynamic anomalies.43 Recent simulations demonstrated that the CS potential can give rise to different types of phase transitions, such as liquid−liquid phase transitions35,37 and isostructural solid− solid transitions,47 which are generally unusual for spherically symmetric one-component fluids. The CS model also captures beautiful stripe formation in a two-dimensional (2D) material.29 Not only that the bulk anomalous properties and unusual phase transitions are exhibited, waterlike solvation thermodynamics48 of nonpolar solutes has also been mimicked by the coarse-grained (CG) core-softened potential. More remarkably, the study48 has been successfully extended to polymeric hardsphere system, for which microscopic density fluctuations do not obey Gaussian statistics. Apart from reproducing bulk properties of waterlike fluids, the CS potential has recently been applied to mimic various nanoscopic hydrophobic solutes in water49−52 as well. Given that the spherically symmetric twolengthscale potential can predict bulk as well as interfacial
2. RESULTS AND DISCUSSION 2.1. Imprints of the C60 Atomic Structure on Its First Solvation Layer. It is quite common to consider the surface of a C60 molecule to be homogeneous53−56 and thus the conventional C 60 −solvent radial distribution function (RDF)56 is generally used to gauge the ordering of the solvent particles around the fullerene molecule. However, C60 has a truncated dodecahedron (soccer ball) structure, in which there are 20 hexagonal and 12 pentagonal faces. It is therefore of interest to know whether the present CS model of fluid can depict the effect of the truncated dodecahedron structure of C60 on its surrounding solvation layers. To account for this effect, 1061
DOI: 10.1021/acsomega.7b01858 ACS Omega 2018, 3, 1060−1068
Article
ACS Omega
Figure 2. (A) Two-dimensional contour map of interaction energy of a CG waterlike solvent particle placed at a radial distance r* = 2.09 with all of the 60 atoms of C60 (see Models and Methods). The black dots are the carbon atoms of C60. (B) Average number ⟨N⟩ of CG waterlike solvent particles in the first hydration shell (up to r* = 2.5) of C60 as a function of bulk density of CG waterlike solvent. (C) Radial distribution functions of CG waterlike solvent particles around the C60 molecule at T* = 0.25 at bulk densities (top panel) ρ0 = 0.07, 0.12, 0.16, and 0.24 and (bottom panel) ρ0 = 0.4, 0.5, 0.6, and 0.8 when an attractive van der Waals interaction is included in the interaction (see eq 2 with λ = 1) between the atoms of C60 and CG waterlike solvent particle. Distance r* in reduced unit (see Models and Methods) is measured from the center of the C60. Two vertical dotted lines are the positions around which a spherical layer of width 0.1 Å is chosen for calculating ADF. (D) Two-dimensional contour map of potential energy calculated by placing a CG waterlike solvent particle at all of the positions on a face of C60 and calculating the potential energy of this CG waterlike solvent particle with all of the 60 atoms of C60 as well as with 31 waterlike solvent particles occupying the centers of the remaining faces of the C60. In the inset, the energy map inside a six-membered ring is enlarged to show three low-energy wells (green spots) forming a triangular shape with each vertex directed toward the nearest six-membered ring. The black dots are the carbon atoms of C60.
we apply a newly developed21 concept of azimuthal distribution function (ADF), which is a modified version of the spatial distribution function (SDF) developed earlier57,58 (see Models and Methods). The ADFs of solvent particles in a thin layer of 0.1 Å thickness around the C60 (at a distance r* = 2.09 (6.625 Å) from the center of the C60) are shown in Figure 1 for different solvent densities at a fixed temperature of T* = 0.25. At the lowest density ρ0 = 0.07 considered here (Figure 1A), remarkable azimuthal ordering of solvent particles is observed with solvent particles preferentially occupying all of the hexagonal faces of C60, keeping all of the pentagonal faces unoccupied (20 bright spots in Figure 1A correspond to highdensity regions at the middle of 20 hexagonal faces). It is important to note that these high-density spots are at the middle of the hexagonal faces. At an increased density of ρ0 = 0.16 of the solvent, because of increased pressure, solvent particles occupy the middle of all of the pentagonal faces as well, but the density at the middle of the pentagonal faces is lesser than that on the hexagonal faces (see 20 brighter spots on the hexagonal faces compared to 12 not so bright spots on pentagonal faces in Figure 1B). Similar results have been observed by Levitt and co-workers21 with their sophisticated polarizable models for water and C60. With further increase in bulk density to ρ0 = 0.24, all of the hexagonal and pentagonal faces are occupied by the solvent molecules with almost equal probability (see 32 equally bright spots in Figure 1C). At ρ0 = 0.4, due to increased pressure, more and more solvent particles occupy positions corresponding to the middle of both the fiveand six-membered rings (Figure 1D). However, the solvent particles on the hexagonal faces are distributed over a region of triangular shape. At even higher densities (Figure 1E), similar picture remains, but now on each of the hexagonal faces there are three vertices (of a triangle) pointing toward other high-
density vertices in the neighboring six-member rings. We observe that the vertices are not directed toward the adjacent five-membered rings. At ρ0 = 0.8, both the pentagonal and hexagonal faces are occupied by the solvent particles. It is interesting to observe that even at this high density, the positions corresponding to all of the carbon atoms and C−C bonding regions of the C60 molecule remain unoccupied (Figure 1F) by the solvent particles. Thus, the imprint of the C60 surface structure is clearly visible to the solvation layer adjacent to the fullerene surface. It transpires from the above observations that at lower densities hexagonal faces are preferred to pentagonal faces by the solvent particles. To examine it, the solute−solvent potential energy map is calculated by placing a CG waterlike solvent molecule over the entire surface of a bare fullerene molecule (no solvent− solvent interaction is taken into account), and it is shown in Figure 2A. Darker and larger spots at all of the hexagonal faces in the energy map indicate (see Figure 2A) that the fullerene− solvent van der Waals attractive well is wider and deeper at the middle of the hexagonal face compared to that at the pentagonal face (see that the dark black spot is present at the middle of the hexagonal face but not at the middle of the pentagonal face). Average number of solvent particles ⟨N⟩ in the first solvation shell (boundary is defined by the first minimum of C60−solvent RDF) as a function of bulk solvent density and C60−solvent RDF gC60−solv(r) for different bulk densities of the solvent are shown in Figure 2B,C respectively. The average number of solvent particles ⟨N⟩ in the first solvation shell increases linearly at lower bulk densities until it reaches a plateau at higher bulk densities. The RDF plots show that the nature of the C60−solvent g(r) is almost the same for ρ0 up to 0.24 (see Figure 2C), and for bulk densities more than 1062
DOI: 10.1021/acsomega.7b01858 ACS Omega 2018, 3, 1060−1068
Article
ACS Omega
Figure 3. ADF for density of waterlike solvent particles in the entire first solvation shell of C60. The solvation shell is defined by the position of the first minimum (r* = 2.5) of the gC60−water(r). The black dots are the position of the carbon atoms of C60. (A−F) ADFs at T* = 0.25 for different bulk solvent densities (ρ0) as indicated in the legends. Inset in (D): enlarged view of the alluvial fan-shaped (triangular) high-density region corresponding to a six-membered carbon ring.
the value of ⟨N⟩ in Figure 2B at ρ0 = 0.4 and above is more than 32. At ρ0 = 0.4, all of the pentagonal sites are of the same solvent density (see Figure 3D) as those in case of ρ0 = 0.24. The solvent particles in excess of 32 are, therefore, have accumulated preferentially in the hexagonal face, but not distributed uniformly. The solvent particles are arranged in an alluvial fan-shaped (Y-shaped) pattern (see Figure 3D) on the hexagonal face. We note that this pattern is formed due to the growth of high-density vertices (cf. Figure 1E) toward the same in the adjacent six-membered rings (see Figure 3D). In the inset of Figure 3D, the enlarged view of the high-density region on top of one of the six-membered carbon rings is shown. It is consistent with the observed triangular shape of the low-energy region on top of a six-membered ring, as shown in the inset of Figure 2D. The additional particles at even higher bulk density of ρ0 = 0.6 are then accumulated along the line joining the centers of two adjacent six-membered rings and thus we observe the high-density lines joining centers of all of the hexagonal faces surrounding a pentagonal face (see red and purple contours in Figure 3E). At ρ0 = 0.8, the solvent particles are spread all over the surfaces except for the positions corresponding to carbon atoms and bonding regions of C60. To gain further insight, we dissect the first solvation shell for ρ0 = 0.4 into three layers (Figure S2 of Supporting Information (SI)) of equal width. The ordering observed in Figure 3D,E is a consequence of accumulation of additional particles in the subsequent layers within the first solvation shell (see Figure S2B,C). It is also important to note that the accumulation of solvent particles at different sublayers of the first solvation shell is complementary to each other (see Figure S2 and the corresponding discussion in the SI). It suggests that the azimuthal ordering in subsequent layers proceeds due to short-
0.4, a new peak develops at r* = 3.0. In Figure 2D, we show the energy map considering both C60−solvent and solvent−solvent interactions obtained from a configuration, in which each of the pentagonal and hexagonal faces (except the one on which we place a test solvent molecule) is occupied by a solvent particle (see Models and Methods). 2.2. Long-Range Azimuthal Ordering. It is interesting to examine whether the azimuthal ordering observed in the close vicinity of the solute surface (at a distance of r* = 2.09, i.e., 6.63 Å from the C60 center) is long-ranged and extends throughout the first solvation shell. In Figure 3A−F, the ADFs of solvent particles in the entire first solvation shell (i.e., considering all of the solvent particles in the first solvation shell) at different solvent densities are shown. A similar pattern to that observed in Figure 1A−C is observed for bulk densities ρ0 = 0.07, 0.16, or 0.24 (see Figure 3A−C), revealing that the same order as that observed in a thin layer surrounding the C60 molecule (at r* = 2.09) persists throughout the first solvation shell. However, the pattern changes for bulk densities of ρ0 = 0.4 and above (Figure 3D−F). Looking at the plots of average number of solvent particles in the first solvation shell ⟨N⟩ as a function of density (see Figure 2B), it is evident that at ρ0 = 0.24, there are almost 32 solvent molecules around the C60. Because there are 32 faces on a C60 surface, one can assume that the middle of all of these 32 faces is filled by the solvent particles. It is corroborated by the existence of almost equally bright spots on all of the 32 faces of the dodecahedron, as shown in Figure 3C at ρ0 = 0.24. However, the pattern changes significantly at and above the bulk density of 0.40. The change in the pattern at and above the bulk density of ρ0 = 0.24 is caused by the preferential accumulation of additional (in excess of 32) solvent particles in the first hydration shell. We see that 1063
DOI: 10.1021/acsomega.7b01858 ACS Omega 2018, 3, 1060−1068
Article
ACS Omega
Figure 4. Long-range ADF for density of waterlike solvent particles present in front of the second/third solvation shell of C60. Waterlike solvent densities are calculated for CG waterlike solvent in a spherical shell of 0.1 Å thickness around the C60. The middle of the layer is at a distance of 11.98 Å (r* = 3.78) from the center of the C60. (A−F) ADFs at T* = 0.25 for different bulk solvent densities (ρ0) as indicated in the legend.
2.3. Does Azimuthal Ordering Extend beyond First Solvation Shell? To examine whether this azimuthal ordering extends beyond the first hydration shell, we have analyzed the density map of a 0.1 Å slice of solvent layer in the beginning of the second/third hydration shell at around r* = 3.78, i.e., 11.98 Å (shown by second dotted line on the RDF of Figure 2C) and found (Figure 4A,B) that for ρ0 = 0.07 and 0.16 the azimuthal ordering persists and the ordering is complimentary (as the pentagonal faces are now denser than the hexagonal faces) to that in the first solvation shell (cf. Figure 1A,B). At ρ0 = 0.24 and 0.4, the long-range order ceases to manifest and a randomized solvent structure with almost homogeneous density over the entire surface is observed (see Figure 4C,D). It will be clear if we look at the RDF plots of Figure 2C. At bulk densities up to ρ0 = 0.16, the region between the first (at r* = 2.12) and second (at r* = 4.25) g(r) peaks is devoid of any solvent particles (see Figure S3 in SI) and thus the effect of the azimuthally ordered structure of the first solvation shell extends up to the second solvation shell, giving rise to a complimentary order at the second solvation shell. At ρ0 = 0.24 and 0.4, some solvent particles occupy the interlayer region, but with no specific ordering and due to these randomized solvent particles, the effect of the ordered first solvation shell does not extend up to the second solvation shell and therefore we cannot have azimuthally ordered second shell. However, at densities ρ0 = 0.6 and above, a new peak (at r* = 3) originates between the first peak and the one at r* = 4.0 in the RDF (cf. Figure 2C). Because of these radially ordered water molecules at r* = 3.0, a significant but different pattern is observed (see Figure 4E,F) in front of the second solvation shell (at a distance r* = 3.78 from the C60 center). The reemergence of the pattern at ρ0 = 0.6 and 0.8 is a result of short-range correlation arising due to the intervening solvation layer at r* = 3 at higher densities. It is
range correlation arising due to excluded volume effect. The triangular density pattern (Figures 1D,E and 3D,E) at the hexagonal faces suggests that this face is not energetically homogeneous; instead, a triangular low-energy region exists on this face. However, the energy plot in Figure 2A obtained by considering only C60−solvent interaction (neglecting solvent− solvent interaction) suggests that the middle of both the pentagonal and hexagonal faces is energetically homogeneous (as we have found all spots circular in Figure 2A). This is true for a low-density system, in which interparticle distances are large and thus the solvent−solvent interaction is negligible. That is why, at lower bulk densities, we observe circular bright spots on both the pentagonal and hexagonal faces, indicating that the solvent particles are homogeneously distributed on the central region of the pentagonal or hexagonal face, if occupied. For a high-density system, in which solvent−solvent correlation is significant, both C60−solvent and solvent−solvent interactions are to be considered and therefore we have calculated a modified energy map by considering only a single layer of 32 solvent particles (one particle is placed on the middle of each face) on the C60 surface. We remove one solvent molecule from a particular face and place a test solvent particle at different positions on that face (keeping other 31 solvent particles in the prefixed positions) and calculate the energy of the test solvent particle with the remaining 31 solvent particles as well as with all of the carbon atoms of the C60 molecule. We repeat the calculation for all of the faces. The two-dimensional energy map is shown in Figure 2D. Now, we observe (see Figure 2D) the triangular region of low energy on the hexagonal faces. Thus, it justifies the triangular pattern on the hexagonal faces in the ADFs at higher densities (cf. Figures 1D,E and 3D,E), where solvent−solvent interaction is important. 1064
DOI: 10.1021/acsomega.7b01858 ACS Omega 2018, 3, 1060−1068
Article
ACS Omega
pattern formation. The core-soft model not only reproduces properties of waterlike hydrogen-bonded (HB) liquids, but also represents a large class of non-HB materials with anomalous behavior, e.g., Te, Ga, S, Bi, Ge15Te85, graphite, silica, phosphorous, silicon, and BeF2. This pattern-formation ability at the nanoscale can therefore be utilized to generate patterned layers of any of these materials using a structurally patterned template like fullerene. The present study will pose a challenge to experimentally observe such patterns at the nanoscale. The nonpolarizable, classical−mechanical description adopted in the present study excludes any possibility of charge transfer as observed60 in a recent DFT study on the C60−water system.
important to note that the long-range order seen by Levitt and co-workers21 was a consequence of long-ranged interaction due to polarization, whereas the long-range azimuthal ordering observed in the present study is due to short-range correlation. 2.4. Role of C60−Water Attractive Dispersion Interaction. It is important to assess the role of small attractive C60−water attractive dispersion interaction, which is generally thought59 to be inconsequential in determining liquid structure, on the pattern formation. We performed a new simulation of the C60−solvent system, in which the interaction between the carbon atom of the C60 and the solvent particle is made purely repulsive (λ = 0 in eq 2; see Models and Methods and SI) and the density of the solvent is considered to be ρ0 = 0.07. The ADF of the solvent molecules in a thin layer of width 0.1 Å at a radial distance of r* = 2.09 from the C60 center is shown in Figure 5. No azimuthal ordering in the ADF in this case (see
4. MODELS AND METHODS Here, we consider an atomistic C60 immersed in a box of solvent molecules interacting with each other through a spherically symmetric two-lengthscale potential. We use a coarse-grained model for the solvent (water or any anomalous liquid), and the solvent−solvent interaction is represented by a spherically symmetric model potential with two lengthscales. This potential was originally designed long back by Hemmer and Stell41 by combining the Lennard-Jones potential with a Gaussian function, and the dimensionless form of this model potential form can be expressed as ⎡⎛ σ ⎞12 ⎛ σ ⎞6 ⎤ ⎡ 1 ⎛ r − r ⎞2 ⎤ 0⎟ ⎥ Uvv(r ) = 4⎢⎜ ⎟ − ⎜ ⎟ ⎥ + a exp⎢ − 2 ⎜ ⎝ ⎠⎦ ⎝ ⎠ ⎝ ⎠ σ r r c ⎣ ⎣ ⎦
Figure 5. Short-range ADF for local density of waterlike solvent particles present in front of the first solvation shell of C60 at a bulk density ρ0 = 0.07 when the interaction between the carbon atoms of C60 and the CG waterlike solvent particle is purely repulsive (λ = 0 in eq 2). No pattern (randomized waterlike solvent particles all over the C60 surface) observed in this case indicates the importance of C60− waterlike solvent particle van der Waals interaction in the pattern formation in the immediate vicinity of the solute.
(1)
where a and c, respectively, are dimensionless parameters associated with the height and width of the Gaussian part of the potential. The parameter σ in the above equation is the usual size parameter. We have considered a = 5.0, c = 1.0, and r0/σ = 0.7, which represents a two-lengthscale potential, consisting of a hard core and a soft corona (see labels in Figure S1a) represented by a repulsive shoulder (see Figure S1a of the Supporting Information (SI)). This potential of water represents most of the structural, dynamic, and thermodynamic anomalies of real water.33−42,44−46,48−52 The fullerene molecule is however described by atomistic model, in which a C60 molecule consists of 60 carbon atoms. Each carbon atom of the C60 molecule is modeled as an uncharged sp2 carbon atom with nonbonded potential parameters taken from the literature.53 The dimensionless interaction potential Uuv(r) between the carbon atom of C60 and the CG solvent particle is modeled by
Figure 5) has been observed. Thus, the pattern on the C60 surface is a consequence of solute−solvent attractive van der Waals interaction, and the same has been shown earlier.21
3. CONCLUSIONS In essence, using a simple, spherically symmetric, twolengthscale model for water, we are able to demonstrate beautiful pattern formation on the solvation layer of C60 and, more importantly, the pattern can be modified by changing the bulk density of the solvent. The effect of solute structure on the surrounding solvent ordering extends up to a distance of 11.98 Å from the C60 center. It is demonstrated that although the attractive solute−solvent van der Waals interaction is responsible for such pattern formation in the immediate vicinity of the solute, the change in the pattern at higher densities is a consequence of both solute−solvent and solvent− solvent interactions. Contrary to earlier study,21 in which the long-range effect of the surface structure of the C60 molecule was shown to be a consequence of the long-range interaction due to polarization, in the present investigation, as we have used nonpolar and nonpolarizable models, the long-range azimuthal ordering observed is due to short-range correlation among the adjacent radially ordered layers. Thus, the coresoftened model not only reproduces anomalous bulk properties33−42,44−46 and solvation thermodynamics,48 but also captures the formation of remarkable patterns on the solvation layers of a fullerene molecule and the long-range nature of the
⎡⎛ σ ⎞12 ⎛ σ ⎞6 ⎤ Uuv(r ) = 4⎢⎜ ⎟ − λ⎜ ⎟ ⎥ ⎝r⎠ ⎦ ⎣⎝ r ⎠
(2)
where the parameter λ tunes the attractive interaction between the atoms of the C60 with the CG waterlike solvent particle and thus dictates hydrophobicity of the C60 molecule. The solute− solvent potential Uuv(r) with λ = 0 represents a purely repulsive C60−water system with maximum hydrophobicity (see red curve in Figure S1b) and that with λ = 1 corresponds to full van der Waals dispersion interaction corresponding to force-field parameters for the carbon atom σC = 3.47 Å and εC = 0.2765 kJ/mol53,61 (see blue curve in Figure S1b). The interparticle potentials Uvv(r) and Uuv(r) with λ = 0 and 1 are shown in Figure S1a,b respectively. All of the quantities are made dimensionless by suitable scaling, e.g., distance r* = r/σvv, temperature T* = kBT/εvv, density ρ0 = ρσvv3, pressure P* = 1065
DOI: 10.1021/acsomega.7b01858 ACS Omega 2018, 3, 1060−1068
Article
ACS Omega
(
Pσvv3/εvv, and t * = t
m vv σvv
of carbon atoms of C60. The neglect of the interwater interaction in calculating the energy map is a good approximation for systems with lower bulk densities. However, at larger bulk densities, it is important to consider water−water interaction along with C60−water interaction when calculating the energy map. Therefore, we consider a representative configuration with water molecules occupying the centers of all of the faces of C60. Then, we remove a water molecule from one face, on which we place a test water molecule at different positions on that face and calculate the interaction energy of the test water molecule with all of the other 31 water molecules and 60 carbon atoms of C60. The same procedure is followed for all of the faces of C60. These energy values are then mapped onto the θ−ϕ plane along with the positions of the carbon atoms of C60. The energy map thus calculated is shown in Figure 2D.
1/2
εvv 2
)
where kB is the Boltzmann
constant and εvv and σvv are the energy and length parameters of the solvent, respectively. For scaling, we use water parameters σvv = 3.17 Å and εvv = 0.6502 kJ/mol corresponding to nonbonded parameters of water oxygen according to SPC/E model.62 The cross parameters are obtained using the Lorentz− Berthelot mixing rule. All of the simulations in this study were performed using an MD simulation code developed by us. A box of solvent molecules corresponding to a particular density ρ0 is simulated in NVE ensemble, and the temperature of the system is fixed by applying velocity scaling technique. During simulation, periodic boundary conditions and minimum image conventions are employed in all three directions. In each simulated system, the velocity scaling is performed in the first 50 000 steps and after that velocity scaling is stopped and additional 150 000 steps running in normal NVE simulation are discarded for equilibration. In each of the systems, a production run of 107 steps with the trajectories saved at an interval of 20 steps is used for final analyses. 4.1. Azimuthal Distribution Function (ADF). A threedimensional (3D) representation of solvent densities around a solute of any arbitrary shape can be obtained from the spatial distribution function (SDF) introduced by Kusalik and Svishchev.57,58 In this method, particles are accumulated in a three-dimensional bin using polar coordinates of solvent particles with respect to the solute and finally density in each bin is calculated with a proper normalization. However, to get a complete spatial arrangement of solvent molecules with respect to the arrangement of surface atoms of the C60, a twodimensional bin at a particular radial distance with azimuthal resolution is essential. Recently, Levitt and co-workers21 have introduced a modification of SDF known as azimuthal distribution function (ADF). They used the concept of Sanson−Flamsteed projection, also known as the Mercator equal-area projection method,63 which is generally used to draw geographical maps. For details, the reader is referred to ref 21 and its Supporting Information. Here, we have used the same technique to account for azimuthal ordering, if any, in the surrounding solvent of a buckminsterfullerene. For calculating the ADF for the front of the first solvation shell, a 0.1 Å spherical shell of solvent molecules situated at a distance of r* = 2.09 (6.625 Å) measured from the center of the C60 is considered. In case of ADF beyond the first solvation shell, a spherical shell of 0.1 Å thickness around a distance of r* = 3.78 (11.98 Å) is considered. Apart from these two cases, we have also considered the entire first solvation shell (defined by the first minimum in the C60−water RDF) for ADF calculation. For comparison among the ADFs at different densities, all azimuthal density values are scaled by the highest value of the same at a particular bulk density. The positions of the carbon atoms of the C60 are also mapped onto the same θ−ϕ plane. 4.2. Energy Maps. The energy map considering only solute−solvent interaction is shown in Figure 2A. It is obtained by placing a test solvent molecule everywhere in a spherical layer (over the C60 surface) at a distance of r* = 2.09 from the center of C60 and calculating the interaction of this solvent molecule with all of the atoms of the C60 (no other solvent molecules than the test molecule are considered). Finally, these energy values are mapped in the θ−ϕ (corresponding to latitude and longitude) plane in the same way as the density values are mapped in the case of ADF along with the positions
■
ASSOCIATED CONTENT
S Supporting Information *
The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acsomega.7b01858. Details of models and methods, results on the dissection of the solvation shell of C60, average number of particles in the region between the first and second solvation shells of the C60 (PDF)
■
AUTHOR INFORMATION
Corresponding Author
*E-mail:
[email protected],
[email protected],
[email protected]. ORCID
Niharendu Choudhury: 0000-0002-9959-8430 Present Address ∥
Department of Physics, UM-DAE Centre for Excellence in Basic Sciences, University of Mumbai Kalina Campus, Mumbai 400098, India (S.P.). Author Contributions
N.C. contributed in designing the research and writing the paper. N.C., D.B., and S.P. performed research, analyzed data, and prepared graphics. Notes
The authors declare no competing financial interest.
■
ACKNOWLEDGMENTS The authors acknowledge Computer Division, BARC, for providing ANUPAM supercomputing facility and support. They thank Drs. T.K. Ghanty, P.D. Naik, and S. Mohan for their kind interest and encouragement.
■
REFERENCES
(1) Dai, H.; Hafner, J. H.; Rinzler, A. G.; Colbert, D. T.; Smalley, R. E. Nanotubes as nanoprobes in scanning probe microscopy. Nature 1996, 384, 147−150. (2) Wong, S. S.; Harper, J. D.; Lansbury, P. T.; Lieber, C. M. Carbon nanotube tips: High resolution probes for imaging biological systems. J. Am. Chem. Soc. 1998, 120, 603−604. (3) Dresselhaus, M. S.; Dresselhaus, G.; Eklund, P. C. Science of Fullerenes and Carbon Nanotubes; Academic Press: New York, 1996. (4) Wilson, L. J.; et al. Metallofullerene drug design. Coord. Chem. Rev. 1999, 190−192, 199−207.
1066
DOI: 10.1021/acsomega.7b01858 ACS Omega 2018, 3, 1060−1068
Article
ACS Omega (5) Tagmatarchis, N.; Shinohara, H. Fullerenes in medicinal chemistry and their biological applications. Mini-Rev. Med. Chem. 2001, 1, 339−348. (6) Nakamura, E.; Isobe, H. Functionalized fullerenes in water. The first 10 years of their chemistry, biology and nanoscience. Acc. Chem. Res. 2003, 36, 807−814. (7) O̅ sawa, E., Ed.; Perspectives of Fullerene Nanotechnology; Kluwer: Dordrecht, The Netherlands, 2002. (8) Szymańska, I.; Radecka, H.; Radecki, J.; Kikut-Ligaj, D. Fullerene modified supported lipid membrane as sensitive element of sensor for odorants. Biosens. Bioelectron. 2001, 16, 911−915. (9) Da Ros, T.; Spalluto, G.; Prato, M. Biological applications of fullerene derivatives: A brief overview. Croat. Chem. Acta 2001, 74, 743−755. (10) Chandler, D. Interfaces and the driving force of hydrophobic assembly. Nature 2005, 437, 640−647. (11) Choudhury, N.; Pettitt, B. M. On the mechanism of hydrophobic association of nanoscopic solutes. J. Am. Chem. Soc. 2005, 127, 3556−3567. (12) Choudhury, N.; Pettitt, B. M. The dewetting transition and the hydrophobic effect. J. Am. Chem. Soc. 2007, 129, 4847−4852. (13) Choudhury, N. On the manifestation of hydrophobicity at the nanoscale. J. Phys. Chem. B 2008, 112, 6296−6300. (14) Huang, X.; Margulis, C. J.; Berne, B. J. Dewetting-induced collapse of hydrophobic particles. Proc. Natl. Acad. Sci. U.S.A. 2003, 100, 11953. (15) Lee, C.-Y.; McCammon, J. A.; Rossky, P. J. The structure of liquid water at an extended hydrophobic surface. J. Chem. Phys. 1984, 80, 4448−4455. (16) Howard, J. J.; Perkyns, J. S.; Choudhury, N.; Pettitt, B. M. Integral equation study of the hydrophobic interaction between graphene plates. J. Chem. Theory Comput. 2008, 4, 1928−1939. (17) Remsing, R. C.; et al. Pathways to dewetting in hydrophobic confinement. Proc. Natl. Acad. Sci. U.S.A. 2015, 112, 8181−8186. (18) Xi, E.; Remsing, R. C.; Patel, A. J. Sparse sampling of water density fluctuations in interfacial environments. J. Chem. Theory Comput. 2016, 12, 706−713. (19) Tiwary, P.; Mondal, J.; Morrone, J. A.; Berne, B. J. Role of water and steric constraints in the kinetics of cavity−ligand unbinding. Proc. Natl. Acad. Sci. U.S.A. 2015, 112, 12015−12019. (20) Weiss, D. R.; Raschke, T. M.; Levitt, M. How hydrophobic buckminsterfullerene affects surrounding water structure. J. Phys. Chem. B 2008, 112, 2981−2990. (21) Chopra, G.; Levitt, M. Remarkable patterns of surface water ordering around polarized buckminsterfullerene. Proc. Natl. Acad. Sci. U.S.A. 2011, 108, 14455−14460. (22) Bowman, C.; Newell, A. C. Natural patterns and wavelets. Rev. Mod. Phys. 1998, 70, 289−301. (23) Ball, P. The Self-Made Tapestry: Pattern Formation in Nature; Oxford University Press, 2001. (24) Möhwald, H. Direct characterization of monolayers in the air water interface. Thin Solid Films 1988, 159, 1−15. (25) Keller, S. L.; McConnell, H. M. Stripe phases in lipid monolayers near a miscibility critical point. Phys. Rev. Lett. 1999, 82, 1602−1605. (26) Maclennan, J.; Seul, M. Novel stripe textures in nonchiral liquidcrystal films. Phys. Rev. Lett. 1992, 69, 2082−2085. (27) Harrison, C.; et al. Mechanisms of ordering in stripe patterns. Science 2000, 290, 1558−1560. (28) Mosbach, K.; Ramström, O. The emerging technique of molecular imprinting and its future impact on biotechnology. Nat. Biotechnol. 1996, 14, 163−170. (29) Malescio, G.; Pellicane, G. Stripe phases from isotropic repulsive interactions. Nat. Mater. 2003, 2, 97−100. (30) Seul, M.; Andelman, D. Domain shapes and patterns: the phenomenology of modulated phases. Science 1995, 267, 476−483. (31) Head-Gordon, T.; Stillinger, F. H. An orientational perturbation theory for pure liquid water. J. Chem. Phys. 1993, 98, 3313−3327.
(32) Garde, S.; Ashbaugh, H. S. Temperature dependence of hydrophobic hydration and entropy convergence in an isotropic model of water. J. Chem. Phys. 2001, 115, 977−982. (33) Jagla, E. A. Core-softened potentials and the anomalous properties of water. J. Chem. Phys. 1999, 111, 8980−8986. (34) Jagla, E. A. Phase behavior of a system of particles with core collapse. Phys. Rev. E 1998, 58, 1478−1486. (35) Buldyrev, S. V.; Franzese, G.; Giovambattista, N.; Malescio, G.; Sadr-Lahijany, M. R.; Scala, A.; Skibinsky, A.; Stanley, H. E. Models for a liquid−liquid phase transition. Phys. A 2002, 304, 23−42. (36) Xu, L.; Buldyrev, S. V.; Angell, C. A.; Stanley, H. E. Thermodynamics and dynamics of the two-scale spherically symmetric Jagla ramp model of anomalous liquids. Phys. Rev. E 2006, 74, No. 031108. (37) Franzese, G.; Malescio, G.; Skibinsky, A.; Buldyrev, S. V.; Stanley, H. E. Generic mechanism for generating a liquid−liquid phase transition. Nature 2001, 409, 692−695. (38) Kumar, P.; Buldyrev, S. V.; Sciortino, F.; Zaccarelli, E.; Stanley, H. E. Static and dynamic anomalies in a repulsive spherical ramp liquid: Theory and simulation. Phys. Rev. E 2005, 72, No. 021501. (39) Yan, Z.; Buldyrev, S. V.; Giovambattista, N.; Stanley, H. E. Structural order for one-scale and two-scale potentials. Phys. Rev. Lett. 2005, 95, No. 130604. (40) Yan, Z.; Buldyrev, S. V.; Giovambattista, N.; Debenedetti, P. G.; Stanley, H. E. Family of tunable spherically symmetric potentials that span the range from hard spheres to waterlike behavior. Phys. Rev. E 2006, 73, No. 051204. (41) Hemmer, P. C.; Stell, G. Fluids with several phase transitions. Phys. Rev. Lett. 1970, 24, 1284−1287. (42) Sadr-Lahijany, M. R.; Scala, A.; Buldyrev, S. V.; Stanley, H. E. Liquid-state anomalies and the Stell-Hemmer core-softened potential. Phys. Rev. Lett. 1998, 81, 4895−4898. (43) Errington, J. R.; Debenedetti, P. G. Relationship between structural order and the anomalies of liquid water. Nature 2001, 409, 318−321. (44) de Oliveira, A. B.; Netz, P. A.; Colla, T.; Barbosa, M. C. Structural anomalies for a three dimensional isotropic core-softened potential. J. Chem. Phys. 2006, 125, No. 124503. (45) Barraz, N. M., Jr.; Salcedo, E.; Barbosa, M. C. Thermodynamic, dynamic, and structural anomalies for shoulderlike potentials. J. Chem. Phys. 2009, 131, No. 094504. (46) Pant, S.; Gera, T.; Choudhury, N. Effect of attractive interactions on the water-like anomalies of a core-softened model potential. J. Chem. Phys. 2013, 139, No. 244505. (47) Bolhuis, P.; Frenkel, D. Prediction of an expanded-to-condensed transition in colloidal crystals. Phys. Rev. Lett. 1994, 72, 2211−2214. (48) Buldyrev, S. V.; Kumar, P.; Debenedetti, P. G.; Rossky, P. J.; Stanley, H. E. Water-like solvation thermodynamics in a spherically symmetric solvent model with two characteristic lengths. Proc. Natl. Acad. Sci. U.S.A. 2007, 104, 20177−20182. (49) Bordin, J. R.; de Oliveira, A. B.; Diehl, A.; Barbosa, M. C. Diffusion enhancement in core-softened fluid confined in nanotubes. J. Chem. Phys. 2012, 137, No. 084504. (50) Krott, L. B.; Barbosa, M. C. Anomalies in a water-like model confined between plates. J. Chem. Phys. 2013, 138, No. 084505. (51) Patawane, S.; Pant, S.; Choudhury, N. Solvation of fullerene in a course grained water: A molecular dynamics simulation study. AIP Conf. Proc. 2015, 1665, No. 040024. (52) Pantawane, S.; Choudhury, N. Molecular dynamics simulation of water in and around carbon nanotubes: A coarse-grained description. AIP Conf. Proc. 2016, 1731, No. 040008. (53) Choudhury, N. A molecular dynamics simulation study of buckyballs in water: Atomistic versus coarse-grained models of C60. J. Chem. Phys. 2006, 125, No. 034502. (54) Choudhury, N. Dynamics of water in solvation shells and intersolute regions of C60: A molecular dynamics simulation study. J. Phys. Chem. C 2007, 111, 2565−2572. 1067
DOI: 10.1021/acsomega.7b01858 ACS Omega 2018, 3, 1060−1068
Article
ACS Omega (55) Choudhury, N. Dynamics of water in the hydration shells of C60: Molecular dynamics simulation using a coarse-grained model. J. Phys. Chem. B 2007, 111, 10474−10480. (56) Allen, M. P.; Tildesley, D. J. Computer Simulation of Liquids; Academic Press: London, 1986. (57) Svishchev, I. M.; Kusalik, P. G. Structure in liquid water: A study of spatial distribution functions. J. Chem. Phys. 1993, 99, 3049−3058. (58) Kusalik, P. G.; Svishchev, I. M. The spatial structure in liquid water. Science 1994, 265, 1219−1221. (59) Weeks, J. D.; Chandler, D.; Andersen, H. C. Role of repulsive forces in determining the equilibrium structure of simple liquids. J. Chem. Phys. 1971, 54, 5237−5247. (60) Choi, J. I.; Snow, S. D.; Kim, J.-H.; Jang, S. S. Interaction of C60 with water: First-principles modeling and environmental implications. Environ. Sci. Technol. 2015, 49, 1529−1536. (61) Girifalco, L. A. Molecular properties of fullerene in the gas and solid phases. J. Phys. Chem. 1992, 96, 858−861. (62) Berendsen, H. J. C.; Grigera, J. R.; Straatsma, T. P. The missing term in effective pair potentials. J. Phys. Chem. 1987, 91, 6269−6271. (63) Steers, J. A. An Introduction to the Study of Map Projections; University London Press: London, 1949.
1068
DOI: 10.1021/acsomega.7b01858 ACS Omega 2018, 3, 1060−1068
