Dynamic Solvation Shell and Solubility of C60 in Organic Solvents

Aug 1, 2014 - The notion of (static) solvation shells has recently proved fruitful in revealing key molecular factors that dictate the solubility and ...
0 downloads 9 Views 2MB Size
Article pubs.acs.org/JPCB

Dynamic Solvation Shell and Solubility of C60 in Organic Solvents Chun I Wang,† Chi C. Hua,*,† and Show A. Chen‡ †

Department of Chemical Engineering, National Chung Cheng University, Chia-Yi 62102, Taiwan, Republic of China Department of Chemical Engineering, National Tsing Hua University, Hsin-Chu 30013, Taiwan, Republic of China



S Supporting Information *

ABSTRACT: The notion of (static) solvation shells has recently proved fruitful in revealing key molecular factors that dictate the solubility and aggregation properties of fullerene species in polar or ionic solvent media. Using molecular dynamics schemes with carefully evaluated force fields, we have scrutinized both the static and the dynamic features of the solvation shells of single C60 particle for three nonpolar organic solvents (i.e., chloroform, toluene, and chlorobenzene) and a range of system temperatures (i.e., T = 250−330 K). The central findings have been that, while the static structures of the solvation shell remain, in general, insensitive to the effects of changing solvent type or system temperature, the dynamic behavior of solvent molecules within the shell exhibits prominent dependence on both factors. Detailed analyses led us to propose the notion of dynamically stable solvation shell, effectiveness of which can be characterized by a new physical parameter defined as the ratio of two fundamental time constants representing, respectively, the solvent relaxation (or residence) time within the first solvation shell and the characteristic time required for the fullerene particle to diffuse a distance comparable to the shell thickness. We show that, for the five (two from the literature) different solvent media and the range of system temperatures examined herein, this parameter bears a value around unity and, in particular, correlates intimately with known trends of solubility for C60 solutions. We also provide evidence revealing that, in addition to fullerene−solvent interactions, solvent−solvent interactions play an important role, too, in shaping the dynamic solvation shell, as implied by recent experimental trends.

1. INTRODUCTION Since the successful synthesis of fullerene (C60) in 1985,1 it has served as a promising material for versatile applications in photovoltaics, cosmetics, biochemistry, biomedicine, and fuel cells.2,3 Often, the success of these applications was found to correlate intimately with the properties of fullerene in solution state. Among the most imperative tasks is, conceivably, to overcome its generally poor solubility as well as to capture the attendant aggregation state. In the past three decades, solubility of C60 in various solvent media has been extensively scrutinized.2,3 In classification, aromatic, nonpolar solvents are typically regarded as relatively good solvents, compared with aliphatic solvents, such as chloroform, and polar protic solvents, such as ethanol and water;4 see, for example, the information provided in Table 1. Theoretically, through thermodynamic models and multiparameter statistical methods, the solubility of C60 has been systematically correlated with solvent’s topological characteristics (e.g., atom types, polarizability, and molecular connectivity), geometric features (e.g., radius of gyration and surface area), and electronic properties (e.g., atomic charge and molecular orbital).4−12 Alternatively, it has been common to observe that the fullerene solubility exhibits a maximum with varying system temperature,13−17 so far explained within the notion of “hypothetical solubility.”7,9,15,18 Although the major factors, as outlined above, affecting the solubility of C60 have been generally identified and well © 2014 American Chemical Society

Table 1. Solubility of C60 (s), Half-life Time (t1/2), Diffusivity of C60 (D), and the Corresponding Physical Parameter (ξ) at 300 K in Various Solvent Media CB T CF ethanolb H2Ob

sa [mg/mL]

t1/2 [ps]

6.45 2.40 0.28 1.2 × 10−3 1.3 × 10−11

261 ± 35 113 ± 17 55 ± 14 50 18

D [10−10 m2/s]

ξ

± ± ± ± ±

1.79 0.74 0.47 0.48 0.23

4.7 4.5 5.8 6.6 8.9

0.5 0.5 0.7 3.8 1.1

a

The solubility of C60 at 298 K as reported in ref 4. bThe half-life times of ethanol and water were obtained from refs 33 and 23, respectively.

explained within classical theories, in many cases the underlying physics remains not entirely transparent. In this respect, computer simulations have proven to be very useful in complementing essential molecular insights. Much computational work, for instance, has been devoted to understanding the behavior of hydrophobic hydration of C60 molecules in aqueous, biological systems using atomistic molecular dynamics (AMD) or Monte Carlo schemes. Specifically, a series of studies on the potential mean force (PMF) of C60 pairs have clearly revealed that the high atomistic surface density Received: July 2, 2014 Revised: July 28, 2014 Published: August 1, 2014 9964

dx.doi.org/10.1021/jp506572p | J. Phys. Chem. B 2014, 118, 9964−9973

The Journal of Physical Chemistry B

Article

associated with a single C60 molecule can bring about strong van der Waals (vdW) interaction with other C60 or solvent molecules.19−22 This, in turn, leads to low diffusivity and slow rotation mobility of water molecules and, consequently, large friction coefficients of C60.20,22−26 However, the presence of solvent (water) molecules may also induce repulsive interactions between two fullerene particles.27−29 In particular, we noticed that for water,21,23,25,26,30,31 ethanol,32,33 and ionic liquids,34−36 the erratic solvent molecules have been suggested to form shell-like structures encompassing a C60usually referred to as the solvation shells. As we elucidate in this work, a better understanding of the very nature of these solvation shells is key to deciphering some of the sophistications associated with fullerene solutions. To date, there have been few computer simulations dedicated to nonaqueous, organic solvent systems,26,37 nor has the nature of the solvation shells been unveiled. Banerjee, for example, had examined the effect of surface charge on the aggregation properties of C60 in three distinct solvents, including two organic solvents (i.e., acetone and toluene).26 There have also been coarse-grained models which offer perspectives of C60 aggregation in organic solvents.37−39 Recently, a significant proliferation of applications with organic photovoltaics imperatively demands a better understanding of the solution properties of fullerene in organic solvent systems. This observation, in part, motivates the present study aimed at resolving the solvation behavior of, as the simplest model, C60 in three organic solvents commonly utilized in the fabrication of photovoltaic devices. Through AMD simulations with carefully evaluated force fields, we have scrutinized single C60 particle suspended in aromatic (toluene and chlorobenzene) or aliphatic (chloroform) solvent medium for a range of system temperatures of practical interest. The (static and dynamic) analyses provided compelling evidence revealing that the properties of solvent molecules residing in the shell-like regions can substantially impact the solvation behavior of C60. For the first time, we elaborate on the notion of dynamically stable solvation shell, effectiveness of which can be well characterized by a new physical parameter defined as the ratio of two time constants representing, respectively, the solvent relaxation time within the first solvation shell and the characteristic time required for the fullerene particle to diffuse a distance comparable to the shell thickness. We show that, for the five different solvent media and a range of system temperatures examined herein, this parameter can be utilized to predict excellently the known trends of solubility for C60 solutions. This paper is organized as follows: in the following section, the force fields and computational methodology are described. Afterward, we discuss in detail both static and dynamic features of single-particle C60 solutions, and systematically correlate the central results with known salvation behavior of C60 whenever relevant data are available. Finally, a summary is provided with some future outlooks.

Figure 1. Chemical structures of fullerene (C60), toluene (T), chloroform (CF), and chlorobenzene (CB), the distribution of partial atomic charges on each solvent molecule, and the resulting dipole moments.

molecule. Preliminary investigations on the system density and self-diffusion coefficient (see the Supporting Information, SI) indicated that united-atom model works well for all molecules under investigation except for CF, for which full atomistic description was found to be necessary. For fullerene solutions, a series of force fields had previously been scrutinized by Maciel et al.42 Noticing that the DREIDING force field had not been assessed for this purposedespite that we had demonstrated its general plausibility for conducting polymer solutions43,44 its validity for C60 solution systems has been carefully inspected, as detailed in the SI as well as in prior work on a condensed system.45 In early AMD simulations of fullerene solutions, we noticed that partial atomic charges were typically considered only for polar solvents, such as water; for nonpolar solvents like those investigated herein, the effect was often ignored. To our surprise, the electrostatic interactions arising from partial atomic charges on an aromatic solvent species, chlorobenzene, have shown notable influences on the dynamic properties of solvation shells, although the corresponding static features seem to be little affected. So far, literature has accommodated various strategies to evaluate partial atomic charges.46 On the basis of a systematic comparison of the predicted dipole moments with the known ones (see the SI), the generalized atomic polar tensor (GAPT) scheme47 appears to describe rather well the solvent systems investigated. In this work, GAPT atomic charges were calculated using Gaussian 09 program,48 with the density functional theory (DFT) performed at the B3LYP/631G(d) level, and the results are illustrated in Figure 1. For comparative purposes, AMD simulations with or without incorporating partial atomic charges were compared at a system temperature of 300 K. For the rest of the work, partial atomic charges were included for a range of system temperatures: T = 250, 260, 270, 280, 290, 310, 320, 330 K. The general procedure of AMD simulations is as follows: A single C60 molecule was fixed at the center of a cubic box consisting of 1916 toluene, 1916 chlorobenzene, or 2185 chloroform. Each system was equilibrated using Nosé−Hoover NPT ensemble49 (with a coupling constant of 0.1 ps for both thermostat and barostat) at 1 atm and the target system temperature, for a duration of 1 ns and a time step of 1 fs. Periodic boundary conditions in all three directions were implemented. Longrange electrostatic interactions were calculated using the

2. SIMULATION PROTOCOL AMD simulations were performed using the software package DL_POLY_440 and the force field from DREIDING.41 Three solvent systems of C60 have been investigated, i.e., toluene (T), chlorobenzene (CB), and chloroform (CF). Figure 1 displays the chemical structures of C60 as well as the solvent molecules, with all hydrogen atoms lumped into the parent carbon atomsthe so-called united-atom modelexcept for the CF 9965

dx.doi.org/10.1021/jp506572p | J. Phys. Chem. B 2014, 118, 9964−9973

The Journal of Physical Chemistry B

Article

can be clearly identified at short distances, 6.0−9.7 Å, for both C60-T and C60-CB pairs, as the geometric center of phenyl ring of a solvent molecule was selected as the mapping center. In contrast, the commonly used mass center was adopted for the C60-CF system, where there is only one main peak in this region. In the former cases, the first peak corresponds to a faceon motif, and the second peak a side-on motif. Similar structural features had been observed in condensed MEH-PPV (poly(2-methoxy-5-(2′-ethylhexyloxy)-1,4-phenylenevinylene))/C60 hybrid system.45 As the mass center was adopted instead as the mapping center for the toluene molecule, we noticed good agreement with the results reported by Banerjee26 and Lambeth et al.31 Aside from the slight disparity between the two different types of solvents in this region, the RDFs in Figure 2 basically reveal no information about the distinct solvation behaviors observed experimentally for the three solvent systems considered. The same observations apply to the effect of varying system temperature or whether partial atomic charges were included (see Figures S5 and S6 of the SI). Still, the RDFs are very useful in revealing three major solvation “shells,” approximately in the ranges of 6.0−9.7, 9.7−14.6, and 14.6− 19.2 Å, as measured from the mass center of C60. Previously, similar shell regions had been identified with water21,23,25,26,30,31 and ethanol32,33 solutions of fullerene. The height of the peak corresponding to the first shell is usually very prominent among various C60 solutions, attributable to a high atomistic surface density of C60 and the resulting strong vdW interactions with the solvent molecules. Given that these common (static) features of C60 solutions seem not especially revealing of the actual solvation behavior (e.g., the RDFs of CB and T are nearly indistinguishable), we should pay particular attention to the corresponding dynamic features of, especially, the first solvation shell. 3.1.2. Orientation Profile. The orientation profile of solvent molecules is here referred to as the angular distribution with varying center-of-mass separation (r) from C60. The angle was measured by the vector defining the dipole moment of the solvent (shown in Figure 1) and that connecting the mass center of C60 and the mapping center of solvent. As can be seen in Figure 3, the orientation profile, unlike the RDFs discussed above, varies substantially among different solvent species. Similar to the case with the RDFs shown in Figure 2, however, the results are insensitive to the effect of partial atomic charges. For the two aromatic solvents, Figure 3a,b, the range of r < 12 Å clearly classifies some anisotropic orientation profile and, perhaps, certain network pattern, which differs markedly from what has been observed for a nonaromatic solvent, chloroform, shown in Figure 3c. It may be expected that the planar, aromatic solvents are especially inclined to forming anisotropic structures,52−56 because of the distinct vdW attractions that may be engendered by the two motifs mentioned above. In interesting contrast, the orientation profile of the nonaromatic solvent CF is dominated by a single molecular layer at r < 9 Å, where the dipole moments point inwardly to the radial axis of C60, as depicted in Figure 3c. The fluctuations of the orientation angle (vertical lines in Figure 3) are generally less than 30° in the vicinity of C60, i.e., r < 12 Å for aromatic solvents and r < 9 Å for CF, implying that solvent molecules in this region (which coincides with the first solvation shell) remain somewhat ordered, as the cartoons in Figure 3 illustrate. As solvent molecules move farther away from C60, a mean orientation angle of 90° is quickly resumed,

smoothed-particle-mesh-Ewald (SPME) technique with a real space cutoff of 12 Å. The same cutoff distance was utilized for the truncation of Lennard−Jones interactions, where standard Lorentz−Berthelot combining rules were adopted for unlike pairs of atoms. After the equilibration, statistical data were collected every 0.02 ps for another duration of 0.4 ns, which results in fairly convergent results for all subsequent analyses. Note that a relatively short time interval for data gathering was found to be necessary to yield accurate dynamic properties. Finally, the C60 was released to conduct free diffusion for a duration of 0.4 ns for evaluating the diffusivity of C60 in each solvent system. The results were obtained via time integration of the velocity autocorrelation function50 and listed in Table 1 (details shown in the SI). A few words about the effect of periodic boundary condition should be added. Yeh and Hummer reported that the predicted self-diffusivity of water molecules reduces with decreasing system size.51 In this work, the solvent self-diffusivities were also noted to systematically fall below the experimental values (see the SI), although the system sizes (70 × 70 × 70 Å3) are, in fact, considerably larger than in the work of Yeh and Hummer. In contrast, the (local) structural and relaxation properties of solvent molecules within a solvation shellthe main focuses of this studyare relatively insensitive to the system size, as might be expected.

3. RESULTS AND DISCUSSION In this section, we show that while the static analyses, including radial distribution function (RDF) and orientation profile, of the solvent molecules in single C60 solution system may be utilized to reveal the existence and geometrical feature of solvation shells, the dynamic features based on the occupation time correlation function (OTCF) and orientation autocorrelation function for the same shell regions are especially illuminating for interpreting the known trends of C60 solubility under varying solvent species or system temperature. 3.1. Static Features. 3.1.1. Radial Distribution Function. Figure 2 shows the RDFs for the C60-solvent pair in the solution system composed by T, CB, or CF. Two major peaks

Figure 2. Radial distribution functions (RDFs) for C60-T pair (blue solid line), C60-CB pair (red dotted line), and C60-CF pair (green dashed line). The black vertical lines mark three solvation shell regions, which are separated by the principal minima, for the later calculations of the occupation time correlation function of each solution system. 9966

dx.doi.org/10.1021/jp506572p | J. Phys. Chem. B 2014, 118, 9964−9973

The Journal of Physical Chemistry B

Article

Figure 3. Orientation profiles of solvent molecules surrounding C60 in (a) toluene, (b) chlorobenzene, and (c) chloroform. Red (blue) lines represent results for which the effect of partial atomic charges is included (excluded), whereas the vertical lines delineate the standard deviations from the mean orientation angle.

Figure 4. (a) Sketch of three solvation shells of C60; (b) the potential mean force of a C60 pair in vacuum (Tref = 300 K).45 The vertical lines in (b) mark the corresponding regions where two fullerenes may be kept separated in solution by the effect of the first or the second solvation shell.

Figure 5. Occupation time correlation functions of solvent molecules in the three solvation shells in (a) toluene, (b) chlorobenzene, and (c) chloroform solution.

and the wild fluctuations associated with it are indicative of a nearly free rotation, nonordered state of solvent molecules. Prior work had revealed analogous configurational ordering for C60 in ionic liquids (polar solvents)35 and for fulleride anions (C60n‑) in metal−ammonia solutions,34,57,58 albeit the individual effects of vdW and electrostatic interactions remain yet to be discriminated. Interestingly, a contrary phenomenon had been noted for C60 in water, where the originally welldefined molecular alignments seem to be greatly disrupted by the breaking of hydrogen bonds upon insertion of a fullerene molecule.59−65 In summary, as is evident from the literature and the present study, the (static) RDF and orientation profile analysis have proved very useful in identifying (three) shell-like regions and

attendant anisotropic structures of solvent molecules wherein. However, definitive solvent quality that is crucially relevant to determining the solvation power and aggregation properties of a fullerene solution remains faintly discernible in general relying only on these, more conventional, static analyses. Furthermore, we have so far noticed that the effect of partial atomic charges remains slim for the three nonpolar solvents investigated, yet the dynamic features to be discussed next seem to suggest otherwise. 3.2. Dynamic Properties. 3.2.1. General Analysis. The local dynamics of solvent molecules encompassing a single C60 can be revealed by the so-called occupation time correlation function (OTCF), which reflects the way a solvent molecule is escaping from a solvation shell where it initially resided. This 9967

dx.doi.org/10.1021/jp506572p | J. Phys. Chem. B 2014, 118, 9964−9973

The Journal of Physical Chemistry B

Article

analyzing scheme had previously been employed to investigate the solvation behavior of nanoparticles in aqueous or ethanol solution.23,24,33,66 A sketch illustrating the three solvation shells is provided in Figure 4a, where local minima of RDFs in Figure 2 have been utilized to define these shells. The OTCF was evaluated as follows: N

R (t ) =

⟨∑i = 1 θi(t0)θi(t + t0)⟩ N

⟨∑i = 1 θi(t0)θi(t0)⟩

(1)

where θi(t) bears a value of unity provided that the mass center of the ith solvent molecule remains in a specific shell at time t (t0 being the starting time), and is set to be zero otherwise; the angular brackets denote ensemble-plus-time average. Note that the above scheme does not preclude the contribution of solvent molecules that had previously escaped the shell and returned to it at some later time. Alternatively, the OTCF may be computed by excluding the contribution from these “doubly visited” molecules, by setting θi(t) to be zero once the molecule was detected to leave the shell for the first (passage) time.33 The last methodology is designated as “OTCF-fpt.” The first strategy is utilized in most of the following discussion, and results based on the second are compared in some cases so as to shed light on the dynamic stability of a solvation shell. As is evident in Figure 5, the relaxation pattern of OTCF for the first solvation shell differs markedly from those of the second and the third in all three solvent systems. Using the “half-life” time, t1/2, when half of the initial solvent population survives, i.e., OTCF = 0.5, as a basis for comparison,67 it is clear that solvent molecules on average stay substantially longer in the first solvation shell than in the other ones. For the sake of comparison and later discussion, a comprehensive gathering of these results is provided in Table 2. The reduced mobility of

Figure 6. Occupation time correlation functions of T, CB, and CF molecules in the first solvation shell, with (solid lines) or without (dashed lines) considering the effect of partial atomic charges in the AMD simulation; the inset show the corresponding results of OTCFfpt.

pair from being in close contact, as this is likely to result in energetically stable aggregation. Thus, Figure 4b shows the potential mean force of a C60 pair in vacuum, where a potential well and binding energy as large as ca. 18 kBTref can be realized once the pair is driven within a separation distance of about 10 Åthe innermost region of the first solvation shell in solution state. It can also been seen that as the separation distance is further enlarged to fall in the (equivalent) range of the second solvation shell, as the cartoons in Figure 4b illustrate, the interaction potential drops abruptly to below 2 k B T ref comparable to the thermal energy of a C60 particle. The central implication is, clearly, that solvent molecules in the first solvation shell may play a crucial role in preventing aggregation of C60 particles. Accordingly, it appears that the longer the solvent molecules stay within the first solvation shell, the more effective they could be to improve the solubility of the fullerene particles. A similar notion seemed to be embedded in recent attempts to create an effective solvation shell by adding charges to the fullerene particles,57,58 using ionic liquid,36,68 or mixing solvents.69 The correlation noted above between the half-lifetime of a solvent and its phenomenological solubility for C60, as listed in Table 1, seems to be well supported in general by noting, in passing, that a considerably smaller half-lifetime of ca. 20 ps in water23 compared with ca. 50 ps in ethanol33 for single C60 solutions corresponds consistently to a much poorer solubility of 1.3 × 10−11 mg/mL for the former than 1.2 × 10−3 mg/mL for the latter.4 The molecular factors practically affecting the (bulk) solubility of C60 can be numerous and quite complex, though, and thus might be difficult to be precisely differentiated one by one. The agreement noted presently, as well as in a later discussion concerning the effect of system temperature, should not be pure coincidence, however. A sound physical argument

Table 2. Occupation Half-life Time (t1/2) of Solvent Molecules Encompassing C60 solvation shell (unit: ps) first second third first (w/o)a first (fpt)a

T 113 47 35 86 30

± ± ± ± ±

CB 17 8 6 15 12

261 90 60 159 39

± ± ± ± ±

H2O23 ethanol33

CF 35 14 8 32 17

55 28 22 48 20

± ± ± ± ±

14 6 3 12 5

18 4 3

50 40 30 15

The abbreviation “fpt” represents the half-lifetime of OTCF-fpt; w/o denotes the cases where the effect of partial atomic charges was not accounted for in the simulation.

a

the solvent molecules within the first solvation shell can be attributed to the strong vdW interaction between C60 and solvent molecule, as Banerjee had pointed out for water, acetone and toluene solutions.26 For reasons to be explained shortly, it is especially instrumental to scrutinize features of this first solvation shell. Figure 6 compares the OTCFs of three solvent systems for the first solvation shell. It can be seen that CB possesses the longest half-lifetime (261 ps), followed by T (113 ps) and CF (55 ps). At first sight, the observed trend appears to correlate well with the relative solubility of C60 as listed in Table 1. Before the significance of this correlation would become even clearer, some additional background should be supplemented. A rational way to perceive the solvation power of a particular solvent species for C60 is to consider its ability to prevent a C60 9968

dx.doi.org/10.1021/jp506572p | J. Phys. Chem. B 2014, 118, 9964−9973

The Journal of Physical Chemistry B

Article

are seen to display the longest orientational relaxation in the first solvation shell, followed by T and CF molecules. This trend is in accord with what has been observed previously for the (purely diffusive) occupation time. Additional remarks on the effect of partial atomic charges can be found in SV of the SI. The above features, along with the static orientation profiles shown in Figure 3, clearly suggest that the two aromatic solvents are capable of forming patterned, and thus more stable, solvation shell in comparison with the nonaromatic solvent CF, which assumes a relatively simple (isotropic) layer structure. A better stabilized solvent shell, in turn, yields prolonged occupation time and, as we inferred earlier, a better solvation power for fullerene. A generally prominent impact of partial atomic charges, as seen in Figures 6 and 7 and Table 2, further suggested that, in addition to fullerene-solvent interactions, solvent−solvent interactions play an important role, too, in determining their solvation power for fullerenea central implication that, to our knowledge, has not captured widespread attention, yet it seems to be supported by recent experimental trends.4−11 3.2.3. Notion of Dynamically Stable Solvation Shell and a Fundamental Parameter. As introduced in a prior section, past research has been fruitful in making use of the renowned (static) features of solvation shells in explaining the highly sophisticated solvation and aggregation properties of fullerene solutions. Despite the fact that our dynamic analyses have clearly revealed the nonpermanence (which has been manifested by, for instance, a complete relaxation in occupation correlation time) of these static shells, their dynamic stability seems to play important role. Below, we show that comparison between OTCF and its counterpart OTCF-fpt sheds new light on how often the solvent molecules simply oscillate between adjacent solvation shells without actually drifting away, evocative of the notion of dynamically stable solvation shell. Recall that for the second scenario, referred to as OTCF-fpt, the solvent molecules will not be counted for the survival probability once they left a certain solvation shell for the first time. The results are illustrated in the inset of Figure 6, and the corresponding half-life times are listed in Table 2 for ease of comparison. In the first solvation shell, the relaxation of OTCFfpt for all solvent systems, including an ethanol solution reported in ref 33, proceeds much more rapidly than their counterpart OTCF. Notably, the half-life times of OTCF-fpt all fall in a range of 15−40 ps, compared with 50−250 ps for the case of OTCF, suggesting that solvent molecules spend a considerable time (about four to five times longer) simply oscillating between adjacent shells. Overall, the oscillatory solvent dynamics so unveiled was clearly suggestive of the existence of dynamically stable solvation shells, effectiveness of which can be anticipated to help stabilize fullerene particles in solution state. Given the dynamic attribute of these solvation shells, their effectiveness in preventing C60 aggregation awaits proper characterization. For this purpose, we propose a new physical parameter defined as the ratio of two fundamental time constants representing, respectively, the relaxation time of the first solvation shell and the characteristic time required for the fullerene particle to diffuse a distance comparable with the shell thickness. The physical significance of this parameter is to quantify the permanence of a solvation shell with respect to the attempting time that fullerene particles might approach one another and form into energetically stable aggregates. Clearly, the larger this parameter value is, the more effective the (first)

will be provided shortly that capitalizes on a more rigorous comparison between two fundamental time constants relevant to the ultimate determination of C60 solubility. 3.2.2. Stability of Solvation Shells. At this point, an intriguing question is to pursue what actually gives rise to the notable difference in the occupation time for various solvent systems. Although high surface atomic density of C60 would ubiquitously render strong vdW interactions between C60 and solvent molecules, especially in the first solvation shell, the detailed structural features of the solvent seem to differ substantially. In Figure 3, it was readily apparent that the two aromatic solvents are comparatively more structured than the aliphatic one in close vicinity of C60. The question is whether the dynamic features follow a similar trend as with the static structure above? More precisely, one may ask does a more structured static state guarantee its dynamic stability as well? The detailed, somewhat lengthy, analysis we shall present shortly seems to lend support to such assertion for C60 solutions. Namely, the capability of forming locally structured pattern for a particular solvent species helps sustain, through complex many-body interactions, the integrity of a solvation shell that consequently holds the solvent molecules within in a relatively prolonged period. Not only can this argument well explain the behavior of the three organic solvent species under investigation, but it also sheds light on the generally low solubility of aqueous solutions of C60 by noting that the collapse of hydrogen-bond network for water molecules adjacent to C60 would help lift the restraining force necessary to hold them any longer in a solvation shell.59−64 In short, it is a similar type of entropy-enthalpy coupling that is believed to contribute significantly to broadly varying the occupation time distributions in the solvation shell that, in turn, substantially diversifies fullerene solution properties. Below, we state that orientation relaxation of the solvent molecules can provide a meaningful measure of the structural stability (or endurance) of a solvation shell, through the following time autocorrelation function of the dipole moment: C(t ) = ⟨μ(t ) ·μ(t0)⟩

(2)

where μ is the unit vector of solvent dipole moment at the initial time (t0) or some later time (t), and the angular brackets denote ensemble-plus-time average. In Figure 7, CB molecules

Figure 7. Autocorrelation functions of orientation relaxation for three different solvent systems of C60 in the first solvation shell with (solid lines) or without (dashed lines) considering the effect of partial atomic charges in the AMD simulation. Results for the second and the third solvation shells follow a similar trend yet exhibit faster relaxations (see Figure S7 in the SI). 9969

dx.doi.org/10.1021/jp506572p | J. Phys. Chem. B 2014, 118, 9964−9973

The Journal of Physical Chemistry B

Article

solvation shell may act to impede any close contact of two fullerene particles. It is f rom this (relative) perspective that an effectively static solvation shell might be f ulf illed in reality. The physical parameter (ξ) is suggested to bear the following form: 2t1/2 ξ= 2 (3) a /D where a is the thickness of the first solvation shell (having in the present case a value of 3.7 Å), and t1/2 is the half-lifetime defined previously; D is the diffusivity of C60 in a specific solvent medium. Shown in Figure 8 and Table 1, this physical

Figure 8. Correlation between the solubility of C60 and the parameter ξ for T, CB, CF, ethanol, and water as the solvent medium; the linear regression so obtained is well described by ξ = 0.23s + 0.34, where s denotes C60 solubility. Figure 9. Predicted parameter values of ξ for (a) toluene, (b) chlorobenzene, and (c) chloroform solution of single C60 as a function of system temperature. The half-life times for chlorobenzene at 280, 270, 260, and 250 K all exceed 400 ps and, therefore, were beyond the reach of this AMD simulation.

parameter correlates intimately with the solubility of C60 for the five solvent media examined. The agreement is somewhat surprising, indeed, considering the significant disparity in solvent attributes. It is illuminating, therefore, to further consult the effect of system temperature for individual solvent species. 3.3. Temperature Effect. In Figure 9, we show the parameter ξ as functions of system temperature for the three solvent species investigated; the associated OTCFs can be found in the SI. For the range of system temperatures examined, the results displayed a profound maximum in the parameter of interest. We note, in particular, that the solubility of C60 in toluene had previously been reported for a wide range of system temperatures that may be compared with the present result. The experimental data were indicative of a maximal fullerene solubility at about 290 K,13,14,16,17 in good agreement with the present finding (see Figure 9a) assuming the correlation between ξ and C60 solubility. Note, however, that at a system temperature of 260 K, ξ was also found to peak with yet unverifiable sources; see discussion in SIII of the SI. Although no related data are available for the cases with CB and CF, it had been suggested that the incongruent melting temperature of C60 solid solvates (in a solvent medium) corresponds to the system temperature at maximal solubility. Accordingly, the maximal solubility of C60 in CB may be expected to occur at about 300 K,9 which is also consistent with the present finding (see Figure 9b).

Recent studies suggested that the solubility of C60 and similar molecules is an outcome of sophisticated compromise of system free energy,69−71 as also implied by this work. Clearly, much work will be needed to establish the correlations among a multitude of theoretical aspects that have so far been invoked to explain the solvation behavior of fullerene solutions. The present findings are, in this respect, aimed to shed some light on the molecular origins that are otherwise difficult to be perceived in the current pursuit of resolving the many puzzles associated with this marvelous material.

4. CONCLUSIONS This study revealed that, albeit no permanent static structures might be assigned to the recently discovered solvation shells, there appear to be dynamically stable ones whose properties correlate intimately with the known solvation behavior of C60 under varying solvent species or system temperature. In particular, we proposed a new physical parameter defined as the ratio of two fundamental time constants representing, respectively, the solvent molecule relaxation time within the first solvation shell and the characteristic time required for the fullerene particle to diffuse a distance comparable to the shell thickness. This parameter, bearing the physical significance of 9970

dx.doi.org/10.1021/jp506572p | J. Phys. Chem. B 2014, 118, 9964−9973

The Journal of Physical Chemistry B

Article

(5) Sivaraman, N.; Dhamodaran, R.; Kaliappan, I.; Srinivasan, T. G.; Rao, P. R. V.; Mathews, C. K. Solubility of C60 in Organic Solvents. J. Org. Chem. 1992, 57, 6077−6079. (6) Ruoff, R. S.; Tse, D. S.; Malhotra, R.; Lorents, D. C. Solubility of C60 in a Variety of Solvents. J. Phys. Chem. 1993, 97, 3379−3383. (7) Marcus, Y. Solubilities of Buckminsterfullerene and Sulfur Hexafluoride in Various Solvents. J. Phys. Chem. B 1997, 101, 8617− 8623. (8) Danauskas, S. M.; Jurs, P. C. Prediction of C60 Solubilities from Solvent Molecular Structures. J. Chem. Inf. Comput. Sci. 2001, 41, 419− 424. (9) Marcus, Y.; Smith, A. L.; Korobov, M. V.; Mirakyan, A. L.; Avramenko, N. V.; Stukalin, E. B. Solubility of C60 Fullerene. J. Phys. Chem. B 2001, 105, 2499−2506. (10) Sivaraman, N.; Srinivasan, T. G.; Vasudeva Rao, P. R.; Natarajan, R. QSPR Modeling for Solubility of Fullerene (C60) in Organic Solvents. J. Chem. Inf. Comput. Sci. 2001, 41, 1067−1074. (11) Stukalin, E. B.; Korobov, M. V.; Avramenko, N. V. Solvation Free Energies of the Fullerenes C60 and C70 in the Framework of Polarizable Continuum Model. J. Phys. Chem. B 2003, 107, 9692− 9700. (12) Van Lier, G.; De Vleeschouwer, F.; De Pril, P.; Geerlings, P. Theoretical Prediction of the Solubility of Fluorinated C60. Phys. Chem. Chem. Phys. 2009, 11, 5175−5179. (13) Ruoff, R. S.; Malhotra, R.; Huestis, D. L.; Tse, D. S.; Lorents, D. C. Anomalous Solubility Behaviour of C60. Nature 1993, 362, 140− 141. (14) Bezmelnitsin, V. N.; Eletskii, A. V.; Stepanov, E. V. Cluster Origin of Fullerene Solubility. J. Phys. Chem. 1994, 98, 6665−6667. (15) Smith, A. L.; Walter, E.; Korobov, M. V.; Gurvich, O. L. Some Enthalpies of Solution of C60 and C70. Thermodynamics of the Temperature Dependence of Fullerene Solubility. J. Phys. Chem. 1996, 100, 6775−6780. (16) Zhou, X.; Liu, J.; Jin, Z.; Gu, Z.; Wu, Y.; Sun, Y. Solubility of Fullerene C60 and C70 in Toluene, o-Xylene and Carbon Disulfide at Various Temperatures. Fullerene Sci. Technol. 1997, 5, 285−290. (17) Sawamura, S.; Fujita, N. High-Pressure Solubility of Fullerene C60 in Toluene. Carbon 2007, 45, 965−970. (18) Korobov, M. V.; Mirakyan, A. L.; Avramenko, N. V.; Olofsson, G.; Smith, A. L.; Ruoff, R. S. Calorimetric Studies of Solvates of C60 and C70 with Aromatic Solvents. J. Phys. Chem. B 1999, 103, 1339− 1346. (19) Walther, J. H.; Jaffe, R. L.; Kotsalis, E. M.; Werder, T.; Halicioglu, T.; Koumoutsakos, P. Hydrophobic Hydration of C60 and Carbon Nanotubes in Water. Carbon 2004, 42, 1185−1194. (20) Hotta, T.; Kimura, A.; Sasai, M. Fluctuating Hydration Structure around Nanometer-Size Hydrophobic Solutes. I. Caging and Drying around C60 and C60H60 Spheres. J. Phys. Chem. B 2005, 109, 18600− 18608. (21) Makowski, M.; Czaplewski, C.; Liwo, A.; Scheraga, H. A. Potential of Mean Force of Association of Large Hydrophobic Particles: Toward the Nanoscale Limit. J. Phys. Chem. B 2010, 114, 993−1003. (22) Morrone, J. A.; Li, J.; Berne, B. J. Interplay between Hydrodynamics and the Free Energy Surface in the Assembly of Nanoscale Hydrophobes. J. Phys. Chem. B 2012, 116, 378−389. (23) 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. (24) 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. (25) Fang, K.-C.; Weng, C.-I. Molecular Dynamics Simulations of Structural Features and Diffusion Properties of Fullerene-in-Water Suspensions. J. Colloid Interface Sci. 2008, 318, 188−194. (26) Banerjee, S. Molecular Dynamics Study of Self-Agglomeration of Charged Fullerenes in Solvents. J. Chem. Phys. 2013, 138, 044318.

the permanence of a solvation shell with respect to the attempting time of fullerene particles to form energetically stable aggregate, was shown to describe excellently known trends of solubility for five different C60 solutions and a range of system temperatures. In the future outlook, the universality of the discovered correlation awaits further testimony based on a wider variety of fullerene species and solutions. In addition to the usual fullerene-solvent interactions, we also found that solvent−solvent interactions play an essential role, too, in shaping the dynamic solvation shell. This fundamental observation, in agreement with recent experimental trends, helps rationalize the curious fact that chlorobenzene (CB) with its more pronounced polar attribute turns out to be a better solvent than toluene (T) for the nonpolar solute C60. Overall, it appears that understanding the nature of dynamic solvation shells could open up a new direction, among existing and yet unexplored ones, to perceive the highly sophisticated solvation behavior of fullerene species.



ASSOCIATED CONTENT

S Supporting Information *

Validation of force fields of the present AMD simulation; detailed comparison of the predicted partial atomic charges from four standard schemes; methodologies for computing solvent self-diffusivity and the C60 diffusivity in various solvent media; a complete gathering of RDFs and OTCFs not shown in the main text; additional remarks on the effect of partial atomic charges. This material is available free of charge via the Internet at http://pubs.acs.org.



AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS This work is sponsored by the National Science Council of ROC. Resources provided by the National Center for HighPerformance Computing of ROC are also acknowledged.



ABBREVIATIONS AMD, atomistic molecular dynamics; PMF, potential mean force; vdW, van der Waals; T, toluene; CB, chlorobenzene; CF, chloroform; GAPT, generalized atomic polar tensor; DFT, density functional theory; SPME, smoothed-particle-meshEwald; RDF, radial distribution function; OTCF, occupation time correlation function; OTCF-fpt, occupation time correlation function subject to the-first-passage-time criterion



REFERENCES

(1) Kroto, H. W.; Heath, J. R.; O’Brien, S. C.; Curl, R. F.; Smalley, R. E. C60: Buckminsterfullerene. Nature 1985, 318, 162−163. (2) McHedlov-Petrossyan, N. O. Fullerenes in Molecular Liquids. Solutions in “Good” Solvents: Another View. J. Mol. Liq. 2011, 161, 1−12. (3) McHedlov-Petrossyan, N. O. Fullerenes in Liquid Media: An Unsettling Intrusion into the Solution Chemistry. Chem. Rev. 2013, 113, 5149−5193. (4) Semenov, K. N.; Charykov, N. A.; Keskinov, V. A.; Piartman, A. K.; Blokhin, A. A.; Kopyrin, A. A. Solubility of Light Fullerenes in Organic Solvents. J. Chem. Eng. Data 2010, 55, 13−36. 9971

dx.doi.org/10.1021/jp506572p | J. Phys. Chem. B 2014, 118, 9964−9973

The Journal of Physical Chemistry B

Article

B.; Petersson, G. A.; et al. Gaussian 09, revision C.01; Gaussian, Inc.: Wallingford, CT, 2009. (49) Hoover, W. G. Constant-Pressure Equations of Motion. Phys. Rev. A 1986, 34, 2499−2500. (50) Allen, M. P.; Tildesley, D. J. Computer Simulation of Liquids; Oxford University Press: New York, USA, 1989. (51) Yeh, I.-C.; Hummer, G. System-Size Dependence of Diffusion Coefficients and Viscosities from Molecular Dynamics Simulations with Periodic Boundary Conditions. J. Phys. Chem. B 2004, 108, 15873−15879. (52) Fioroni, M.; Vogt, D. Toluene Model for Molecular Dynamics Simulations in the Ranges 298 < T (K) < 350 and 0.1 < P (MPa) < 10. J. Phys. Chem. B 2004, 108, 11774−11781. (53) Baker, C. M.; Grant, G. H. The Structure of Liquid Benzene. J. Chem. Theory Comput. 2006, 2, 947−955. (54) Headen, T. F.; Howard, C. A.; Skipper, N. T.; Wilkinson, M. A.; Bowron, D. T.; Soper, A. K. Structure of π−π Interactions in Aromatic Liquids. J. Am. Chem. Soc. 2010, 132, 5735−5742. (55) Gamieldien, M. R.; Strümpfer, J.; Naidoo, K. J. HydrationDetermined Orientational Preferences in Aromatic Association from Benzene Dimer Free Energy Volumes. J. Phys. Chem. B 2012, 116, 324−331. (56) Takeuchi, H. Structural Features of Small Benzene Clusters (C6H6)n (n ≤ 30) As Investigated with the All-Atom OPLS Potential. J. Phys. Chem. A 2012, 116, 10172−10181. (57) Howard, C. A.; Wasse, J. C.; Skipper, N. T.; Thompson, H.; Soper, A. K. The Solvation Structure of Fulleride C605− Anions in Potassium Ammonia Solution. J. Phys. Chem. C 2007, 111, 5640− 5647. (58) Howard, C. A.; Thompson, H.; Wasse, J. C.; Skipper, N. T. Formation of Giant Solvation Shells around Fulleride Anions in Liquid Ammonia. J. Am. Chem. Soc. 2004, 126, 13228−13229. (59) Scharff, P.; Risch, K.; Carta-Abelmann, L.; Dmytruk, I. M.; Bilyi, M. M.; Golub, O. A.; Khavryuchenko, A. V.; Buzaneva, E. V.; Aksenov, V. L.; Avdeev, M. V.; et al. Structure of C60 Fullerene in Water: Spectroscopic Data. Carbon 2004, 42, 1203−1206. (60) Hernández-Rojas, J.; Bretón, J.; Gomez Llorente, J. M.; Wales, D. J. Global Potential Energy Minima of C60(H2O)n Clusters. J. Phys. Chem. B 2006, 110, 13357−13362. (61) Weiss, D. R.; Raschke, T. M.; Levitt, M. How Hydrophobic Buckminsterfullerene Affects Surrounding Water Structure. J. Phys. Chem. B 2008, 112, 2981−2990. (62) Labille, J.; Masion, A.; Ziarelli, F.; Rose, J.; Brant, J.; Villiéras, F.; Pelletier, M.; Borschneck, D.; Wiesner, M. R.; Bottero, J.-Y. Hydration and Dispersion of C60 in Aqueous Systems: The Nature of Water− Fullerene Interactions. Langmuir 2009, 25, 11232−11235. (63) Maciel, C.; Fileti, E. E.; Rivelino, R. Assessing the Solvation Mechanism of C60(OH)24 in Aqueous Solution. Chem. Phys. Lett. 2011, 507, 244−247. (64) Colherinhas, G.; Fonseca, T. L.; Fileti, E. E. Theoretical Analysis of the Hydration of C60 in Normal and Supercritical Conditions. Carbon 2011, 49, 187−192. (65) Chopra, G.; Levitt, M. Remarkable Patterns of Surface Water Ordering around Polarized Buckminsterfullerene. Proc. Natl. Acad. Sci. U. S. A. 2011, 108, 14455−14460. (66) Choudhury, N.; Pettitt, B. M. Dynamics of Water Trapped between Hydrophobic Solutes. J. Phys. Chem. B 2005, 109, 6422− 6429. (67) In refs 23, 24, and 66, Choudhury showed that the relaxation time constants may be fitted from a stretched biexponential function for the OTCFs of an aqueous system. For the organic solvent systems investigated herein, we noticed that the biexponential model is unable to precisely capture our OTCF results. Therefore, the half-lifetime defined in the text was chosen as an alternative way to characterize the OTCF behavior. (68) Fileti, E. E.; Chaban, V. V. Structure and Supersaturation of Highly Concentrated Solutions of Buckyball in 1-Butyl-3-Methylimidazolium Tetrafluoroborate. J. Phys. Chem. B 2014, 118, 7376− 7382.

(27) Li, L.; Bedrov, D.; Smith, G. D. A Molecular-Dynamics Simulation Study of Solvent-Induced Repulsion between C60 Fullerenes in Water. J. Chem. Phys. 2005, 123, 204504. (28) Li, L.; Bedrov, D.; Smith, G. D. Repulsive Solvent-Induced Interaction between C60 Fullerenes in Water. Phys. Rev. E 2005, 71, 011502. (29) Kim, H.; Bedrov, D.; Smith, G. D. Molecular Dynamics Simulation Study of the Influence of Cluster Geometry on Formation of C60 Fullerene Clusters in Aqueous Solution. J. Chem. Theory Comput. 2008, 4, 335−340. (30) Rivelino, R.; Maniero, A. M.; Prudente, F. V.; Costa, L. S. Theoretical Calculations of the Structure and UV−vis Absorption Spectra of Hydrated C60 Fullerene. Carbon 2006, 44, 2925−2930. (31) Lambeth, B. P.; Junghans, C.; Kremer, K.; Clementi, C.; Site, L. D. Communication: On the Locality of Hydrogen Bond Networks at Hydrophobic Interfaces. J. Chem. Phys. 2010, 133, 221101. (32) Malaspina, T.; Fileti, E. E.; Rivelino, R. Structure and UV−vis Spectrum of C60 Fullerene in Ethanol: A Sequential Molecular Dynamics/Quantum Mechanics Study. J. Phys. Chem. B 2007, 111, 11935−11939. (33) Cao, Z.; Peng, Y.; Li, S.; Liu, L.; Yan, T. Molecular Dynamics Simulation of Fullerene C60 in Ethanol Solution. J. Phys. Chem. C 2009, 113, 3096−3104. (34) Howard, C. A.; Skipper, N. T. Computer Simulations of Fulleride Anions in Metal-Ammonia Solutions. J. Phys. Chem. B 2009, 113, 3324−3332. (35) Maciel, C.; Fileti, E. E. Molecular Interactions between Fullerene C60 and Ionic Liquids. Chem. Phys. Lett. 2013, 568−569, 75−79. (36) Chaban, V.; Maciel, C.; Fileti, E. Does the Like Dissolves Like Rule Hold for Fullerene and Ionic Liquids? J. Solution Chem. 2014, 43, 1019−1031. (37) Fritsch, S.; Junghans, C.; Kremer, K. Structure Formation of Toluene around C60: Implementation of the Adaptive Resolution Scheme (AdResS) into GROMACS. J. Chem. Theory Comput. 2012, 8, 398−403. (38) Chiu, C.-C.; DeVane, R.; Klein, M. L.; Shinoda, W.; Moore, P. B.; Nielsen, S. O. Coarse-Grained Potential Models for Phenyl-Based Molecules: II. Application to Fullerenes. J. Phys. Chem. B 2010, 114, 6394−6400. (39) Monticelli, L. On Atomistic and Coarse-Grained Models for C60 Fullerene. J. Chem. Theory Comput. 2012, 8, 1370−1378. (40) Todorov, I. T.; Smith, W.; Trachenko, K.; Dove, M. T. DL_POLY_3: New Dimensions in Molecular Dynamics Simulations via Massive Parallelism. J. Mater. Chem. 2006, 16, 1911−1918. (41) Mayo, S. L.; Olafson, B. D.; Goddard, W. A., III DREIDING: A Generic Force Field for Molecular Simulations. J. Phys. Chem. 1990, 94, 8897−8909. (42) Maciel, C.; Fileti, E. E.; Rivelino, R. Note on the Free Energy of Transfer of Fullerene C60 Simulated by Using Classical Potentials. J. Phys. Chem. B 2009, 113, 7045−7048. (43) Lee, C. K.; Hua, C. C.; Chen, S. A. Single-Chain and Aggregation Properties of Semiconducting Polymer Solutions Investigated by Coarse-Grained Langevin Dynamics Simulation. J. Phys. Chem. B 2008, 112, 11479−11489. (44) Lee, C. K.; Hua, C. C.; Chen, S. A. Multiscale Simulation for Conducting Conjugated Polymers from Solution to the Quenching State. J. Phys. Chem. B 2009, 113, 15937−15948. (45) Wang, C. I; Hsu, C. H.; Hua, C. C.; Chen, S. A. Molecular Dynamics Study of Pair Interactions, Interfacial Microstructure, and Nanomorphology of C60/MEH-PPV Hybrids. J. Polym. Res. 2013, 20, 188. (46) Cramer, C. J. Essentials of Computational Chemistry: Theories and Models; John Wiley & Sons Inc.: West Sussex, U.K., 2002. (47) Cioslowski, J. A New Population Analysis Based on Atomic Polar Tensors. J. Am. Chem. Soc. 1989, 111, 8333−8336. (48) Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Scuseria, G. E.; Robb, M. A.; Cheeseman, J. R.; Scalmani, G.; Barone, V.; Mennucci, 9972

dx.doi.org/10.1021/jp506572p | J. Phys. Chem. B 2014, 118, 9964−9973

The Journal of Physical Chemistry B

Article

(69) Chaban, V. V.; Maciel, C.; Fileti, E. E. Solvent Polarity Considerations Are Unable to Describe Fullerene Solvation Behavior. J. Phys. Chem. B 2014, 118, 3378−3384. (70) Herbst, M. H.; Dias, G. H. M.; Magalhães, J. G.; Tôrres, R. B.; Volpe, P. L. O. Enthalpy of Solution of Fullerene[60] in Some Aromatic Solvents. J. Mol. Liq. 2005, 118, 9−13. (71) Redmill, P. S.; Capps, S. L.; Cummings, P. T.; McCabe, C. A Molecular Dynamics Study of the Gibbs Free Energy of Solvation of Fullerene Particles in Octanol and Water. Carbon 2009, 47, 2865− 2874.

9973

dx.doi.org/10.1021/jp506572p | J. Phys. Chem. B 2014, 118, 9964−9973