Coarse-Grained Molecular Dynamics Simulations of the Sphere to

ARTICLE pubs.acs.org/Langmuir

Coarse-Grained Molecular Dynamics Simulations of the Sphere to Rod Transition in Surfactant Micelles Ashish V. Sangwai† and Radhakrishna Sureshkumar*,†,‡ †

Department of Biomedical and Chemical Engineering and ‡Department of Physics, Syracuse University, Syracuse, New York 13244, United States

bS Supporting Information ABSTRACT: Surfactant molecules self-assemble in aqueous solutions to form various micellar structures such as spheres, rods, or lamellae. Although phase transitions in surfactant solutions have been studied experimentally, their molecular mechanisms are still not well understood. In this work, we show that molecular dynamics (MD) simulations using the coarsegrained (CG) MARTINI force field and explicit CG solvent, validated against atomistic MD studies, can accurately represent micellar assemblies of cetyltrimethylammonium chloride (CTAC). The effect of salt on micellar structures is studied for aromatic anionic salts, e.g., sodium salicylate, and simple inorganic salts, e.g., sodium chloride. Above a threshold concentration, sodium salicylate induces a sphere to rod transition in the micelle. CG MD simulations are shown to capture the dynamics of this shape transition and support a mechanism based on the reduction in the micelle
water interfacial tension induced by the adsorption of the amphiphilic salicylate ions. At the threshold salt concentration, the interface is nearly saturated with adsorbed salicylate ions. Predictions of the effect of salt on the micelle structure in different CG solvent models, namely, single-site standard water and threesite polarizable water, show qualitative agreement. This suggests that phase transitions in aqueous micelle solutions could be investigated by using standard CG water models which allow for 3 orders of magnitude reduction in the computational time as compared to that required for atomistic MD simulations.

I. INTRODUCTION Surfactants are amphiphilic molecules with hydrophobic ‘tail’ and hydrophilic ‘head’ groups, which self-assemble in aqueous solutions into a wide variety of micellar structures such as spherical, rodlike, lamellar, vesicular, or crystalline phases.1
4 The structural phase of micelles depends on various parameters such as concentration and molecular structure of surfactant, counterions, salt, and solution temperature/pH. Certain micellar phases, for example cylindrical micelles at high surfactant and salt concentrations, exhibit rich rheological properties including viscoelasticity.5
11 Micelles find a wide range of applications in drug delivery,12
15 in turbulent friction drag reduction,16 as dispersant and detergents, and as templates for the synthesis of active nanomaterials.11,17
24 Applications of detergent micelles are found also in the nuclear magnetic resonance spectroscopy for structural determination of membrane proteins, where the difficulty to solubilize membrane proteins in polar solutions is overcome by their reconstitution into the micellar assemblies.25
29 Despite vast theoretical and experimental studies, understanding the molecular mechanisms of the formation of micelle structures of various shapes and how they are influenced by the chemical environment remains a challenge. In the past decade, many computational approaches have been utilized to study self-assembly in surfactant solutions. However, probing shape transitions in molecular detail requires relatively large systems consisting of tens of thousands of atoms to be sampled over several r 2011 American Chemical Society

hundreds of nanoseconds. This is an onerous task using atomistic molecular dynamics (MD) simulations even on modern parallel computing platforms.30
32 In this work, we carry out MD simulations of cetyltrimethylammonium chloride (CTAC) surfactant micelles modeled with coarse-grained (CG) potentials, validated against atomistic MD studies, and exploit the computational advantages of CG simulations to probe the dynamics of sphere to rod shape transitions. Several studies have reported that simple inorganic salts such as sodium chloride (NaCl) or sodium bromide (NaBr) induce sphere to rod transition in the micelles of various cationic surfactants above a certain threshold concentration.10,33
38 This transition is known to occur at much lower concentrations of salts with aromatic anions, e.g., sodium salicylate (NaSal).9,11,34,38
40 It is believed that anions such as Cl
or Br
are weakly bound to the cationic surfactant head groups constituting an electrical double layer around the micelle
water interface and effectively screen the intermicelle electrostatic repulsion. This results in a decrease in the surface area per surfactant and causes a sphere to rod transition of the micelle. Alternatively, it is proposed that aromatic anions such as Sal
penetrate into the micellar shell, increasing intramicellar ordering, which provides an Received: February 17, 2011 Revised: April 7, 2011 Published: April 27, 2011 6628

dx.doi.org/10.1021/la2006315 | Langmuir 2011, 27, 6628–6638

Langmuir entropic driving force for the shape transition.41,42 Many MD simulation studies have attempted to investigate the effects of different salts on the sphere to rod transition of micelles.30,41
46 Maillet et al. performed short (∼3 ns) atomistic MD simulations of preassembled spherical and cylindrical micelles of n-nonyltrimethylammonium chloride (C9 TAC) and erucyl bis[2-hydroxyethyl]methylammonium chloride (EMAC).42 The cylindrical micelles were shown to be stabilized with NaSal, while at high concentrations of NaCl the cylindrical micelle broke into smaller cylindrical fragments. In a notable work by Marrink et al., atomistic MD simulations of spontaneous formation of micelles of dodecylphosphocholine (DPC) surfactant molecules in the absence of salt were carried out.30 Formation of a long threadlike micelle was preferred through the periodic boundaries when the size of the surfactant
aggregate was comparable to the size of the simulation box. It was pointed out that a small simulation box induces artifacts in the dynamics of micellar aggregate formation due to spatial correlations or long hydrophobic interactions at the periodic boundaries. Piotrovskaya et al. analyzed the penetration of different additives such as acetone, 2-propanol, and sodium benzoate into the threadlike micelles of CTAC chloride to estimate their effects on the stability of structures.43 Wang and Larson used atomistic MD simulations to study the stability of an infinitely long threadlike micelle of CTAC in the presence of sodium salicylate and sodium chloride salts.45 The effect of different salts on the kinetics of shape fluctuations and transitions was estimated using relatively long (∼10 ns) simulations. While the above-mentioned studies render important molecular details to the mechanisms of the sphere to rod transition, a direct dynamical manifestation of the sphere to rod shape transition has not yet been captured in MD simulations. In a recent work, the structural properties of both spontaneously formed and preassembled micelles of nonionic dimethyldodecylamine-Noxide (DDAO) upon addition of ethyl butyrate-like oil molecules were studied.46 The oil molecules induced a near-spherical shape in DDAO micelles which were otherwise somewhat elongated. To overcome the limitations of time and length scales in atomistic MD, different coarse-grained schemes have been used to model surfactant self-assemblies.47
57 Lattice models with Monte Carlo simulations have been developed to study the phase behavior of surfactants and determine their critical micelle concentrations.55,57
59 Further, CG models of dimer surfactants represented by two types of particles (head or tail) with Lennard
Jones interactions60 have been utilized to predict the structural phase behavior and dynamical equilibrium of micelles.47,54 Implicit solvation has been widely used in both atomistic and CG systems to improve the sampling of micellar structures for longer time and length scales.56,61
67 Finitely extensible nonlinear elastic (FENE) potentials used with dissipative particle dynamics or even molecular dynamics have been explored for large-scale simulations of micellar self-assemblies to study rheological properties of wormlike micelles.48
50 While these studies are capable of making qualitative predictions of the phase behavior and rheological properties of micelles, they do not account for the effect of the chemical specificity of the intermolecular interactions among the surfactant, salt, and water molecules. As evidenced by experiments, the chemical environment plays a central role in determining the phase behavior and rheology of surfactant solutions. For instance, hydrophobic salts such as NaSal can not only greatly promote the sphere to rod transition but also abet formation of persistent or permanent flow-induced structures.68,69 The present study is an important first step toward realizing sufficiently long time scale CG MD simulations

ARTICLE

Table 1. Ratio of the Salt to CTAC Concentration in the 12 Different Simulations Employed in Set I, II, IIIa no salt

[NaSal]/[CTAC]

[NaCl]/[CTAC]

0

0.5, 1, 1.5, 2, 3, 5

1, 2, 3, 4, 10

Set I: Atomistic GROMOS-96 force field. Set II: MARTINI potentials with standard CG water. Set III: MARTINI potentials with polarizable CG water. a

of micellar fluids which allow for studying the effect of chemical environment on phase transitions. More recently, CG MARTINI models which reproduce the partitioning free energies between polar and apolar phases of a large number of chemical compounds have been developed.70,71 The advantage of the MARTINI force field is that it assigns interaction parameters to CG beads, in the varying degrees of polarity, charge, and van der Waals interactions, thereby preserving the chemical and physical characteristics of the molecule such as solubility/hydrophobicity as well as the dispersion interactions of the underlying atomistic structure. Coarse graining allows for the use of larger time steps of several 10 fs for the dynamics as compared to 1
2 fs in atomistic MD since the time scale of atomistic vibrations is neglected in the CG descriptions. Further, mapping of many (3
4) heavy atoms into a single bead significantly reduces the number of degrees of freedom in the simulations. In this work we used two explicit CG models for water, standard CG water and polarizable CG water.72,73 Both models lump 4 water molecules to a single CG entity. The standard CG water model is a single-site electroneutral Lennard
Jones particle in a constant dielectric field. Polarizable CG water is also electroneutral, although with a total of three sites, two of which have partial charges to maintain orientational polarizability. As a result, standard CG water effectively reduces a total of 12 atoms (4 water molecules) to 1 site, while polarizable CG water reduces it to 3. Overall, compared to atomistic MD, the computational time for CG MD is 3
4 orders of magnitude smaller. Hence, dynamics over time scales on the order of microseconds can be studied for systems with a linear dimension of O(10) nm. In order to validate the CG molecular models, we performed extensive comparisons of the structural features such as the distribution functions of surfactant, water, and ions in the micellar assemblies predicted by the CG and atomistic simulations. The atomistic and CG simulation systems are similar in terms of box dimensions, concentration of surfactant, molar ratios of the constituent chemical species, and ionic conditions. In section II, we describe the computational methods in more detail. Our findings are discussed in section III. Conclusions are offered in section IV.

II. COMPUTATIONAL METHODS Gromacs 4.0 molecular dynamics package74 was used to carry out a total of 36 simulations grouped into three different sets corresponding to different force fields. Specifically, sets I, II, and III correspond to the atomistic GROMOS-96 force field,75 MARTINI potentials in standard CG water,70 and MARTINI potentials with polarizable CG water73, respectively. Each set consists of 12 different simulations at varying molar ratios of NaSal or NaCl: see Table 1. CG MD simulation at a particular salt concentration was matched with the corresponding atomistic simulation to retain a cubic system size of approximately 9
10 nm on each side, keeping the number of surfactants and ions the same. A typical simulation system contained a preassembled micelle 6629

dx.doi.org/10.1021/la2006315 |Langmuir 2011, 27, 6628–6638

Langmuir

Figure 1. MARTINI coarse-graining (CG) scheme: (green) uncharged CG bead, (purple) CG bead with þ1 unit charge, and (red) CG bead with
1 unit charge. (a) Mapping of atomistic CTAC onto CG beads. C1: nonpolar bead with 4:1 mapping of C atoms. C2: nonpolar bead with 3:1 mapping of C atoms. Q0 = þ1 charged head group particle. (b) CG structure of CTAC with þ1 charge. (c) Mapping of atomistic Sal
anion onto CG beads. AC: aromatic carbon. SC4: slightly polar CG aromatic bead with 2:1 mapping. Qa:
1 charged particle CG bead with hydrogen bond acceptor capability. (d) CG structure of Sal
anion with
1 charge. Preassembled micelle structure: (e) atomistic and (f) CG. made up of 81 CTAþ surfactants with Cl
as counterions surrounded by positive and negative salt ions at a desired concentration in explicit water. This corresponds to 60 000
65 000 atoms with GROMOS-96 (set I), 8000
8300 particles with MARTINI in single-site standard CG water (set II), and 20 000
23 000 particles with MARTINI in three-site polarizable CG water (set III). Methods implemented in each set of simulations are as follows. II.A. Set I: Atomistic Simulations. GROMOS96 45a375 was used for modeling CTAþ, Sal
, Cl
, and Naþ, while water was modeled with SPC/E potential.76 The partial charges on CTAþ were taken from the literature on molecular dynamics of threadlike micelles of CTAC modeled with the same GROMOS96 force fields.45 The CH2 and CH3 groups on CTAþ as well as the aromatic CH on Sal
have been assigned the GROMOS96 united atom description, as shown in Figure 1a and 1c. The partial charges on the Sal
molecule were calculated by a geometry optimization with basis set 6-31G* and hybrid functional B3LYP77
80 using Gaussian 03.81 The optimized geometry was subjected to the restrained electrostatic potential (RESP) method to evaluate the atom-centered partial charges using RED-III.82,83 A constraint of zero charge was used on hydrogen atoms of aromatic CH on Sal
to obtain the partial charges for the united atom sites. The initial structures of the monomer CTAþ and Sal
were generated using the PRODRG server.84 The Gromacs topology files for atomistic CTAþ and Sal
are provided in the Supporting Information. For surfactants with long hydrocarbon chain such as C16 in CTAC, the equilibration times in the atomistic simulations for spontaneous micelle formation could be large. We therefore use a preassembled micelle structure which was generated using the PACKMOL program.85 The spherical structure was randomly generated from 81 surfactant molecules by constraining the head atom Nþ (Figure 1a) outside a sphere of a radius 1.9 nm and the other end C16 tail atom to the inside of a concentric sphere of radius 0.4 nm. The aggregation number 81 is close to the values observed in experiments.86
90 The initial structure generated by PACKMOL was solvated in roughly 22 000 pre-equilibrated SPC/E water molecules to

ARTICLE

make the box vector comparable to twice the diameter of the micelle (∼9.0 nm) and to avoid interactions among the periodic images of the micelle with itself. Counterions Cl
were added by replacing the overlapping water molecules. Subsequently, MD simulation was run for 5 ns in an isothermal
isobaric (NPT) ensemble to relax the micelle structure as well as to reach equilibrium density and energy. Systems in set I were subsequently prepared at different molar ratios of NaSal or NaCl (Table 1) by adding appropriate numbers of Naþ, Cl
, and/or Sal
into the equilibrated micelle box at energetically favorable coordinates. Equilibration at each molar ratio was performed for 3 ns, followed by a 15 ns production run. The simulations were maintained at 300 K and 1 atm in a NPT ensemble using Nose
Hoover91
94 temperature and Parrinello
Rahman isotropic pressure coupling. Newton’s equations of motion were integrated with a time step of 2 fs, and trajectory snapshots were collected at intervals of 0.4 ps. Periodic boundary conditions were applied with the minimum image convention.95 Bonds involving hydrogen were constrained with the LINCS algorithm.96 SPC/E water geometry was constrained using the SETTLE algorithm.97 Lennard
Jones interactions were terminated beyond 1.4 nm with the continuum correction for longer range contributions. Long-range electrostatic interactions were evaluated with the particle mesh Ewald (PME) method98 with a real space cutoff of 1.0 nm.

II.B. Sets II and III: Coarse-Grained (CG) Simulations. Simulations were performed with two water model schemes: set II standard CG70 and set III polarizable CG water.73 Standard CG water is a Lennard
Jones particle, mapped from 4 atomistic water molecules in a constant dielectric medium. This is an oversimplification of explicit water since it neglects the effect of explicit polarization in the solvent which can be significant in our simulation systems with charged particles or charged interfacial boundaries with the water phase. Further, the screening of intermolecular interactions by the solvent cannot be accurately represented. Alternatively, the polarizable CG water model accounts for the orientational polarizability of real water. To compare for the trade off between the two different CG water schemes in terms of the computational speed and accurateness of predictions, we performed MD simulations of standard and polarizable CG water models under identical conditions. Methodologically, the major difference between the two models is the incorporation of long-range electrostatic interactions in the latter. Such long-range interactions between charged surfactants and partially charged sites on polarizable water could become significant in determining the stability and dynamics of micelle structures. Consequently, PME was employed in the simulations of set III. The mapping of atomistic CTAþ into CG CTAþ is depicted in Figure 1a and 1b. The resultant coarse-grained CTAþ is a structure represented by 5 connected beads at a bond length of 0.47 nm. Detailed definitions of the nature of CG particle types are reported by Marrink et al.70 Particle types C1 and C2 represent 4 and 3 hydrocarbon groups, respectively. The trimethylammonium head group along with two adjoining hydrocarbon groups from the tail are represented by particle type Q0. Particle types C1 and C2 represent butane- and propane-type hydrocarbon analogs, respectively, while particle type Q0 has a unit positive charge with no hydrogen bond donating capability. The choice of CG particles for CTAC was based on previously modeled CG dipalmitoylphosphatidylcholine using the MARTINI force field.70 In Figure 1c and 1d, CG Sal
is represented by mapping the aromatic ring of Sal
onto an equilateral triangular MARTINI ring with a bond length of 0.27 nm made up of particle type SC4. Charged carboxyl group OdC
O
is represented by particle type Qa with unit negative charge which also possesses hydrogen bond acceptance capacity. Although CG Sal
could not explicitly incorporate the additional hydroxyl group at the ortho position in the aromatic ring of Sal
, the criterion of creating an aromatic anion with a slightly polar nature was satisfied. The Gromacs 6630

dx.doi.org/10.1021/la2006315 |Langmuir 2011, 27, 6628–6638

Langmuir

ARTICLE

Figure 2. Radial distribution function of surfactant sites from the center of mass (COM) of the micelle, gmic-surf(r), at different concentrations of NaSal for (a) atomistic and (b) CG simulations. Radial distribution function of water molecules from the center of mass of the micelle, gmic-wat(r), at different concentrations of NaSal for (c) atomistic and (d) CG simulations. Error bars are smaller than the symbol size. [NaSal]/[CTAC] molar ratio: 0.0 (0), 0.5 (4), 1.0 (), 1.5 (b), 2.0 (O). topology files for CG CTAþ and Sal
are provided in the Supporting Information. The preassembled spherical micelle structure of CG CTAþ was obtained from PACKMOL85 using 81 surfactants. The end group on the tail, C1, was constrained to an inner sphere of radius 0.5 nm and head group Q0 to a concentric outer sphere of radius 2.0 nm. The structure was solvated in a pre-equilibrated cubic box of CG water of roughly 10 nm on each side by removing overlapping water molecules and subsequently adding CG Cl
counterions at energetically favorable positions. Standard CG water has a freezing point in the range of 280
300 K. To avoid freezing, 10% of water particles were modeled with the antifreeze CG water provided in the MARTINI force field database. The structure was subjected to energy minimization and an equilibration run of 50 ns in a NPT ensemble to reach the correct density and energy. The density-equilibrated system was used to create initial structures at different concentrations (Table 1) by replacing water molecules with appropriate numbers of Sal
, Naþ, and Cl
ions. Systems were further allowed to relax for 10 ns, followed by production runs of 0.5 μs. Snapshots were collected every 7.5 ps. The time step for integrating Newton’s equations of motion was 30 fs. In the CG MD simulations, the Coulombic potential between charged groups was calculated as jij ¼

qi qj 4πε0 εr rij

ð1Þ

where qi and qj are the charges on two beads i and j, ε0 is the permittivity of the vacuum, and εr is the relative dielectric constant of the medium which is set to 15.0 for standard CG water. In polarizable CG water, εr = 2.5. For standard CG water, the Coulombic potential was truncated within a cutoff of 1.2 nm and long-range electrostatics forces were neglected, while PME summation with a real-space cutoff of 1.2 nm was used for polarizable water. Lennard
Jones interactions were reduced to zero gradually in the range of 0.9
1.2 nm using shifted potentials for both models. The temperature was maintained at 300 K with improved

Berendsen coupling,99 while pressure was maintained at 1 atm by employing the Parrinello
Rahman scheme.93,94

III. RESULTS AND DISCUSSION This section is divided into three subsections. In section III.A, a comparison of atomistic and CG micelle structure is provided at various concentrations of NaSal and NaCl salts. Section III.B describes the sphere to rod shape transition of the micelle. A comparison between micelle structures predicted by standard (set II) and polarizable (set III) CG water models is provided in section III.C. III.A. Micelle Structure Comparison: Atomistic vs CG MD Simulations. To compare the structure of the micelle obtained

from atomistic and CG simulations (Figure 1 e and 1f) using polarizable water, the radial distribution function gmic-surf(r) of surfactant atoms from the center of mass (COM) of the micelle is shown in Figure 2a and 2b, respectively. The distributions are generated from equilibrium sampling of 15 ns and 0.5 μs in the atomistic and CG simulations, respectively. Head group sites were also included in the calculation of gmic-surf(r) to estimate the size of the whole micelle. In the case of atomistic simulations at zero salt concentration, gmic-surf(r) has a value of about 10 at the COM of the micelle, which attains its maximum at r ≈ 1.6 nm and subsequently decreases monotonically to zero at r ≈ 3 nm. This implies that the micelle has a very well-defined hydrophobic core of radius 1.6 nm and the micelle
water interface extends up to a distance of 3 nm from the micelle COM. These features of the atomistic gmic-surf(r) are also exhibited by the CG structure despite the coarse graining of CTAþ in a simple 5 bead representation. The difference in peak heights in gmic-surf(r) originates from the evaluation of gmic-surf(r) with a different number of sites on CTAþ (i.e., atomistic, 20 sites; CG, 5 sites). Atomistic gmic-surf(r) from our study is in good agreement with 6631

dx.doi.org/10.1021/la2006315 |Langmuir 2011, 27, 6628–6638

Langmuir

ARTICLE

Figure 3. Radial distribution function g(r) of surfactant head groups from the center of mass (COM) of the micelle for CG (


) and atomistic simulations (). Radial distribution function g(r) of Cl
ions from CG (—) and atomistic simulations (O). (a) [NaSal]/[CTAC] = 0 and (b) [NaSal]/[CTAC] = 1.

Figure 4. Radial distribution function of surfactant sites from the center of mass (COM) of the micelle, gmic-surf(r), at different concentrations of NaCl for (a) atomistic and (b) CG simulations. Radial distribution function of surfactant head groups from the center of mass of micelle, gmic-head(r), at different concentrations of NaCl for (c) atomistic and (d) CG simulations. Error bars are smaller than the symbol size. [NaCl]/[CTAC] molar ratio: 0 (0), 1 (2), 2 (), 3 (b), 4 (O), 10 ((). [NaCl]/[CTAC] = 10 is shown only for the CG case.

the micelle COM
surfactant g(r) previously reported by Wang and Larson.45 Due to significantly enhanced sampling times, CG gmic-surf(r) profiles are smoother compared to their atomistic counterparts. As NaSal salt concentration is increased, gmic-surf(r) shows a decrease in the peak height. This effect is more evident in the CG micelles where the peak in gmic-surf(r) decreases in height while shifting toward smaller distances from the COM of the micelle. Simultaneously, the micelle
water interface broadens, indicating enhanced shape/size fluctuations, and extends past a distance of 3.5 nm from the micellar COM at [NaSal]/[CTAC] molar ratios of 1.5 and 2.0. This can be attributed to a sphere to rod transition of the micelle as discussed in section III.B The radial distribution functions of water from the COM of micelle gmic-wat(r) obtained from atomistic and CG simulations are shown in Figure 2c and 2d, respectively. At zero salt concentration,

the micelle completely expels water up to 1.3 nm from the COM beyond which gmic-wat(r) is nonzero, suggesting that the micelle
water interface starts at r ≈ 1.3 nm. As NaSal is added, water is gradually pushed out and the excluded volume of the whole micelle increases by ∼0.3 nm in radius. The increase in the size of micelle is due to the penetration and accommodation of aromatic anion Sal
into the micellar shell. The features of the CG gmic-wat(r) are very similar to that of the atomistic gmic-wat(r) except for the broadening of the interface observed at higher [NaSal]/[CTAC] ratios of 1.5 and 2.0 due to shape transitions. Figure 3a and 3b shows the distribution of surfactant head groups gmic-head(r) around the COM of micelle for atomistic and CG simulations at two different NaSal salt ratios. At zero NaSal concentration, the maximum in gmic-head(r) predicted by the two simulations is in good agreement, located at r ≈ 2.2 nm. 6632

dx.doi.org/10.1021/la2006315 |Langmuir 2011, 27, 6628–6638

Langmuir

ARTICLE

Figure 5. Radial distribution function of (a) surfactant head, (b) Cl
, (c) Sal
, and (d) Naþ ions from the center of mass (COM) of the micelle at different concentrations of NaSal in the solutions for CG system. [NaSal]/[CTAC] molar ratios are 0.0 (0), 0.5(4), 1.0 (), 1.5 (b), 2.0 (O). Insets in a indicate a change in the shape of the micelle corresponding to the trends in g(r), where surfactant head groups are purple beads and tail groups are green beads.

However, the profile of gmic-head(r) in the CG micelle at [NaSal]/[CTAC] = 1 is broader as compared to the atomistic one with the peak shifting by 0.2 nm toward a larger r. This is a result of shape fluctuations sampled by the longer CG simulations and the slight swelling of the micelle due to penetration of Sal
ions into the micelle corona. In the absence of salt, both atomistic and CG display magnitudes in gmic-Cl-(r) peaked to a value to 8 near the micelle
water interface. The maximum in gmic-Cl-(r) for the CG micelle is located at a r = 2.5 nm, which is a shift of 0.3 nm to the right. This can be attributed to the larger size of CG Cl
ions due to inclusion of coordinated water molecules in determining its van der Waals radius. At [NaSal]/ [CTAC] = 1, the peak in gmic-Cl-(r) disappears for atomistic and CG micelles, suggesting that Cl
ions dissociate from the interface into the bulk water as they are replaced by Sal
ions. Hydrophobic salts such as NaSal have a pronounced effect on the phase behavior of micellar solutions as compared to that of simple inorganic salts such as NaCl.35,36 We investigated whether CG potentials are capable of differentiating the effects of NaSal on the micelle structure from that of NaCl. As shown in Figure 4, radial distribution functions for the CG and atomistic micelles are practically identical for all cases studied. Further, no shape transitions were observed in the presence of NaCl even at salt concentrations as high as [NaCl]/[CTAC] = 10. The maximum of gmic-head(r) for molar ratio 10 was in fact somewhat elevated compared to the cases of lower NaCl molar ratios, implying that the structure is more spherical and head groups are more ordered. Note that the Debye length κ
1 estimates change significantly from 7.6 Å in the absence of added salt to 1.7 Å at [NaCl]/[CTAC] = 10. At high concentrations of NaCl, the lowering of κ
1 is signified by an electrical double layer of weakly bound Naþ and Cl
ions, which is expected to screen the electrostatic repulsion among head groups to promote shape

transitions. It is plausible that we are either in a regime below the threshold concentration of salt or surfactant or both to observe structural changes induced solely by electrostatic screening by simple inorganic salts. Overall, the above results establish that the response of CG micelles to changes in the chemical environment is quantitatively similar to that of atomistic ones. III.B. Effect of NaSal Salt: Sphere to Rod Transition. The distribution of surfactant head groups at the micelle
water interface showed peculiar variations with increasing NaSal molar ratios. At lower ratios ([NaSal]/[CTAC] e 1), g mic-head (r) was a relatively narrow Gaussian centered at ∼2.4 nm, indicating a well-defined spherically shaped micelle (Figure 5a). As we increased [NaSal]/[CTAC] from 0.0 to 1.0, there was a minor shift of the Gaussian peak to larger r values, suggesting slight swelling of the micelle. The initial spherical shape was preserved after long times (500 ns) even at a moderate salt ratio of [NaSal]/[CTAC] = 1 (see Figure 6a for micelle image). However, at concentrations of [NaSal]/[CTAC] > 1, the gmic-head(r) profile became skewed with lower maxima closer to the COM of the micelle. Further, the tail of the distribution extended beyond r = 3.7
3.9 nm. These features can be attributed to a sphere to rod shape transition. The shape transition to cylindrical structure at [NaSal]/[CTAC] = 2 is visually shown in Figure 6b and 6c, where an initially spherical micelle (at ∼6 ns) transformed into a cylindrical structure after 60 ns of simulation, respectively. The effect of NaSal on shape is also manifested in the distribution of ionic species Cl
, Sal
, and Naþ in solution. For example, the maximum in gmic-Cl-(r) gradually disappears with increasing ratios of added NaSal (Figure 5b). At [NaSal]/[CTAC] > 1, Cl
ions were completely dissociated from the micelle
water interface and are solvated in bulk water. This effect was also confirmed in atomistic simulations (Figure 3). The effect on gmic-Sal-(r) with increasing NaSal was quite similar to the trends 6633

dx.doi.org/10.1021/la2006315 |Langmuir 2011, 27, 6628–6638

Langmuir seen in gmic-head(r); see Figure 5c. The gmic-Sal-(r) distribution at [NaSal]/[CTAC] e 1 was broader and slightly shifted to larger r compared to the case of zero salt. On the other hand, features of cylindrical shape were evident in gmic-Sal-(r) at ratios of [NaSal]/[CTAC] > 1. The qualitative similarity in the trends of gmic-head(r) and gmic-Sal-(r) also suggests that the Sal
ions were significantly bound to the surfactant head groups. An increase in

Figure 6. Snapshots of micelle shape in the CG system (a) after 500 ns in [NaSal]/[CTAC] = 1 (spherical), (b) after 6 ns in [NaSal]/[CTAC] = 2 (spherical), and (c) after 60 ns in [NaSal]/[CTAC] = 2 (cylindrical). Surfactant head groups are purple beads and tail groups green beads.

ARTICLE

Sal
binding increased the number of Naþ ions distributed closer to the micelle (Figure 5d). Although this effect was not very pronounced, gmic-Naþ(r) gradually shifted closer to the COM of micelle with added NaSal. In order to investigate the dynamics of the spherical to rodlike shape transition, the variation in the radius of gyration rg of the CG micelle as a function of time is shown in Figure 7a. Only surfactant head groups were used in the calculation of rg, and the molar ratio of [NaSal]/[CTAC] was 2. A dynamic visualization of this transition is provided in a video file in the Supporting Information. Prior to collecting data for the calculation of rg, the system was simulated for about 10 ns to reach equilibrium density and energy, during which it acquired a slightly elliptical shape. After about 12 ns in the production run, rg started increasing from ∼2.45 to ∼2.75 nm during a span of 50 ns. The fluctuations in rg even after completely transitioning to the rodlike structure (t > 100 ns) were considerably large. The aggregation number of cylindrical micelles is larger than that in spherical ones. However, the preassembled structures used as initial conditions contain a fixed number of surfactants, implying

Figure 7. (a) Radius of gyration (nm) of CG micellar heads as a function of simulation time (ns) in production run for the case of [NaSal]/[CTAC] = 2. (b) Maximum radial distance (—) and minimum radial distance (O) of the CG surfactant head from the center of mass of the micelle as a function of simulation time.

Figure 8. (a) Radial distribution functions of surfactant beads (black line), water (blue dashed line), Sal
aromatic group (green line), and Sal
charge group (red dashed line) to show the distribution and orientation of Sal
at the micelle
water interface. (b) Number of bound Sal
ions to the micelle corona for spherical (b) and cylindrical (O) shapes and total number of Sal
in the system (). Snapshots showing the distribution and orientation of Sal
ions in a cylindrical micelle: (c) cross-sectional view and (d) lateral view in a system of molar ratio [NaSal]/[CTAC] = 2: aromatic group (cyan color) and charged carboxyl group (red). 6634

dx.doi.org/10.1021/la2006315 |Langmuir 2011, 27, 6628–6638

Langmuir

ARTICLE

Figure 9. Radial distribution function evaluated from the center of mass (COM) of the micelle for (a) surfactant beads, (b) water, (c) surfactant head group, and (d) Cl
ion. CG system in polarizable (set III) and standard water (set II) is denoted by thick and thin lines, respectively, for [NaSal]/[CTAC] = 0. CG system in polarizable (set III) and standard water (set II) is denoted by short thick dash and long thin dash, respectively, for [NaSal]/[CTAC] = 2.

a relatively undersaturated cylindrical micelle. In some snapshots, the micelle was seen in a transition state structure which was neither a well-defined rod nor a sphere. As an additional measure for quantifying the shape dynamics, the distance between the farthest and the closest head group from the COM of the micelle is plotted in Figure 7b. The longest radial distance from the center of a perfect hemispherically capped rod would be one-half the length of its major axis, while the shortest would be one-half the length of the minor axis. In a micelle transitioning to a rod from a sphere, the two distances (the farthest head group and the closest one from the COM of micelle) should increase with time. There was an initial difference of ∼1.3 nm among two extreme distances, implying an imperfect spherical structure, which increased to ∼3 nm after 70
100 ns of simulation. During this period, the aspect ratio of the cylindrical micelle ranged from 2.5 to 4. The partitioning of Sal
ions near the micelle
water interface is provided in Figure 8a by plotting distributions of surfactant, water, charged acid group on Sal
, and hydrophobic aromatic ring from the COM of micelle ([NaSal]/[CTAC] = 1). The distribution of water is plotted on a differently scaled right-hand Y axis to produce a comparison in g(r) of the surfactant and water. The distributions showed a distinct micelle
water interface with Sal
ions mainly distributed at the interface. From the g(r) of Sal
, the surface excess concentration of Sal
at this interface can be seen to be fairly high. In addition, the difference in the peak locations of the hydrophobic aromatic ring and charged group was ∼0.4 nm with the peak of the aromatic ring g(r) closer to the micelle core. This suggests that the aromatic rings were penetrating into the micelle core and the charge group remained oriented radially outward into the water phase, as evidenced by Figure 8c and 8d. The above-mentioned behavior of Sal
is very similar to that in surfactants at the oil
water interface in emulsions where the hydrophobic tails point into the

oil phase and the heads are situated in the aqueous phase. The degree of adsorption of Sal
on the micelle is dependent on [NaSal]/[CTAC] molar ratios. Average Sal
ions bound to the micelle surface as well as the total Sal
ions in the system are plotted in Figure 8b against the molar ratio of NaSal. At lower molar ratios of 0.5 and 1.0 of [NaSal]/[CTAC], almost all the Sal
ions (>95%) were bound to the micelle corona. As the salt ratio was increased beyond 1, a significant number of Sal
ions dissociated into the bulk of water, indicating that the saturation of Sal
binding to the micelle interface occurred in the range of [NaSal]/[CTAC] of 1
1.5. The lowest [NaSal]/[CTAC] ratio at which the shape transition from a sphere to a rod occurred was 1.5, which suggests that saturation in micelle-bound Sal
is a necessary condition to induce shape transition. The mechanism by which Sal
stimulates the shape transitions is reduction of interfacial tension at the micelle
water interface due to the surface active amphiphilic nature of Sal
. Our findings of Sal
orientation are consistent with the predictions of atomistic simulations. However, this is the first dynamical simulation of a sphere to a rod shape transition obtained using CG models which capture the influence of solvation and chemical environment at the molecular level. Previously, atomistic MD and Monte Carlo studies combined with free energy models of amphiphile aggregation showed that CTAþ cylindrical micelles without Sal
have a more positive head group free energy than spherical CTAþ micelles.41,45 Sal
shields repulsive interactions among the heads to cause a reduction in the free energy. With appropriate ordering of surfactants and Sal
, a sharp decrease of free energy of the surfactant head groups was observed in cylindrical micelles. In addition, the normal vector to the planar aromatic π rings of Sal
was found to be parallel to the axis of cylindrical micelles. Signatures of such ordering were not manifested during the shape transition predicted by the CG MD simulations. 6635

dx.doi.org/10.1021/la2006315 |Langmuir 2011, 27, 6628–6638

Langmuir III.C. MARTINI Solvents: Polarizable versus Standard CG Water. As discussed in the previous two sections, CG MARTINI

potentials implemented with polarizable CG water were capable of reproducing the structure of atomistic spherical micelle and predicted the shape changes as a function of NaSal salt. However, polarizable CG water is not computationally as advantageous as standard CG water, given the three-site nature of polarizable versus single-site standard CG water model. In addition, longrange electrostatic interactions present in polarizable CG water significantly reduced the simulation speed. On the other hand, system sizes are considerably small with standard CG water, i.e., system sizes in set III were roughly 20 000 CG particles, while set II was typically comprised of approximately 8000 solvent particles. In addition, calculation of long-range electrostatics was omitted to increase simulation speeds by roughly 15 times in the case of standard CG water. The distributions of surfactant from the COM of micelle gmicsurf(r) obtained using polarizable and standard CG water models are compared at two different ionic conditions of NaSal in Figure 9. In the absence of NaSal, the profiles of gmic-surf(r) remained qualitatively similar although the peak value of gmicsurf(r) was slightly lower in polarizable CG water compared to that of the standard CG model (Figure 9a). After addition of NaSal and subsequent shape transition, gmic-surf(r) flattened out for both models with their peaks shifting closer to the COM of the micelle. Similarly, as seen from gmic-wat(r) profiles, the micellar core is well reproduced at no salt condition (Figure 9b). The micellar core expands moderately after addition of NaSal with polarizable CG water accompanied by a shift in gmic-wat (r) toward larger r values. This effect, however, is more pronounced with standard CG water where the gmic-wat(r) is broader, suggesting a larger propensity to the cylindrical shape. These observations confirm that although the shape fluctuations in standard CG water are relatively larger in the presence of salt, the structure compares qualitatively well with that predicted by the polarizable CG water. Further, predictions of the profiles of both gmic-head(r) and gmic-Cl-(r) are also in qualitative agreement (Figure 9c and 9d). Hence, the computationally advantageous standard CG water model could be used to explore boundaries of shape/phase transitions. However, properties with strong dependence on electrostatics and solvent polarizability could get affected by such oversimplification of explicit solvent. Hence, more expensive polarizable CG water simulations may be necessary to fine tune the qualitative predictions of the simpler model.

IV. CONCLUSIONS The MARTINI CG model offers orders of magnitude speedup in MD simulations of aqueous micellar solutions as compared to those based on atomistic descriptions. The CG model is capable of predicting the structural features of spherical micelles of CTAC accurately. Specifically, the predictions for the distributions of different ionic species such as counterions, salt ions, and head groups and their variations with respect to aromatic and simple salts are quantitatively consistent with those obtained from atomistic MD simulations. CG simulations enable dynamical explorations using larger systems over time scales sufficiently long enough to capture shape fluctuations and shape transition of the micelle, achieved by addition of NaSal. Further, they allow for probing the molecular mechanisms driving such transitions. The simulations show that the amphiphilic nature of the Sal
ions result in their partitioning at the micelle
water

ARTICLE

interface in such a way that the hydrophobic aromatic groups penetrate into the micelle core, whereas the charged carboxylic moieties point radially outward into the aqueous phase. When the salt to surfactant molar ratio is >1, the micelle corona is almost saturated with Sal
. The reduction in micelle
water interfacial tension generated by the adsorption of the amphiphilic salt is seen to induce a sphere to rod shape transition. Predictions of micelle shape transitions obtained from different explicit solvent models, namely, standard CG water and polarizable CG water, are in qualitative agreement with each other. Standard CG water, which allows for an order of increase in computational speeds, can therefore be used to study the molecular dynamics of systems where polarizability of solvent is not a dominant factor.

’ ASSOCIATED CONTENT Supporting Information. Gromacs force field files of CTAþ cation and Sal
anion are provided for the atomistic and CG potentials; video showing sphere to rod transition of the CG surfactant micelle at [NaSal]/[CTAC] = 2.0 from simulation set III is included. This material is available free of charge via the Internet at http://pubs.acs.org.

bS

’ AUTHOR INFORMATION Corresponding Author

*E-mail: [email protected].

’ ACKNOWLEDGMENT The authors acknowledge the National Science Foundation grant CBET-1049454 for partial support of this research and Dr. Shikha Nangia for helpful discussions and input. This research utilized resources at the New York Center for Computational Sciences at Stony Brook University/Brookhaven National Laboratory which is supported by the U.S. Department of Energy under Contract No. DE-AC02-98CH10886 and by the State of New York. ’ REFERENCES (1) Israelachvili, J. N.; Mitchell, D. J.; Ninham, B. W. J. Chem. Soc., Faraday Trans. II 1976, 72, 1525. (2) Israelachvili, J. Intermolecular and Surface Forces, 2nd ed.; Academic Press: New York, 1992. (3) Larson, R. The Structure and Rheology of Complex Fluids (Topics in Chemical Engineering); Oxford University Press: New York, 1998. (4) Cates, M. E.; Candau, S. J. J. Physi.: Condens. Matter 1990, 2, 6869. (5) Turner, M. S.; Cates, M. E. J. Phys.: Condens. Matter 1992, 4, 3719. (6) Rehage, H.; Hoffmann, H. J. Phys. Chem. 1988, 92, 4712. (7) Wunderlich, I.; Hoffmann, H.; Rehage, H. Rheol. Acta 1987, 26, 532. (8) Hu, Y. T.; Boltenhagen, P.; Pine, D. J. J. Rheol. 1998, 42, 1185. (9) Hu, Y. T.; Wang, S. Q.; Jamieson, A. M. J. Rheol. 1993, 37, 531. (10) Hartmann, V.; Cressely, R. Europhys. Lett. 1997, 40, 691. (11) Vasudevan, M.; Shen, A.; Khomami, B.; Sureshkumar, R. J. Rheol. 2008, 52, 527. (12) Templeton, N. S.; Lasic, D. D.; Frederik, P. M.; Strey, H. H.; Roberts, D. D.; Pavlakis, G. N. Nat. Biotechnol. 1997, 15, 647. (13) Gao, X.; Huang, L. Gene Ther. 1995, 2, 710. (14) Morrissey, D. V.; Lockridge, J. A.; Shaw, L.; Blanchard, K.; Jensen, K.; Breen, W.; Hartsough, K.; Machemer, L.; Radka, S.; Jadhav, V.; Vaish, N.; Zinnen, S.; Vargeese, C.; Bowman, K.; Shaffer, C. S.; Jeffs, 6636

dx.doi.org/10.1021/la2006315 |Langmuir 2011, 27, 6628–6638

Langmuir L. B.; Judge, A.; MacLachlan, I.; Polisky, B. Nat. Biotechnol. 2005, 23, 1002. (15) Gregoriadis, G. Trends Biotechnol. 1995, 13, 527. (16) Yang, J. Curr. Opin. Colloid Interface Sci. 2002, 7, 276. (17) Tseng, T.-C.; McGarrity, E. S.; Kiel, J. W.; Duxbury, P. M.; Mackay, M. E.; Frischknecht, A. L.; Asokan, S.; Wong, M. S. Soft Matter 2010, 6, 1533. (18) Duque, J. G.; Eukel, J. A.; Pasquali, M.; Schmid, H. K. J. Phys. Chem. C 2009, 113, 18863. (19) Yang, S. M.; Kim, W. J. Adv. Mater. 2001, 13, 1191. (20) Li, M.; Schnablegger, H.; Mann, S. Nature 1999, 402, 393. (21) Taylor, K. M. L.; Kim, J. S.; Rieter, W. J.; An, H.; Lin, W. L.; Lin, W. B. J. Am. Chem. Soc. 2008, 130, 2154. (22) Nettesheim, F.; Liberatore, M. W.; Hodgdon, T. K.; Wagner, N. J.; Kaler, E. W.; Vethamuthu, M. Langmuir 2008, 24, 7718. (23) Helgeson, M. E.; Hodgdon, T. K.; Kaler, E. W.; Wagner, N. J.; Vethamuthu, M.; Ananthapadmanabhan, K. P. Langmuir 2010, 26, 8049. (24) Bagaria, H. G.; Kini, G. C.; Wong, M. S. J. Phys. Chem. C 2010, 114, 19901. (25) Fern
andez, C.; W€uthrich, K. FEBS Lett. 2003, 555, 144. (26) Bond, P. J.; Sansom, M. S. P. J. Am. Chem. Soc. 2006, 128, 2697. (27) Bond, P. J.; Cuthbertson, J. M.; Deol, S. S.; Sansom, M. S. P. J. Am. Chem. Soc. 2004, 126, 15948. (28) Inooka, H.; Ohtaki, T.; Kitahara, O.; Ikegami, T.; Endo, S.; Kitada, C.; Ogi, K.; Onda, H.; Fujino, M.; Shirakawa, M. Nat. Struct. Biol. 2001, 8, 161. (29) Coles, M.; Bicknell, W.; Watson, A. A.; Fairlie, D. P.; Craik, D. J. Biochemistry 1998, 37, 11064. (30) Marrink, S. J.; Tieleman, D. P.; Mark, A. E. J. Phys. Chem. B 2000, 104, 12165. (31) Tieleman, D. P.; van der Spoel, D.; Berendsen, H. J. C. J. Phys. Chem. B 2000, 104, 6380. (32) Bruce, C. D.; Senapati, S.; Berkowitz, M. L.; Perera, L.; Forbes, M. D. E. J. Phys. Chem. B 2002, 106, 10902. (33) Hayter, J. B.; Penfold, J. Colloid Polym. Sci. 1983, 261, 1022. (34) Raghavan, S. R.; Kaler, E. W. Langmuir 2001, 17, 300. (35) Imae, T.; Ikeda, S. Colloid Polym. Sci. 1987, 265, 1090. (36) Ikeda, S. Colloid Polym. Sci. 1991, 269, 49. (37) Johnson, S. B.; Drummond, C. J.; Scales, P. J.; Nishimura, S. Colloids Surf. A: Physicochem. Eng. Aspects 1995, 103, 195. (38) Magid, L. J.; Han, Z.; Warr, G. G.; Cassidy, M. A.; Butler, P. D.; Hamilton, W. A. J. Phys. Chem. B 1997, 101, 7919. (39) Clausen, T. M.; Vinson, P. K.; Minter, J. R.; Davis, H. T.; Talmon, Y.; Miller, W. G. J. Phys. Chem. 1992, 96, 474. (40) Kim, W. J.; Yang, S. M. Langmuir 2000, 16, 4761. (41) Mohanty, S.; Davis, H. T.; McCormick, A. V. Langmuir 2001, 17, 7160. (42) Maillet, J. B.; Lachet, V.; Coveney, P. V. Phys. Chem. Chem. Phys. 1999, 1, 5277. (43) Piotrovskaya, E. M.; Vanin, A. A.; Smirnova, N. A. Mol. Phys. 2006, 104, 3645. (44) Yakovlev, D. S.; Boek, E. S. Langmuir 2007, 23, 6588. (45) Wang, Z. W.; Larson, R. G. J. Phys. Chem. B 2009, 113, 13697. (46) Lorenz, C. D.; Hsieh, C.-M.; Dreiss, C. C. A.; Lawrence, M. J. Langmuir 2010, 27, 546. (47) Smit, B.; Hilbers, P. A. J.; Esselink, K.; Rupert, L. A. M.; Vanos, N. M.; Schlijper, A. G. Nature 1990, 348, 624. (48) Padding, J. T.; Boek, E. S.; Briels, W. J. J. Phys.: Condens. Matter 2005, 17, S3347. (49) Padding, J. T.; Boek, E. S.; Briels, W. J. J. Chem. Phys. 2008, 129, 074903. (50) Padding, J. T.; Briels, W. J.; Stukan, M. R.; Boek, E. S. Soft Matter 2009, 5, 4367. (51) Boek, E. S.; Jusufi, A.; Lowen, H.; Maitland, G. C. J. Phys.: Condens. Matter 2002, 14, 9413. (52) Boek, E. S.; Padding, J. T.; Anderson, V. J.; Briels, W. J.; Crawshaw, J. P. J. Non-Newtonian Fluid Mech. 2007, 146, 11.

ARTICLE

(53) Sanders, S. A.; Panagiotopoulos, A. Z. J. Chem. Phys. 2010, 132, 114902. (54) Smit, B.; Hilbers, P. A. J.; Esselink, K.; Rupert, L. A. M.; Vanos, N. M.; Schlijper, A. G. J. Phys. Chem. 1991, 95, 6361. (55) Cheong, D. W.; Panagiotopoulos, A. Z. Langmuir 2006, 22, 4076. (56) Davis, J. R.; Panagiotopoulos, A. Z. Mol. Phys. 2009, 107, 2359. (57) Mackie, A. D.; Onur, K.; Panagiotopoulos, A. Z. J. Chem. Phys. 1996, 104, 3718. (58) Floriano, M. A.; Caponetti, E.; Panagiotopoulos, A. Z. Langmuir 1999, 15, 3143. (59) Larson, R. G. J. Chem. Phys. 1988, 89, 1642. (60) Jones, J. E. Proc. R. Soc. London, Ser. A 1924, 106, 463. (61) Lazaridis, T.; Mallik, B.; Chen, Y. J. Phys. Chem. B 2005, 109, 15098. (62) Allen, E. C.; Rutledge, G. C. J. Chem. Phys. 2009, 130. (63) Angelikopoulos, P.; Bock, H. J. Phys. Chem. B 2008, 112, 13793. (64) Morisada, S.; Shinto, H. J. Phys. Chem. B 2010, 114, 6337. (65) Morisada, S.; Shinto, H.; Higashitani, K. J. Phys. Chem. B 2005, 109, 11762. (66) Shinto, H.; Morisada, S.; Miyahara, M.; Higashitani, K. Langmuir 2004, 20, 2017. (67) Verde, A. V.; Frenkel, D. Soft Matter 2010, 6, 3815. (68) Pasquali, M. Nat. Mater. 2010, 9, 381. (69) Ketner, A. M.; Kumar, R.; Davies, T. S.; Elder, P. W.; Raghavan, S. R. J. Am. Chem. Soc. 2007, 129, 1553. (70) Marrink, S. J.; Risselada, H. J.; Yefimov, S.; Tieleman, D. P.; de Vries, A. H. J. Phys. Chem. B 2007, 111, 7812. (71) Baron, R.; de Vries, A. H.; H€unenberger, P. H.; van Gunsteren, W. F. J. Phys. Chem. B 2006, 110, 15602. (72) Monticelli, L.; Kandasamy, S. K.; Periole, X.; Larson, R. G.; Tieleman, D. P.; Marrink, S. J. J. Chem. Theory Comput. 2008, 4, 819. (73) Wu, Z.; Cui, Q. A.; Yethiraj, A. J. Phys. Chem. B 2010, 114, 10524. (74) Hess, B.; Kutzner, C.; van der Spoel, D.; Lindahl, E. J. Chem. Theory Comput. 2008, 4, 435. (75) Schuler, L. D.; Daura, X.; Van Gunsteren, W. F. J. Comput. Chem. 2001, 22, 1205. (76) Berendsen, H. J. C.; Grigera, J. R.; Straatsma, T. P. J. Phys. Chem. 1987, 91, 6269. (77) Becke, A. D. J. Chem. Phys. 1993, 98, 1372. (78) Vosko, S. H.; Wilk, L.; Nusair, M. Can. J. Phys. 1980, 58, 1200. (79) Stephens, P. J.; Devlin, F. J.; Chabalowski, C. F.; Frisch, M. J. J. Phys. Chem. 1994, 98, 11623. (80) Lee, C. T.; Yang, W. T.; Parr, R. G. Phys. Rev. B 1988, 37, 785. (81) Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Scuseria, G. E.; Robb, M. A.; Cheeseman, J. R.; Montgomery, J. A.; Vreven, T.; Kudin, K. N.; Burant, J. C.; Millam, J. M.; Iyengar, S. S.; Tomasi, J.; Barone, V.; Mennucci, B.; Cossi, M.; Scalmani, G.; Rega, N.; Petersson, G. A.; Nakatsuji, H.; Hada, M.; Ehara, M.; Toyota, K.; Fukuda, R.; Hasegawa, J.; Ishida, M.; Nakajima, T.; Honda, Y.; Kitao, O.; Nakai, H.; Klene, M.; Li, X.; Knox, J. E.; Hratchian, H. P.; Cross, J. B.; Bakken, V.; Adamo, C.; Jaramillo, J.; Gomperts, R.; Stratmann, R. E.; Yazyev, O.; Austin, A. J.; Cammi, R.; Pomelli, C.; Ochterski, J. W.; Ayala, P. Y.; Morokuma, K.; Voth, G. A.; Salvador, P.; Dannenberg, J. J.; Zakrzewski, V. G.; Dapprich, S.; Daniels, A. D.; Strain, M. C.; Farkas, O.; Malick, D. K.; Rabuck, A. D.; Raghavachari, K.; Foresman, J. B.; Ortiz, J. V.; Cui, Q.; Baboul, A. G.; Clifford, S.; Cioslowski, J.; Stefanov, B. B.; Liu, G.; Liashenko, A.; Piskorz, P.; Komaromi, I.; Martin, R. L.; Fox, D. J.; Keith, T.; Laham, A.; Peng, C. Y.; Nanayakkara, A.; Challacombe, M.; Gill, P. M. W.; Johnson, B.; Chen, W.; Wong, M. W.; Gonzalez, C.; Pople, J. A. Gaussian 03, revision C.02; Guassian, Inc.: Wallingford, CT, 2004. (82) Bayly, C. I.; Cieplak, P.; Cornell, W. D.; Kollman, P. A. J. Phys. Chem. 1993, 97, 10269. (83) Dupradeau, F.-Y.; Pigache, A.; Zaffran, T.; Savineau, C.; Lelong, R.; Grivel, N.; Lelong, D.; Rosanski, W.; Cieplak, P. Phys. Chem. Chem. Phys. 2010, 12, 7821. (84) Schuttelkopf, A. W.; van Aalten, D. M. F. Acta Crystallogr., Sect. D: Biol. Crystallogr. 2004, 60, 1355. 6637

dx.doi.org/10.1021/la2006315 |Langmuir 2011, 27, 6628–6638

Langmuir

ARTICLE

(85) Martinez, L.; Andrade, R.; Birgin, E. G.; Martinez, J. M. J. Comput. Chem. 2009, 30, 2157. (86) Reekmans, S.; Bernik, D.; Gehlen, M.; Vanstam, J.; Vanderauweraer, M.; Deschryver, F. C. Langmuir 1993, 9, 2289. (87) Roelants, E.; De Schryver, F. C. Langmuir 1987, 3, 209. (88) Roelants, E.; Gelade, E.; Vanderauweraer, M.; Croonen, Y.; Deschryver, F. C. J. Colloid Interface Sci. 1983, 96, 288. (89) Almgren, M.; Garamus, V. M. J. Phys. Chem. B 2005, 109, 11348. (90) Daraio, M. E.; V€olker, A.; Aramendia, P. F.; Rom
an, E. S. Photochem. Photobiol. 1998, 67, 371. (91) Nos
e, S. Mol. Phys. 1984, 52, 255. (92) Hoover, W. G. Phys. Rev. A 1985, 31, 1695. (93) Parrinello, M.; Rahman, A. Phys. Rev. Lett. 1980, 45, 1196. (94) Parrinello, M.; Rahman, A. J. Appl. Phys. 1981, 52, 7182. (95) Allen, M. P.; Tildesley, D. J. Computer Simulation of Liquids; Oxford University Press: New York, 1989. (96) Hess, B.; Bekker, H.; Berendsen, H. J. C.; Fraaije, J. J. Comput. Chem. 1997, 18, 1463. (97) Miyamoto, S.; Kollman, P. A. J. Comput. Chem. 1992, 13, 952. (98) Darden, T.; York, D.; Pedersen, L. J. Chem. Phys. 1993, 98, 10089. (99) Bussi, G.; Donadio, D.; Parrinello, M. J. Chem. Phys. 2007, 126, 014101.

6638

dx.doi.org/10.1021/la2006315 |Langmuir 2011, 27, 6628–6638