Molecular Dynamics Simulations of Threadlike ... - ACS Publications

May 28, 2009 - This effect indicates that the electrolyte ions drive the micellar shape transition by screening the electrostatic repulsions between t...
27 downloads 0 Views 592KB Size
J. Phys. Chem. B 2009, 113, 13697–13710

13697

Molecular Dynamics Simulations of Threadlike Cetyltrimethylammonium Chloride Micelles: Effects of Sodium Chloride and Sodium Salicylate Salts† Zuowei Wang‡ and Ronald G. Larson* Department of Chemical Engineering, UniVersity of Michigan, Ann Arbor, Michigan 48109-2136 ReceiVed: February 20, 2009; ReVised Manuscript ReceiVed: April 30, 2009

We use atomistic molecular dynamics simulations to probe the effects of added sodium chloride (NaCl) and sodium salicylate (NaSal) salts on the spherical-to-threadlike micelle shape transition in aqueous solutions of cetyltrimethylammonium chloride (CTAC) surfactants. Long threadlike micelles are found to be unstable and break into spherical micelles at low concentrations of NaCl, but remain stable for 20 ns above a threshold value of [NaCl] ≈ 3.0 M, which is about 2.5 times larger than the experimental salt concentration at which the transition between spherical and rodlike micelles occurs. The chloride counterions associate weakly on the surface of the CTAC micelles with the degree of counterion dissociation decreasing slightly with increasing [NaCl] on spherical micelles, but dropping significantly on the threadlike micelles at high [NaCl]. This effect indicates that the electrolyte ions drive the micellar shape transition by screening the electrostatic repulsions between the micellar headgroups. The aromatic salicylate counterions, on the other hand, penetrate inside the micelle with their hydrophilic groups staying in the surfactant headgroup region and the hydrophobic groups partially embedded into the hydrophobic core of the micelle. The strong association of the salicylate ions with the surfactant headgroups leads to dense packing of the surfactant molecules, which effectively reduces the surface area per surfactant, and increases intramicellar ordering of the surfactant headgroups, favoring the formation of long threadlike micelles. Simulation predictions of the geometric and electrostatic properties of the spherical and threadlike micelles are in good agreement with experiments. 1. Introduction Surfactant molecules in aqueous solution self-assemble into a wide variety of morphological structures ranging from spherical and cylindrical micelles to lamellar sheets and exotic ordered bicontinuous phases.1-7 The size, shape, and shape transitions of the aggregates depend sensitively on the molecular structure and concentration of surfactant, nature of counterion, presence of added salts, and temperature. Adding electrolyte or aromatic salts into solutions of cationic surfactants, such as alkyltrimethylammonium halide and alkylpridinium halide, can induce a spherical-to-cylindrical micellar shape transition.8-16 Long equilibrium threadlike or wormlike micelles are sometimes formed above a critical salt concentration, and the solutions are accordingly transformed from low-viscosity Newtonian liquids into highly viscoelastic solutions.4,11,13-15,17-20 At surfactant concentrations below the overlap threshold of the threadlike micelles, such solutions exhibit shear thickening with a large increase in viscosity under shear or extensional flows, attributed to a flow-induced gel-like structure.17,21-27 At high enough surfactant concentrations, these micelles can entangle with each other, leading to pronounced viscoelastic behavior with a spectrum of relaxation times, analogous to that observed in semidilute solutions of long synthetic polymers.3,4,10,13,18 A further increase in the salt concentration can lead to a regime in which the relaxation spectrum is dominated by a single relaxation time.4,6,13,28 The rich phase and rheological behavior of self-assembling threadlike micellar solutions lend them practical applications as drag reducing agents, hydraulic fractur†

Part of the “H. Ted Davis Special Section”. * To whom correspondence should be addressed. E-mail: rlarson@ umich.edu. ‡ E-mail address: [email protected].

ing fluids for use in the oilfield,29 and as templates for structuring colloids and ceramics.30 Theoretically, constitutive models based on a combination of the tube theory for entangled polymers with the reversible kinetics of self-assembly have been developed by Cates and co-workers to describe the linear and steady-state nonlinear viscoelastic properties of threadlike micelle solutions.3,14,31,32 Several distinctive features of the rheology in the strongly entangled regime, including extreme shear thinning, can be well explained by these models.32 Analytical models17 and mesoscopic simulations33-35 have also been proposed to correlate the micelle breakage and recombination kinetics to the flow-induced structures. Despite extensive experimental and theoretical work on overall structure and rheological behavior, less progress has been made in understanding the microscopic mechanisms of the structural formation and transition of surfactant micelles under various salt and flow conditions, which are necessary for building a microscopic foundation for the theoretical models and linking their input parameters to the molecular properties of real surfactant systems. In this work, we investigate the key role of added salts in driving the spherical-to-cylindrical micelle shape transition and in stabilizing the threadlike micelles by performing atomistic molecular dynamics (MD) simulations. The effects of salt species and salinity are examined systematically by adding two types of salts, sodium chloride (NaCl) and sodium salicylate (NaSal), at various concentrations into the aqueous solutions of cetyltrimethylammonium chloride (CTAC) surfactants. Experiments on solutions of CTAC10,13,16,36-38 and cetyltrimethylammonium bromide (CTAB)8,9,11,15,19,23,24,27,39-42 indicate that aromatic anions like salicylate (Sal-) are much more effective in inducing threadlike micelle formation than are

10.1021/jp901576e CCC: $40.75  2009 American Chemical Society Published on Web 05/28/2009

13698

J. Phys. Chem. B, Vol. 113, No. 42, 2009

simple inorganic anions, such as chloride (Cl-) and bromide (Br-). The threshold salt concentration for the appearance of rodlike or threadlike micelles was found to be about 1.18 M for NaCl10 and 0.06 M for NaBr,8 but as low as 0.15 mM for NaSal.11,13,43 Different ion binding mechanisms have been suggested to elucidate these phenomena.8,10,11,13,15,16,19,39,43-45 The halide ions, such as Cl- or Br-, are believed to be adsorbed on the surface of the micelles due to their weakly or moderately association with the surfactant cations. The screening of the electrostatic repulsion between the surfactant head groups produced by anion binding leads to a reduction of the surface area per surfactant, and consequently drives the formation of rod- or threadlike micelles. On the other hand, nuclear magnetic resonance (NMR) experiments suggest that aromatic counterions like Sal- penetrate into the micelle and bind strongly with the surfactant head groups, promoting the micelle shape transition at much lower salt concentrations.11,13,39,43,44,46 Some of the most important theoretical and experimental work in this area have been carried out by Ted Davis and co-workers.13,15,47 For instance, Mohanty, Davis, and McCormick combined free energy models with atomistic-level Monte Carlo (MC) simulations to examine the interaction mechanism between the CTAB micelles and NaSal salt.47 They modeled the micelle as a compressible spherical or cylindrical hydrocarbon fluid core surrounded by a shell representing the rest of the micelle, including the CTA+ head groups, the Sal- ions, and part of the hydrocarbon tail. The interactions of the headgroups with counterions and water molecules were modeled in a mean-field manner. Penetration of Sal- ions into the micellar shell was found to enhance the intramicellar ordering of the headgroups by shielding the electrostatic repulsions between them, which in turn cause a sharper decrease of the free energy of cylindrical micelles than of the spherical ones, hence favoring the formation of wormlike micelles. These results qualitatively concur with the experimental observations. The goal of our work is to augment the model of Mohanty et al. by examining the kinetics of the geometric shape fluctuations and transitions of the micelles using more realistic atomistic simulation algorithms and running for much longer times. A few atomistic simulations have already contributed to the understanding of the structural and dynamic properties of wormlike or threadlike surfactant micelles. In one of the first fully atomistic MD simulations, Maillet et al. investigated n-nonyltrimethylammonium chloride (C9TAC) and erucyl bis[2-hydroxyethyl] methylammonium chloride (EMAC) micelles in aqueous solution.48 In the absence of salt, the preconstructed infinite cylindrical micelle (in a periodic box) of short-tail surfactant C9TAC was found to break up into smaller spherical micelles after 700 ps, but that of the long-chained EMAC retained its cylindrical shape throughout the simulation run of 1.85 ns. At a NaCl concentration of 1.11 M, the cylindrical C9TAC micelle still broke into finite-size segments each with a roughly cylindrical shape after a run of 2.5 ns. These short segments were taken as examples of the formation of finite, rodlike micelles, but their transition into spherical micelles might have been prevented by the finite size of the simulation box containing only 50 surfactants and the limited simulation time. Marrink et al. have pointed out in their MD simulations of the spontaneous micelle formation of dodecylphosphocholine surfactants in water that if the side length of the simulation box is comparable to the size of a spherical micelle, the surfactants would prefer to aggregate into a periodic cylindrical micelle rather than a spherical one to reduce the surface area per surfactant.49 Thus a small simulation box could impose an

Wang and Larson

Figure 1. Structure of cetyltrimethylammonium chloride (CTAC) molecule.

Figure 2. Structure of salicylate ion.

artificial energy barrier to the dissociation of the wormlike micelle. Yakovlev and Boek used a large simulation box to study the stability of wormlike micelles consisting of up to 320 cationic octyltrimethylammonium chloride (C8TAC) mixed with up to 320 anionic sodium oleate (NaOA) at different molar ratios.50 As with C9TAC,48 the infinite wormlike micelles formed by pure C8TAC were unstable and transformed into smaller spherical micelles with and without added NaCl salt. But the mixed cationic/anionic worms were quite stable because of the electrostatic attractions between the oppositely charged headgroups. The associations of various types of additives, like sodium benzoate, 2-propanol, and acetone, with the CTAC micelles have also been analyzed by Piotrovskaya et al. using MD simulations.51 Although the length and time scales accessible to atomistic simulations are still not large enough to address the macroscopic rheological properties, they do provide a microscopic basis for the development of more coarse-grained simulation models. Along this line, Padding, Boek, Briels and co-workers have obtained from atomistic simulations input parameters for a mesoscopic Brownian dynamics model which they used to study the dynamic and rheological behavior of wormlike micelle fluids.52-55 We started our MD simulations from preassembled periodic and therefore infinitely long cylindrical micelles.50,51,53 The effect of the NaCl and NaSal salts on inducing the wormlike micelle formation were investigated from the stabilization and shape transition of these preassembled micelles at various salt concentrations. Experimentally conjectured binding mechanisms for the Cl- and Sal- counterions were examined in detail by analyzing their radial density distributions relative to the micelle center, the degree of counterion dissociation, and the penetration depth and orientation of the aromatic salicylate molecules. The details of the simulation method and the system setup can be found in Section 2. Simulation results are presented and discussed in Section 3, and we draw conclusions in Section 4. 2. Methods Force-field parameters. We employed the GROMOS96 45a3 force field for all our simulations.56 In our recent MD simulation work, the use of this GROMOS forcefield in combination with the SPC (simple point charge) water model has been shown to provide structural properties of spherical sodium dodecyl sulfate (SDS) micelles in good agreement with other simulations and experimental data.57 The initial coordinates and the forcefield parameters for the CTA+ surfactant monomers (Figure 1) and salicylate ions (Figure 2) were generated using the ProDRG server.58 The aromatic hydrogens on the salicylate molecules are included explicitly. Table 1 lists the atomic charges on these

MD Simulations of Threadlike CTAC Micelles

J. Phys. Chem. B, Vol. 113, No. 42, 2009 13699

TABLE 1: Charge Assignments for the CTAC Molecule (Upper Two Rows) and the Salicylate Ion61 (Lower Rows), Respectively atom

charge (e)

atom

charge (e)

atom

charge (e)

CN1 CN2 O1 O2 O3 C1 C2

0.25 0.25 -0.843 -0.658 -0.731 -0.292 0.425

CN3 C1 C3 C4 C5 C6 C7

0.25 0.25 -0.307 -0.087 -0.237 -0.079 0.803

Cl N H2 H3 H4 H5 H6

-1.0 0 0.458 0.167 0.131 0.131 0.119

molecules. The charge assignment for CTA+ is the same as that has been used in the MD simulations of decyltrimethylammonium bromide (C10TAB) micelles,59,60 where the three methyl groups in the headgroup and the methylene group (C1) adjacent to the nitrogen carry a positive charge of + 0.25e each and the central nitrogen atom is set to be neutral. The partial charges on the salicylate atoms were taken from Song et al.61 as estimated from their quantum mechanical/molecular (QM/MM) calculations. Similar charge distribution on salicylate was obtained by Mohanty et al. using the commercial software package INSIGHT.47 Simulation predictions for the interactions of salicylate ions with dipalmitoylphosphatidylcholine (DPPC) bilayers61 and with CTAB micelles47 made by utilizing these rather similar partial charge values were shown to be in good agreement with experimental observations, even though different types of force fields have been applied in various simulation studies. In the current work, the choice of the partial charge distributions on the CTAC and salicylate molecules in combination with GROMOS96 45a3 force field is validated by the good agreement of our simulation results presented in the next section with available experimental and simulation data. More details about the forcefield and MD simulation parameters can be found in the Gromacs files (.itp, .top and .mdp files) provided in the Supporting Information. In our study the chloride counterion was used instead of the experimentally more generally studied bromide ion.50,51 Experiments showed that the degree of counterion binding increases with the radius of the hydrated ions;62 thus a higher Cl- ion concentration is required to induce the sphere-to-threadlike micelle transition or to stabilize a threadlike micelle of a given aggregation number than for Br- ions at the same concentration of CTA+ surfactants.8,10 But the general physical behavior proclaimed by the two types of counterions should be qualitatively the same. They both drive micelle formation and shape transitions by associating with the cationic surfactant headgroups and shielding the electrostatic repulsions between them. The use of the weakly binding chloride counterions can also emphasize the contrast with effect of the aromatic Sal- ions. Preparation of Initial Configuration. Since it could be very time-consuming to obtain the spontaneous micellization of the surfactants with long hydrocarbon tails from a random initial configuration, we started our simulations from preassmbled cylindrical micelles, as have been done in most other MD simulations of threadlike surfactant micelles.48,50-53,63 The cylindrical micelle was assembled by stacking 20 round slices of thickness 0.5 nm along the z-direction with each slice containing 9 CTA+ molecules with their principal axes aligned along the radial direction and evenly separated by an angle of 40°. The slices were rotated by an angle of 20° around the centerline of the micelle relative to their neighbors. In this way the average number of surfactants per nanometer was 18, which is very close to the experimental value of 16-18 in aqueous CTAC/NaCl and CTAB/NaBr solutions.8-10 The radius of the

inner void of the micelle, which measures the radial distance of the end carbon atom (C16) of the CTA hydrocarbon tail from the centerline, was set to 3 Å. A small change of the void radius has no effect on the final simulation results, because the hydrophobic micellar interior quickly squeezes out any water molecules that were initially introduced into the void space. The preassembled cylindrical micelle consisting of 180 CTA+ surfactants was then solvated with SPC water64 in a simulation box of dimensions 7.5 nm × 7.5 nm × 10 nm with its central axis along the z-direction. After equilibration, the surfactant molar concentration CD was around 0.53-0.55 M, depending on the amount of salt and water molecules added, which is well above the cmc (2.5 mM) of CTAC in water.10 Periodic boundary conditions were applied in all three directions. We are thus effectively considering an infinitely long threadlike micelle. The simple ions (Cl- and Na+) were added into the system by randomly replacing them for water molecules. All Sal- ions were initially placed outside the micelle at random orientations and positions. Energy minimization was performed before and after the addition of water and ions to keep the maximum force on any atom below 1000 KJ mol-1 nm-1. Further position restraint runs of duration 20 ps were employed by harmonically restraining the headgroup atoms to allow water molecules and ions to relax around the surfactant molecules. The resulting cylindrical micelle is shown in Figure 3 in a side-view (a) and a z-axial view (b). Simulation Details. All MD simulations were performed using the GROMACS package (version 3.3.1).65 A constant temperature of 300 K and a constant pressure of 1 bar were maintained by the Berendsen bath coupling scheme.66 For the convenience of simulating threadlike micelles, a semi-isotropic pressure coupling was applied, where the simulation cell was rescaled independently in the z-direction from that in the x- and y-directions. All bond lengths were constrained by the LINCS algorithm.67 The Ryckaert-Bellemans potential68 for dihedrals was used for the hydrocarbon tail of CTA+, and the aromatic ring of the salicylate ion was kept planar by using improper dihedrals. The SPC water geometry was constrained using the SETTLE algorithm.69 The Lennard-Jones interactions were subject to a cutoff distance of 1.4 nm. The long-range electrostatic interactions were computed using the Particle Mesh Ewald technique with a real space cutoff of 0.9 nm.70 A time step of 2 fs was employed. The trajectories were saved every 1 ps. The steady state of the systems were monitored by the total potential energy and the dimensions of the simulation cell. Table 2 summarizes all the systems we have studied with the details of the compositions. The total simulation time trun and the time duration tana at the end of the run, used for the data analysis after the system had entered a steady state, are also given. The range of the NaCl salt concentrations we studied is similar to that used by Imae and Ikeda in the static light scattering measurements of CTAC/NaCl solutions (CsNaCl ) 0-4.0 M).10 In cases of added NaSal, a Sal-/CTA+ molar ratio of 0-1.0 was also selected because of its relevance to rheological13,27,39 and NMR experiments,11,43 although our surfactant concentration CD is higher. 3. Results and Discussions 3.1. CTAC in Water: Spherical Micelles. The snapshots in Figure 3 illustrate the structural evolution of the CTAC micelle in water with no salt added. Although the surfactant concentration is relatively high, the initial long threadlike micelle is unstable.36 At early times, it exhibits strong fluctuations both laterally in its global conformation, similar to the thermal

13700

J. Phys. Chem. B, Vol. 113, No. 42, 2009

Wang and Larson

Figure 3. Snapshots of the CTAC-No-Salt system at t ) 0 ns [(a) side view and (b) z-axial view], 2 ns (c), and 6 ns (d). Only the CTA+ surfactant molecules are shown for clarity. In panels a, c, and d, the central simulation box is shown together with half of its periodic image in the z-direction.

TABLE 2: Systems Parameters: Number of CTAC Surfactants, SPC Water Molecules, Chloride, Salicylate, and Sodium Ions in the Central Simulation Box, the Molar Concentration of the Added Salt Cs, Total Simulation Time trun and Time Period for Doing Analysis tana, Respectively System

CTA+

Cl-

Sal-

Na+

SPC

Cs (M)

trun (ns)

tana (ns)

CTAC-No-Salt CTAC-NaCl-45 CTAC-NaCl-90 CTAC-NaCl-180 CTAC-NaCl-360 CTAC-NaCl-720 CTAC-NaCl-1080 CTAC-NaSal-45 CTAC-NaSal-90 CTAC-NaSal-180

180 180 180 180 180 180 180 180 180 180

180 225 270 360 540 900 1260 180 180 180

0 0 0 0 0 0 0 45 90 180

0 45 90 180 360 720 1080 45 90 180

14977 14887 14797 14617 14257 13537 12817 14511 14054 13386

0 0.13 0.27 0.54 1.08 2.17 3.28 0.13 0.27 0.55

10 10 10 16 22 20 20 10 10 10

4 2 2.5 6 3 5 5 4 4 4

fluctuation of a linear polymer chain, and longitudinally in the surfactant molecule density along its contour. To quantify the density fluctuation, we plot the number distribution of the CTA+ surfactant atoms along the z-axis in Figure 4. In these calculations the box was divided into 40 slices evenly spaced in the z-direction. The distribution curve obtained at t ) 2 ns shows clearly that the surfactants are grouping into two large clusters away from the initial even distribution (t ) 0 ns). At locations where the atom densities are low, the micelle surface bends inward due to the tilting of surfactant molecules; see Figure 3c. These tilted surfactants could form end-caps for the micellar segments, and leading to the breakage of the threadlike micelle.

Figure 4. Number distribution of the CTA+ surfactant atoms along the z-direction in the CTAC-No-Salt system. The numbers are binned into 40 slices in each case. The curve for t ) 0 ns was obtained by averaging over the first 25 ps, and the two others were averaged over 50 ps around each specified time (2 ns ( 25 ps and 6 ns ( 25 ps, respectively).

Defining the breakage time, τbreak, as the moment when the surfaces of the two end-caps formed at the first breaking point completely separated from each other, we obtained τbreak ) 5.0 ( 0.2 ns for the system studied in Figure 3. In fact, the CTAC thread broke almost simultaneously at two locations, resulting in two individual segments that quickly evolved into spherical micelles; see Figure 3d and the distribution curve at 6 ns in Figure 4. The two spherical micelles remained stable up to the end of the 10 ns simulation run. The dissociation of the CTAC thread into spherical micelles indicates that our periodic simulation cell is large enough to avoid the artificial stabilization of the threadlike micelles.49,50 In addition, the average number of surfactant molecules per unit length nL ) 18 nm-1 in the initial threadlike micelle is larger than that (nL ≈ 16 nm-1) of stabilized threadlike micelles at high NaCl concentrations, as will be seen below. This indicates that the breakage of this thread is not a result of an overstretched initial configuration. We have also monitored the variation of the box length in the z-direction, Lz, with time and found that during and after the breakage of the thread, Lz slowly increases with time owing to the repulsion between the surface charges of the two spherical micelles. In a 10 ns simulation run, it goes from 10 nm to around 13 nm for the CATC-No-Salt system and from 10 to 12 nm for the CTACNaCl-45 system. This implies that the system prefers to have a smaller value of nL, or in other words, favors the formation of spherical micelles, at low NaCl salt concentrations. On the other hand, when the NaSal salt is added our simulations successfully capture the decrease of Lz from its initial value of 10 to 7-8 nm with increased salt concentration in accordance with experiments. Thus the breakage of the threadlike micelle in our simulation results is not an artifact of the initial configuration.

MD Simulations of Threadlike CTAC Micelles

J. Phys. Chem. B, Vol. 113, No. 42, 2009 13701

TABLE 3: Properties of the CTAC Micelles in Aqueous NaCl Solutions: the Breakage Time τbreak of the Initial Infinitely Long Cylindrical Micelle, the Aggregation Number Nagg of the Spherical Micelles Resulting from the Breakage, and the Degree of Counterion Dissociation r, Respectively system CTAC-NaCl-45 CTAC-NaCl-90 CTAC-NaCl-180 CTAC-NaCl-360 CTAC-NaCl-720 CTAC-NaCl-1080

rod τbreak (ns)

4.7 7.2 8.6 13.4

sphere τbreak (ns)

Nagg

R

5.0 7.7 7.1 9.1 18.5 >20

100, 80 94, 86 94, 86 98, 82 105, 75 180

0.37 0.32 0.30 0.25 0.21 0.13 0.11

The aggregation numbers Nagg of the two individual micelles formed in this case are Nagg ) 100 and 80 (Table 3), falling in the experimental range of Nagg ) 81-115 for CTAC.10,36,37 The geometric properties of the larger micelle (Nagg ) 100) averaged over the last 4 ns are listed in Table 4. Very similar results are obtained for the smaller one. The aspect ratios of Imax/Imin (≈ 1.33-1.39) measured by the moments of inertia along the principal axes indicate that these micelles are of modestly ellipsoidal shape. Nonetheless, in the rest of this paper we will still call them spherical micelles to distinguish from the threadlike ones. The radius of a spherical micelle, Rsp, is defined as the root-mean-square (rms) distance of the nitrogen atoms from the micelle center of mass (COM). Our simulation results for Rsp (2.53 and 2.37 nm, respectively) are consistent with the experimental data of 2.5-3 nm for spherical CTAB micelles obtained from electron micrograph images8 and with another MD simulation value of 2.39 nm for a CTAC micelle of Nagg ) 90.51 The radius of the hydrocarbon core, Rc, measured by the rms distance from the C1 atoms to the micelle COM,57 is around 2.44 and 2.26 nm for the two spherical micelles, respectively. Using a two-shell model, Hayter and Penfold estimated from their small-angle neutron scattering experiments that the radius of the inner dry (hydrophobic) core of a CTAC spherical micelle, which contains the terminal methyl group and a fraction Rm ≈ 0.8 of the methylenes, is R1 ) 2.17 nm at high surfactant concentrations (CD ∼ 0.6) in the presence of zero and 0.1 M NaCl salt.36 To relate the experimental value of R1 to the simulation data, we plot in Figure 5a the radial density distributions of the surfactant methylene groups (hydrocarbon atoms) and water molecules as a function of distance from the micelle COM for the micelle with aggregation number Nagg ) 100. It can be seen that the interior of the micelle is water-free. The density profiles of water and the surfactant hydrocarbon tails cross at Rcrs ≈ 2.22 nm. A fraction of about 0.78 methylenes are found to be within this distance to COM. Similar result was obtained for the other micelle of Nagg ) 80 where the crossover between the two distribution curves occurs at Rcrs ≈ 2.07 nm with 81% methylenes within this radius. Thus the Rcrs value gives a better description of the radius of the “dry” core discussed in experiments,36 although the calculated local density of water molecules still shows a nonzero value in certain parts of this region presumably owing to the micelle surface fluctuations and the rotation of the ellipsoidal-shape micelle as a whole. Figure 5a also shows the radial density distributions of the surfactant nitrogen atoms and the Cl- counterions. The density profile of the nitrogen atoms represents the shell region of the surfactant headgroups. The peak of the Cldistribution is located outside that of the N atoms, showing that the chloride counterions are mostly associated on the surface of the micelle. The counterion association behavior

could be manifested more clearly from the radial distribution functions (RDF) between the nitrogen atoms and the Cl- ions, gN-Cl(r), and between the headgroup methyl groups and these ions, gCN-Cl(r) (Figure 6a,b). The results in Figure 6 were obtained by averaging over the whole system (over the two individual micelles). The electrostatic attraction between the positively charged methyl groups and the Cl- counterions leads to a well-structured gCN-Cl(r) which has the first peak around 0.39 nm, a first minimum at 0.51 nm, and a second peak at 0.59 nm. Since each nitrogen atom is surrounded by four charged methyl (methylene) groups, the N-Cl RDF demonstrates a broader peak around 0.49 nm. The gradual decrease of gN-Cl(r) after the peak is featured by a shoulderlike structure in between r ) 0.6-0.8 nm where no clear minimum could be distinguished. The gN-Cl(r) and gCN-Cl(r) functions demonstrate the formation of the counterion Stern layer around the micelle surface. The hydration shell formation surrounding the surfactant head groups is also elucidated in the RDFs between the nitrogen-water oxygen, gNOW, and the headgroup methyl-water oxygen, gCN-OW, pairs (Figure 6c,d). The N-OW RDF has the first peak at 0.45 nm, a minimum at 0.62 nm and the second peak at 0.75 nm, while those values for the CN-OW RDF are 0.35, 0.46, and 0.55 nm, respectively. The thickness of the hydration shells are normally taken to be the distances to the first minimum locations which are 0.62 and 0.46 nm, respectively. Our RDF results on the CTAC spherical micelles in water (Figure 6) are quite close to that obtained in MD simulations of the decyltrimethylammonium bromide (DTAB or C10TAB) spherical micelles.59,60 Those studies also showed a very shallow minimum in the nitrogen-bromide RDF whose precise location was difficult to determine. The degree of counterion dissociation, R, provides a quantitative description of the counterion binding on the micellar suface. The value of R could be estimated from the average number of chloride ions bound to each surfactant headgroup, nbound, which then gives R ) 1 - nbound. In simulations, nbound was calculated by counting the total number of Cl- ions located within a certain distance from any nitrogen atom and then dividing it by the total number of surfactant molecules. A natural choice for the cutoff distance would be the first minimum in the N-Cl RDF, corresponding to an integration to the first shell of the absorbed counterions. Although this first minimum is absent in the cases of zero or low salt concentrations,60 it appears at high NaCl concentrations (CsNaCl g 2.17 M) and was found to be located at around 0.68 nm; see below. We thus took 0.68 nm as the cutoff distance for calculating R throughout the paper in order to compare the degree of counterion dissociation on the same basis for various salt concentrations. In Figure 6a, this cutoff distance is located in the middle of the shoulder region in the gN-Cl(r) curve. It is also in between the two different first minimum values (0.65 nm and 0.715 nm) of gN-Br(r) estimated from the DTAB simulations by two different groups.59,60 Our calculations rendered nbound ≈ 0.66 and thus R ≈ 0.34 for the spherical micelle with Nagg ) 100, and nbound ≈ 0.60 and R ≈ 0.40 for the other micelle with Nagg ) 80, respectively. These R values are consistent with the experimental range of R ≈ 0.20-0.38 in which R generally decreases with increasing surfactant concentrations.12,36,37,71,72 The effective charges, Q ) RNagg, on the two individual spherical micelles were found to be 34.5 and 31.7e, respectively. In dilute aqueous CTAC solutions, Ikeda proposed that the net charge on a spherical micelle is around 31e, independent of the salt concentration.12 It is worthwhile to note that the use of a different cutoff distance

13702

J. Phys. Chem. B, Vol. 113, No. 42, 2009

Wang and Larson

TABLE 4: Geometric and Electrostatic Properties of the Spherical Micelles Formed in the Aqueous NaCl Solutions of CTAC After the Breakage of the Initial Threadlike Micelles: Aggregation Number Nagg, Radius of Gyration Rg, Micellar Radius Rsp Measured as the rms Distance from the Nitrogen Atom to COM of the Micelle, Aspect Ratio Imax/Imin, Average Volume W and Surface Area a0 Per Surfactant, Average Solvent Accessible Surface Area aSASA Per Surfactant Headgroup, and Effective Micelle Charge Q ) rNagg, Respectively system

Nagg

Rg (nm)

Rsp (nm)

Imax/Imin

V (Å3)

a0 (Å2)

aSASA (Å2)

Q (e)

CTAC-No-Salt CTAC-NaCl-45 CTAC-NaCl-90 CTAC-NaCl-180 CTAC-NaCl-360

100 94 94 98 105

2.02 1.96 1.96 1.99 2.03

2.53 2.48 2.47 2.50 2.55

1.39 1.30 1.25 1.25 1.32

678.3 679.7 671.5 667.9 661.5

80.4 82.2 81.6 80.1 77.8

231.3 228.7 225.5 223.2 215.7

34.5 29.9 27.3 24.1 18.5

will change the quantitative values of R and Q, for example, larger cutoffs leading to smaller R and Q. But the qualitative dependence of these values on the salt concentrations, as will be addressed below, should not be affected. 3.2. CTAC in Aqueous NaCl Solutions: Spherical-toThreadlike Micelle Shape Transition. Static light scattering experiments of Imae and Ikeda revealed that the spherical-torodlike micelle transition in aqueous NaCl solutions of CTAC occurs at a critical salt concentration of CsNaCl ) 1.18 M.10 The aggregation number of the rodlike micelles increases remarkably with increase of the salt concentration. In simulations starting with a single infinite threadlike micelle, the effect of salt on the micellar shape transition is reflected in a stabilization of the thread at various NaCl concentrations. The infinite thread is expected to break up at low salt concentrations where in experiments only spherical or short rodlike micelles are formed, significantly higher than the critical but to remain stable at CNaCl s value (1.18 M), corresponding to the formation of long threadlike micelles in experiments.

Figure 5. (a) Density distribution of different species in the spherical micelle of aggregation number Nagg ) 100 formed in the CTAC-NoSalt system as a function of the three-dimensional radial distance from the micelle center of mass; (b) radial density distribution of different species in the threadlike micelle in the CTAC-NaCl-1080 system as a function of the two-dimensional radial distance from the central axis of the micelle. C1-C16 refers to all carbon atoms in the hydrocarbon tails of the CTA+ surfactants. Values for nitrogen (N), chloride (Cl), and sodium (Na) atoms have been multiplied by 5 for clarity.

Breakage of Threadlike Micelle. Table 3 summarizes the CTAC/NaCl systems where the preassembled threadlike micelle broke up within a simulation run up to 20 ns. In these cases, the break-up kinetics of the threadlike micelle involved density fluctuations of CTA+ surfactants along the micellar contour, similar to those observed in the no-salt case and shown in Figure 3. The breakage in general occurs in a stepwise manner. The rod into a short periodically extended thread breaks at time τbreak + rodlike fragment containing all 180 CTA surfactant molecules in the simulation box. This rodlike micellar fragment has two closed end-caps, but its middle section does not have a welldefined cylindrical geometry. Instead its shape fluctuates with the central part getting thinner and thinner. The rod finally splits at a later time τsphere break into two smaller segments that quickly relax into spherical micelles. In the CTAC-No-Salt and CTAC-NaClrod sphere ≈ τbreak . At a NaCl concentration of CsNaCl ) 90 systems, τbreak 2.17 M, which is above the experimental threshold value of rod kept 1.18 M, the rodlike micellar fragment formed at τbreak fluctuating and touching its periodic neighbors. But the energy barrier between the end-caps prevented their recombination into a periodic infinite threadlike micelle. Thus no further breakage or recombination was observed in the CTAC-NaCl-720 system up to 20 ns. This suggests the existence of rodlike micelles of finite sizes in the solutions at these high salt concentrations. Since the systems we simulated are still rather small, the breakage times obtained only have qualitative meaning. The τbreak values for each system could also vary in different simulation runs starting from different initial configurations. A quantitative determination of the breakage times and their distributions would require a much larger simulation box with many micelles. Nonetheless, the simulation results in Table 3 demonstrate a trend of delayed breakage with increasing salt

Figure 6. Radial distribution functions for (a) nitrogen-chloride ion, (b) headgroup methyl-chloride ion, (c) nitrogen-water oxygen (OW), and (d) headgroup methyl-water oxygen pairs in the CTAC-No-Salt system.

MD Simulations of Threadlike CTAC Micelles

J. Phys. Chem. B, Vol. 113, No. 42, 2009 13703

concentration, indicating the effect of electrolyte salt in driving the formation of larger size micelles. Recently Sammalkorpi et al. presented MD simulations of the fission or splitting pathway of ionic SDS micelles induced by varying the concentration of CaCl2 salt.73 When the salt concentration was suddenly reduced from 0.7 M to zero, a rodlike micelle consisting of about 180 SDS molecules split into two spherical micelles of nearly equal sizes. In this fission process, the initial deformation of the rodlike micelle into a dumbbell geometry and the final rupture of the dumbbell into two spherical daughter micelles are very similar to what we observed in the breakage of CTAC rodlike micelles. But in the simulations of Sammalkorpi et al., there is also an intermediate state in which a long neck or stalk structure is formed between the two micellar heads. The SDS molecules inside the stalk almost fully interdigitated, giving this rodlike structure a diameter slightly larger than the length of a surfactant molecule. The length of the stalk is about 3.5-4 nm, which is comparable to the size of the micellar heads. Such an intermediate-state structure was not directly observed in our simulations. It would be interesting to explore whether the cationic CTAC/CTAB rodlike micelles would follow a similar fission pathway if the salt concentration is abruptly dropped to zero. Structure of Spherical Micelles. The geometric properties of the spherical micelles formed after the threadlike micelle breakage are presented in Table 4. Since the same number of CTA+ surfactants are used in all simulations, the average aggregation number Nagg () 90) does not show the experimentally observed slight increase (by a factor of ∼1.4) with the salt concentration.10 This is acceptable considering the distribution of the micellar sizes in real experimental systems. Since the two spherical micelles formed in each system are very similar in both their sizes and geometric aspect ratios, only the properties of the larger ones, whose aggregation numbers (94s105) are well within the experimental range of Nagg ) 81-115,10,36,37 are given in each case for simplicity. Very similar results were obtained from analysis of the smaller-sized micelles. The average volume, V, and average surface area, a0, per surfactant were estimated from a simple spherical geometry argument 3 V ) (4/3)πRsp /Nagg

(1)

2 a0 ) 4πRsp /Nagg

(2)

although the larger than unity values of the aspect ratios Imax/ Imin reflect the ellipsoidal shape of these micelles. The simulation results in Table 4 show that the values of V and a0 decrease , while a0 remains nearly constant. slightly with increasing CNaCl s The insensitivity of these quantities to the salt concentration indicates that the size and shape of the spherical micelles formed values are mostly determined by geometric packing at low CNaCl s of the surfactant molecules. The nearly identical packing feature ofthesphericalmicellescanalsobeseenfromtheirnitrogen-nitrogen RDFs, gNN(r), which almost sit on top of each other for the systems with CNaCl e 1.08 M; see Figure 7. The average surface s area per surfactant headgroup, a0 ≈ 81 Å2, is very close to that estimated in experiments (78 Å2)9 and in another MD simulation of spherical CTAC micelle in water (79.8 Å2).51 This is expected because both our work and theirs found similar values of the micellar radius Rsp and employed eq 2 to calculate a0. The use of eq 2 to estimate a0 assumes a perfect spherical geometry and a smooth surface of the micelle. A more realistic

Figure 7. Radial distribution function between the nitrogen atoms in the aqueous NaCl solutions of CTAC at various salt concentrations (Table 2).

description of the micellar surface could be obtained from the solvent accessible surface area (SASA). A probe “water” molecule of radius 1.4 Å is moved around the rough surface of the micelle to calculate all area accessible to it. The simulation results for the SASA per surfactant molecule, aSASA, are given in Table 4. The aSASA values are larger than the a0 values, but both show very weak surface area dependence on the NaCl concentration. The contributions to SASA from the surfactant hydrophilic headgroups and from the hydrocarbon tails are found to be nearly equal, indicating that a considerable number of methelyne and methyl end groups are exposed to water due to the bending of the hydrocarbon tails toward the micelle surface. Stabilization of Threadlike Micelle. At the highest salt concentration (CsNaCl ) 3.28 M) investigated in this work, the periodic threadlike micelle remained stable for the entire simulation run of 20 ns. In this case (CTAC-NaCl-1080), the box length in the z-direction increased slowly at early times and reached a plateau value of Lz ≈ 11.3 nm after about 15 ns. This gives an average aggregation number of surfactant molecules per unit length nL ) 16 nm-1, which is in good agreement with experimental values (16-18 nm-1).9,10 The pressure difference between the axial (z-) and lateral directions, ∆P ) Pzz - (Pxx + Pyy)/2, could be used to monitor the tensional state of the micelle, where Pii are the diagonal components of the pressure tensor. The small value of ∆P ≈ 0.64 bar in this system reflects that the threadlike micelle is in a nearly tensionless state.50 Figure 8a,b shows the snapshots of the threadlike micelle at t ) 20 ns. Density fluctuations along the micellar contour still exist (Figure 8a), and the interdigitation of the hydrocarbon tails of the surfactants in the micellar core can be seen from the crosssection view (Figure 8b). Experimentally, the average contour length of threadlike or rodlike micelles formed at NaCl concentrations CsNaCl g 3 M was found to be above 170 nm.10 Thus the breakage of the threadlike micelle in Figure 8a might not be anticipated even if the simulation run was extended for tens of nanoseconds. Structure of Threadlike Micelle. The radial density distributions of various surfactant atoms, water molecules, and Cl- and Na+ ions are plotted in Figure 5b as a function of the twodimensional radial distance of these species from the central axis of the micelle. The peak of the Cl- counterion distribution is again located beyond that of the surfactant headgroups, and the Na+ co-ions are distributed further away from the micelle surface. This reflects the formation of the electrostatic double layer around the micelle surface, which shields the inter- and intramicellar electrostatic repulsions. The radius of a cylindrical micelle, Rth, can be measured as the radial rms distance of the

13704

J. Phys. Chem. B, Vol. 113, No. 42, 2009

Wang and Larson

Figure 8. Snapshots of the threadlike micelle formed in the CTAC-NaCl-1080 system at t ) 20 ns [(a) side view and (b) axial view] and in the CTAC-NaSal-180 system at t ) 10 ns [(c) side view and (d) axial view]. Only the CTA+ surfactant (cyan) and salicylate (red) molecules are shown for clarity. The central simulation box is shown together with 0.5 (a) and 1.5 (c) of its periodic images in the z-direction, respectively, to make the length scales of the two images comparable.

TABLE 5: Geometric Properties of the Threadlike Micelles: Radius of the Cylindrical Micelle Rth, Side Lengths Lx,y,z of the Simulation Box, Average Number of Surfactant Molecules Per Unit Length nL, Average Volume W and Surface Area a0 Per Surfactant, and Average Solvent Accessible Surface Area aSASA Per Surfactant Headgroup, Respectively system

Lx,y (nm)

Lz (nm)

nL (nm-1)

Rth (nm)

V (Å3)

a0 (Å2)

aSASA (Å2)

CTAC-NaCl-1080 CTAC-NaSal-45 CTAC-NaSal-90 CTAC-NaSal-180

6.95 8.18 9.02 8.67

11.3 8.29 6.73 7.17

15.9 21.7 26.7 25.1

1.84 2.09 2.26 2.21

667.7 632.0 599.9 611.2

72.6 60.5 53.1 55.3

214.9 200.2 181.1 163.7

nitrogen atoms to the central axis of the micelle. This value was evaluated to be Rth ≈ 1.84 nm for the threadlike micelle in the CTAC-NaCl-1080 system. It is somewhat smaller than the radii of the spherical micelles (Table 4) because of the tilting and interdigitation of the surfactant molecules in the thread. The average volume and average surface area per surfactant molecule were calculated for the idealized cylindrical geometry 2 V ) πRth /nL

(3)

a0 ) 2πRth /nL

(4)

where nL is the average aggregation number of surfactant molecules per unit length. Simulation results for V and a0 are given in Table 5. Since the hydrocarbon core of the micelles are water-free, the average volume per surfactant in the cylindrical micelle is very close to that in the spherical micelles. As expected for the spherical-to-cylindrical shape transition, the average surface area occupied by a surfactant headgroup in the threadlike micelle (a0 ≈ 72.5 Å2) is smaller than that of spherical micelles (81 Å2). The solvent accessible surface area per surfactant, aSASA, is also reduced in the threadlike micelle. From the simple geometric packing argument that assumes equal radii and average volumes per surfactant V for both spherical and cylindrical micelles,1,2,6 the ratio of the a0 values in the two

cases is expected to be 1.5. In simulations, this ratio was found to be about 1.12, because the radius of the threadlike micelle is smaller than that of the spherical ones. The denser packing of the surfactant headgroups in the threadlike micelle is also reflected in the increased heights of the first peak and the shoulder at shorter distance in the nitrogen-nitrogen RDF [gN-N(r)] in Figure 7. A similar variation trend of the gN-N(r) functions with respect to the spherical-to-rodlike micellar shape transition has been observed in the MD simulations of C9TAC micelles.48 Counterion association. Association of counterions on the micellar surface, which effectively screens the intra- and intermicellar electrostatic repulsions, has been generally accepted as the microscopic mechanism for the spherical-to-cylindrical micelle shape transition induced by adding simple electrolyte salts. There have been a large number of experimental and theoretical works on the counterion binding or association, which were usually aimed at charactering the degree of counterion dissociation R or fractional charge of the micelles.12,36-38,72,74,75 For spherical micelles, R was found to decrease with increasing salt concentration, while the effective micelle charge Q () RNagg) remained roughly constant due to the slight increase of the aggregation number.12,76 Results on the R values of rodlike or threadlike micelles, however, seem to be controversial. On the basis of a linear double logarithmic relation between the micelle molecular weight and the counterion concentration,

MD Simulations of Threadlike CTAC Micelles

Figure 9. Radial density distributions of the chloride (thick lines) and sodium (thin lines) ions relative to the nitrogen atoms in the aqueous NaCl solutions of CTAC.

Ikeda proposed that for CTAC R for rodlike micelles (0.535) is higher than for spherical micelles (0.269s0.383).12 On the other hand, theoretical calculations of Ruckenstein and Beunen predicted that the degree of counterion dissociation is lower for the threadlike micelles than for the spherical ones.74 In atomistic simulations, the counterion association behavior can be examined both in terms of the “macroscopic” R-value and the microscopic spatial distribution of the counterions around the surfactant headgroups. The absorption of Cl- counterions to the micelle surface has already been indicated in Figure 5 in the radial density distribution of these ions relative to the COM of the spherical micelle or to the central axis of the threadlike micelle. The binding of the counterions to the surfactant headgroups could be further elucidated from the radial distribution of Cl- ions around the nitrogen atoms. Figure 9 presents the radial density distributions of the Cl- and Na+ ions relative to the nitrogen atoms, FN-Cl(r) and FN-Na(r). The use of the density distributions F(r), rather than the normalized RDFs g(r), allows the quantitative comparison of the degree of the counterion association. The first peak in the FN-Cl(r) (and also in gN-Cl(r)) function is located at r ) 0.49(5) nm for all these systems, but its height grows significantly with the increase of the salt concentration, indicating that more and more Cl- ions are absorbed into the Stern layer surrounding the surfactant head groups. At high salt g 2.17 M) where the rodlike or threadlike concentrations (CNaCl s micelles are formed, FN-Cl(r) [and gN-Cl(r)] becomes more structured with the appearance of the first minimum at 0.68(5) nm and a second peak at 0.78(5) nm. This reflects the formation of two counterion shells around the nitrogen atoms. In the mean time, more Na+ co-ions are also attracted to the region with high counterion density. The FN-Na(r) [and gN-Na(r)] curves show a peak around 0.56(5) nm located mostly inside the inner counterion shell.48 The density of the Na+ ions at this peak is still lower than that in the bulk because of the repulsions from the cationic head groups. The degree of counterion dissociation R was calculated in the same way as discussed for the no-salt case by counting the counterions inside the first binding shell with r e 0.68 nm. But in the presence of NaCl salt, the contributions from the Na+ ions, which were located inside this binding shell and thus effectively neutralized the same number of counterions, should be deducted when evaluating R and the effective micelle charge Q. The simulation results for R are given in Table 3. For the systems with spherical micelle formation (CsNaCl e 1.08), our R data and their slight decrease with increasing salt concentration are in good agreement with the results given by Ikeda.12 The effective charges Q of the spherical micelles (Table 4) also show

J. Phys. Chem. B, Vol. 113, No. 42, 2009 13705

Figure 10. Radial density distributions of different CTA+ surfactant atoms, water molecules, chloride and sodium ions, and all (SAL), hydrophilic (SAL, hyphil) and hydrophobic (SAL, hyprob) atoms of the salicylate ions in the CTAC-NaSal-180 system as a function of the two-dimensional radial distance from the central axis of the threadlike micelles. Values for nitrogen (N), chloride (Cl), and sodium (Na) atoms have been multiplied by 5 for clarity.

a decrease with increase of CsNaCl, because our simulations with fixed number of surfactant molecules did not capture the increase of the micellar aggregation number with increased salt concentration. For threadlike micelles, the simulation result, R ≈ 0.11, obtained in the CTAC-NaCl-1080 system is much smaller than the value of 0.535 estimated by Ikeda.12 It is also lower than that of the spherical micelles, thus qualitatively supporting the theoretical predictions of Ruckenstein and Beunen.74 The lower value of R on the threadlike micelles is physically reasonable, because a higher degree of counterion binding is needed to effectively screen the electrostatic repulsions between the cationic headgroups in order to induce the micellar shape transition from a highly curved spherical geometry to a moderately curved cylindrical geometry. Experiments on the self-assembled CTAC bilayers have also shown a significantly smaller degree of counterion dissociation for the planar bilayers than for spherical micelles.38 3.3. CTAC in Aqueous NaSal Solutions: Threadlike Micelles. Formation of CTA+ rodlike or threadlike micelles was experimentally observed at very low Sal- concentrations (∼0.15 mM).11,13,43,44,47,77 For the relatively high NaSal concentrations we studied (Table 2), the infinite threadlike micelle stayed stable in all three systems. Figure 8c,d shows the snapshots of the threadlike micelle in the equimolar CTA+/Sal- (CTAC-NaSal180) system at t ) 10 ns. The density distributions of various CTA+ atoms, water molecules, and salicylate ions as a function of the two-dimensional radial distance from the central axis of the threadlike micelle are shown in Figure 10 for this system. The thread in the NaSal solution (Figure 8c) displays a much smaller amplitude of density fluctuations than does that in the NaCl solution (Figure 8a). A comparison of the radial density distribution functions in these two cases (Figure 5b and Figure 10) also illustrates a shift of the micelle-water interface to a larger radial distance, marked by the peak position of the nitrogen atom distribution and the crossover between the density profiles of the surfactant hydrocarbon tails and water molecules. Furthermore the thread in NaSal shows a smaller surface roughness, measured by the height and width of the distribution function of the nitrogen atoms. These differences can be seen more clearly in Figure 11 and Table 5 and are discussed in the following subsection. Ordering of Surfactant Molecules. The geometric properties of the threadlike micelles formed in various salt solutions are summarized in Table 5. In the presence of NaSal salt, the steady-

13706

J. Phys. Chem. B, Vol. 113, No. 42, 2009

Wang and Larson

Figure 12. Radial distribution functions between the nitrogen atoms, gN-N(r), and between nitrogen atoms and headgroup methyl groups, gN-CN(r), in the CTAC-NaSal-180 system.

Figure 11. Enhanced ordering of the surfactant headgroups with increasing salicylate concentration (Table 2). (a) Radial density distribution functions, FN(r), of nitrogen atoms as a function of the two-dimensional radial distance from the central axis of the threadlike micelles; (b) radial distribution function, gN-N(r), between the nitrogen atoms.

state simulation box length in the z-direction shrank from the original value of Lz ) 10 nm to around 6.7-8.3 nm, indicating the more dense packing of the surfactant molecules along the micellar contour than in the NaCl solutions (CTAC-NaCl-1080). The corresponding average aggregation number per unit length is in the range of nL ) 21-27 nm-1 for Sal-/CTA+ molar ratios of 0.25-1.0. This agrees well with the experimental value of nL ≈ 27.8 nm-1 obtained by Imae in a dilute aqueous NaSal 77 19 solution of CTASal where CSal s /CD > 1.0. Shikata et al. have also estimated a higher molecular weight per unit length (ML ≈ 8900 nm-1 and correspondingly nL ≈ 21 nm-1) of the wormlike micelles formed in aqueous CTAB/NaSal solutions than that in CTAB/NaBr solutions (ML ≈ 7030 nm-1 and nL ≈ 17 nm-1).9 The linear packing density of the threadlike micelles increases with increasing NaSal concentration. This is evident when comparing the nL value for the CTAC-NaSal-45 system with those for the CTAC-NaSal-90 and CTAC-NaSal-180 systems. In the two latter cases, the slightly higher nL value or smaller box length Lz in the CTAC-NaSal-90 system could be thought to be within the simulation error bar due to the limited simulation box size and simulation time. It could also be related to the insertion of large number of salicylate ions inside the micelle which prevents a further increase of the surfactant packing density at high Sal- concentrations. In these CTAC/ NaSal systems, the pressure differences in the directions parallel and perpendicular to the micelle axis are also relatively small, given as ∆P ≈ -0.04, -2.37, and -1.06 bar, respectively. The threadlike micelles are basically in the tensionless state. The increase in the surfactant packing density is accompanied by an increase in the micellar radius Rth. The radii of the threadlike micelles stabilized in NaSal solutions are more close to, but still lower than, the radii Rsp of the spherical micelles; see Tables 4 and 5. Free energy calculations of Mohanty et al.

have shown a similar increase of the radius of cylindrical micelles with the Sal- concentration in CTAB/NaSal systems.47 The larger Rth value at higher NaSal density results from the enhanced ordering of the surfactant molecules due to the insertion of the aromatic salicylate ions. The axial-view of the threadlike micelle in the CTAC-NaSal-180 system (Figure 8d) demonstrates that the surfactant hydrocarbon tails are oriented radially with their molecular axes perpendicular to the central axis of the threadlike micelle and possess almost no interdigitation.19 The roughness of the micelle surface has also been reduced dramatically in comparison with the CTAC-NaCl-1080 system (Figure 8b). This corresponds to a narrowed headgroup shell region characterized in Figure 11a by the change of the radial density distribution function FN(r) of the nitrogen atoms in which the width of FN(r) decreases and its peak height increases with increasing NaSal concentration. The shift of the peak location to higher radial distance reflects the increase in Rth. Figure 11b presents the radial distribution functions between the nitrogen atoms, gN-N(r), in the threadlike micelles formed at various salt conditions. Three pronounced peaks appear in gN-N(r) of the CTAC/NaSal systems, while only the second peak is evident in that of the CTAC-NaCl-1080 system. The sharpening of these peaks reveals the increased ordering of the surfactant headgroups. The first peak in gN-N(r) indicates a strong coordination of the surfactant headgroups induced by the inserted Sal- ions. The negatively charged salicylate ion could bind two or more cationic surfactant headgroups to form a temporary cluster, which results in a high probability for a headgroup to find another one in its proximity. Since the direct contact of the nitrogen atoms are sterically hindered by the surrounding methyl groups, the separation between the nitrogen atoms at this peak (r ≈ 0.63 nm) is larger than their van der Waals diameter (∼0.31 nm). The RDFs in Figure 12 further illustrate that the first peak of the gN-N(r) function sits in between the two shells of methyl groups. The strong association of the surfactant headgroups with the Sal- ions leads to their close packing, which in turn reduces the average surface area a0 and the average solvent accessible surface area aSASA per surfactant, as shown in Table 5. The value of aSASA was evaluated by dividing the total SASA of the theadlike micelle in the central simulation box by the total number of CTA+ surfactants (180). Since certain parts of the solubilized salicylate ions are also in contact with water, the results for aSASA have thus included the contribution from the Sal- ions. This contribution is about 28% in the CTAC-NaSal-180 system and is smaller in the systems with lower Sal- concentrations. The simulation results for a0 and aSAS validate the hypothesis of Lin et al.15

MD Simulations of Threadlike CTAC Micelles that the headgroup area is reduced by the favorable interactions of salicylic anions with the surfactant cations. The second and third peaks in gN-N(r) imply the formation of a crystal-like structure. For a confined single layer of likely charged colloidal particles, the ground-state structure is a twodimensional (2D) hexagonal lattice. Supposing that the surfactant headgroups arrange into a quasi-2D hexagonal lattice on the micellar surface with a lattice constant b, the area of the Wigner-Seitz cell around each headgroup would be a0 ) √3b2/2. The value of b, which features the distance of a surfactant headgroup with its nearest neighbors on the lattice, could then be estimated by substituting the simulation result on a0 into this expression. The second nearest neighbors of the reference headgroup should be located at √3b. Using the CTAC-NaSal-180 system as an example, the lattice constant was found to be b ≈ 0.8 nm, which subsequently gives √3b ) 1.39 nm. In gN-N(r) of this system, the second and third peaks are located at r ≈ 0.83 and 1.42 nm, respectively. Similar results have been achieved for the CTAC-NaSal-45 and CTAC-NaSal90 systems. Considering the approximations made in estimating the surface area per surfactant headgroup a0 and in treating the curved cylindrical surface as a planar surface, the simulation results for gN-N(r) of the CTAC/NaSal systems show clear evidence for a 2D hexagonal arrangement of the surfactant headgroups on the surface of the threadlike micelle. Such an ordered structure was not observed on the spherical micelles47 or on the threadlike micelles formed in aqueous NaCl solutions we studied (CTAC-NaCl-1080). The gN-N(r) functions in Figure 11b have no apparent fourth (at r ) 2b) or higher order peaks expected for an ideal hexagonal lattice, implying that the quasi2D lattice structure is still rather local. The ordering of the surfactant headgroups at the intermediate length scales can effectively lower the repulsion energy between them and favors the formation of longer threadlike micelles. Mohanty et al. have pointed out the importance of considering the intramicellar ordering of the surfactant headgroups in estimating the free energy of surfactant micelles.47 Penetration of Salicylate Ions Inside the Micelles. The enhanced ordering of the surfactant molecules in the threadlike micelle is caused by the penetration of salicylate ions into the micelle. The snapshot in Figure 8d and the radial density distribution of the salicylate ions in Figure 10 show clearly that the Sal- ions are inserted into the interior of the micelle and form a well-defined shell spanning the micellar headgroup region and the outer part of the hydrophobic core. Similar penetration behavior of salicylate ions into the CTA+ micelles formed in equimolar CTAB/NaSal mixtures47 and into the dipalmitoylphosphatidylcholine (DPPC) bilayers61 has also observed in Monte Carlo and MD simulations, respectively. For the Sal-/ CTA+ molar ratios we studied, 95 to 100% of the Sal- ions are absorbed into the micelles. This is consistent with the NMR study of dilute CTASal solutions where the degree of the Salcounterion binding was found to be around 92%.43 In the sketch of the salicylate ion in Figure 2, the hydrophilic group of Sal- is defined as atoms O1, O2, H2, O3, and C7, and the hydrophobic group is defined as atoms C1-C6 and H3-H6.61 Figure 10 shows that the hydrophilic group is located in the micellar headgroup region, while the hydrophobic group penetrates more deeply into the hydrophobic core of the micelle. On the basis of their NMR experiments, Rao et al. have suggested that the salicylate ions are oriented radially in the micelle with their carboxylic (COO-) groups pointing out.11 This picture has been visualized in the Monte Carlo simulations of Mohanty et al.47 and can be seen clearly in the colored version

J. Phys. Chem. B, Vol. 113, No. 42, 2009 13707

Figure 13. Snapshots of the salicylate ions in the CTAC-NaSal-180 system at t ) 10 ns: (a) axial view and (b) side view. In (a), only half of the simulation box in the z-direction is shown for clarity. The red, gray, and cyan bonds correspond to the oxygen, hydrogen, and carbon atoms, respectively.

Figure 14. Radial density distributions of individual salicylate atoms and nitrogen atoms in the CTAC-NaSal-180 system as a function of the two-dimensional radial distance from the central axis of the threadlike micelle. The results for the surfactant hydrocarbon tails and water molecules are also included for reference.

of the snapshot in Figure 13a. A more quantitative examination of the Sal- ion orientation and penetration depth is shown in Figure 14, which plots the radial density distribution of the individual Sal- atoms with respect to the central axis of the threadlike micelle. These results for the locations of the

13708

J. Phys. Chem. B, Vol. 113, No. 42, 2009

individual salicylate atoms relative to different regions of the surfactant micelle could be compared directly with potential NMR experiments. The para protons H4 are shown to the closest to the micellar central axis. The radial distances of the other atoms from the micellar axis increase in the order of meta protons H5, H3, and ortho protons H6 for the hydrophobic group, and then from H2, O3, to O1 atoms for the hydrophilic group. The well-ordered arrangement of the Sal- atoms provides direct evidence for the radial orientation of the salicylate ions with their hydrophilic groups pointing outward.11,47 The hydrophilic atoms (O1, O3, H2) are distributed in the surfactant headgroup region with their F(r) extensively overlapping with the density profile of the nitrogen atoms. These atoms are in direct contact with water. The hydrophobic atoms on the π-ring spread among the headgroup region and the outer part of the water-free hydrophobic core. If we assume that the probability for finding a salicylate hydrophobic atom to be in the surfactant headgroup region is determined by the overlapping area between its radial density distribution function and that of the nitrogen atoms, this probability is about 59% for the ortho proton H6 and it decreases systematically for the meta (51% for H3 and 37% for H5) and para (34% for H4) protons. Our simulation results thus confirm the long-standing experimental hypothesis that the salicylate ions penetrate inside the CTA+ micelles with their hydrophilic groups located in the headgroup region (micelle-water interface) and their hydrophobic group (aromatic ring) partially embedded in the hydrophobic core.11,15,43,44,46 Calculation of the SASA contribution from the salicylate ions has also shown that their hydrophilic groups have a higher probability to contact water than the hydrophobic groups. In the CTAC-NaSal-180 system, aSASA per Sal- ion is about 27 Å2 for its hydrophilic group and 18 Å2 for the hydrophobic group. The orientation and positioning of the Sal- ions inside the micelle are determined by the amphiphilic character of their chemical structures. Rao et al. proposed a mechanism for the threadlike micelle formation where the carboxylic (COO-) groups of the salicylate ions protrude out the micellar surface and work as bridges to link a set of spherical micellar beads into long rodlike micelle.11 This phenomenon is neither supported by our simulations, nor by the cryo-TEM studies of Lin et al.15 The insertion of the Sal- ions into the headgroup region effectively reduces the net micelle charge and subsequently increases the degree of dissociation R of the Cl- counterions. Using the same cutoff distance of 0.68 nm as in the aqueous NaCl solutions, the estimated R value of Cl- increases with Sal- concentration from R ≈ 0.46 at CsNaSal/CD ) 0.25 (CTACNaSal-45), to 0.65 at CsNaSal/CD ) 0.5 (CTAC-NaSal-90) and finally to a fully dissociated state of R ) 1.0 in the equimolar Sal-/CTA+ system (CTAC-NaSal-180) where the micelle charge is completely neutralized by the inserted Sal- ions. Figure 10 also shows that in the CTAC-NaSal-180 system the density of Cl- counterions near the micellar surface region is very low and almost coincides with the distribution of the Na+ co-ions. Experimental measurements of Shikata et al. on the electrostatic features of threadlike micelles formed in the CTAB/NaSal solutions found exactly the same trend: the degree of dissociation of Br- increased with NaSal concentration from R ) 0.25-0.5 (depending on the surfactant concentration) at CsNaSal ) 0 up to 1.0 in the equimolar mixture.41 Thus at high Sal- concentrations, the Cl- or Br- counterions make negligible contributions to shielding the intramicellar electrostatic repulsions. Figure 15 shows the probability of finding the normal vector of the π-ring plane of a salicylate ion at an angle θ to the z-direction (the central axis of the micelle). The π-ring is

Wang and Larson

Figure 15. Probability of finding the normal vector of the π-ring plane of a salicylate ion at an angle θ in respect to the z-direction.

perpendicular to the micellar axis when θ ) 0° and in plane with it at θ ) 90°. The results in Figure 15 demonstrate that this probability has a very broad distribution and exhibits a rather weak dependence on the Sal- concentration. Averaging over the whole range of the probability distribution functions, the resulted mean values of θ were found to be 49.0, 47.0 and 47.5° for the CTAC-NaSal-45, 90, and 180 systems, respectively. There is no strong preference for the π-rings of salicylate ions to be perpendicular to the micellar axis, as could be seen in the snapshot in Figure 13b.47 The screening of the electrostatic repulsions between the headgroups due to strong Sal- ion association reduces the stiffness of the threadlike micelle, but the related increases in the surfactant packing density and the micellar diameter are likely more effective in enhancing the micellar stiffness. As a result, the persistence length of the threadlike micelles formed in the presence of Sal- should be larger than that formed in the simple electrolyte salts. Experimentally the persistence length of the CTA+ rodlike micelles formed in aqueous NaSal solutions is about 110-140 nm,77 and that in NaBr or NaCl solutions are in the range of 40-50 nm.9 4. Conclusions Atomistic molecule dynamics simulations have been performed to study the effects of different types of salts on the spherical-to-threadlike micelle shape transition in aqueous solutions of cationic cetyltrimethylammonium chloride surfactants. Two types of salts were examined. One is the simple electrolyte salt, sodium chloride, and the other is an aromatic salt, sodium salicylate. Experiments have suggested that different binding mechanisms of the Cl- and Sal- counterions with the surfactant headgroups lead to the distinct efficiencies of these ions in driving the formation of long threadlike micelles. All our simulations were started with a preassembled infinitely long cylindrical micelle in a periodic box. When no salt is added, the initial threadlike micelle breaks apart quickly into smaller segments, which then evolve into spherical micelles. The aggregation number and geometric properties of these spherical micelles obtained from the simulations are in good agreement with experimental data, thus validating the simulation model employed. In the presence of NaCl salt, the initial threadlike micelles break into spherical micelles at low salt concentrations, but the breakage time grows with the increasing salinity, and above a critical salt concentration (∼3.0 M) we see no breakage over a 20 ns simulation run. The Cl- counterions are found to associate weakly on the surface of these surfactant micelles. The degree of counterion dissociation R shows a slight decrease with increasing salt concentration in the systems where spherical

MD Simulations of Threadlike CTAC Micelles micelles are formed, while the R value for the threadlike micelles is considerably smaller than for the spherical ones. The higher degree of counterion binding at high salt concentrations effectively shields the electrostatic repulsion between the surface headgroups, leading to a reduction in the average surface area per surfactant that consequently drives the micelle shape transition toward the cylindrical geometry. The aromatic salicylate ions are more efficient than the weakly associated electrolyte counertions in driving the formation of long threadlike micelles by penetrating into the micellar interior. The hydrophilic groups of the Sal- ions are found to locate in the headgroup region of the micelle, while their hydrophobic groups (the π-rings) are partially embedded in the hydrophobic core. The simulation results thus provide clear evidence supporting the experimentally based suggestions that the salicylate ions orient radially with hydrophilic groups pointing outward and the π-ring positioned among the surfactant headgroup region and the outer part of the hydrophobic core. Long threadlike micelles remain stable in the CTAC/NaSal systems at the Sal-/ CTA+ molar ratios of 0.25-1.0 in agreement with experimental observations. The strong association of Sal- ions with the surfactant headgroups leads to a dense packing of the surfactant molecules along the micellar contour, which yields a higher aggregation number per unit length and smaller surface area per surfactant headgroup than for the threadlike micelles formed in aqueous NaCl solutions. The ordering of the surfactant headgroups is also enhanced by the insertion of the Sal- ions, represented by the narrowed width of the headgroup shell region and the formation of a local quasi-2D hexagonal lattice of the headgroups on the micellar surface. The addition of Sal- helps to lower the total free energy of micelles by both effectively shielding the electrostatic repulsion between the surfactant headgroups via their penetration into the headgroup shell and enhancing the intramicellar ordering of the headgroups. Long threadlike micelles can thus be formed at very low NaSal concentrations. Acknowledgment. Helpful discussions with Dr. Radhakrishna Sureshkumar are gratefully acknowledged. Supporting Information Available: Gromacs input files providing the forcefield parameters of the CTAC and salicylate molecules, as well as the MD run parameters. This material is available free of charge via the Internet at http://pubs.acs.org. References and Notes (1) Israelachvili, J. N.; Mitchell, D. J.; Ninham, B. W. J. Chem. Soc., Faraday Trans. 2 1976, 72, 1525–1568. (2) Israelachvili, J. N. Intermolecular and surface forces with applications to colloidal and biological systems; Academic Press: New York, 1985. (3) Cates, M. E.; Candau, S. J. J. Phys.: Condens. Matter 1990, 2, 6869–6892. (4) Rehage, H.; Hoffmann, H. Mol. Phys. 1991, 74, 933–973. (5) Laughlin, R. G. The aqueous phase behaVior of surfactants; Academic Press: New York, 1994. (6) Larson, R. G. Structure and Rheology of Complex Fluids; Oxford University Press: New York, 1999. (7) Abdel-Rahem, R. AdV. Colloid Interface Sci. 2008, 141, 24–36. (8) Imae, T.; Kamiya, R.; Ikeda, S. J. Colloid Interface Sci. 1985, 108, 215–225. (9) Imae, T.; Ikeda, S. J. Phys. Chem. 1986, 90, 5216–5223. (10) Imae, T.; Ikeda, S. Colloid Polym. Sci. 1987, 265, 1090–1098. (11) Rao, U. R. K.; Manohar, C.; Valaulikar, B. S.; Iyer, R. M. J. Phys. Chem. 1987, 91, 3286–3291. (12) Ikeda, S. Colloid Polym. Sci. 1991, 269, 49–61. (13) Clausen, T. M.; Vinson, P. K.; Minter, J. R.; Davis, H. T.; Talmon, Y.; Miller, W. G. J. Phys. Chem. 1992, 96, 474–484. (14) Spenley, N. A.; Cates, M. E.; McLeish, T. C. B. Phys. ReV. Lett. 1993, 71, 939–942.

J. Phys. Chem. B, Vol. 113, No. 42, 2009 13709 (15) Lin, Z.; Cai, J. J.; Scriven, L. E.; Davis, H. T. J. Phys. Chem. 1994, 98, 5984–5993. (16) Magid, L. J.; Han, Z.; Warr, G. G.; Cassidy, M. A.; Butler, P. D.; Hamilton, W. A. J. Phys. Chem. B 1997, 101, 7919–7927. (17) Turner, M. S.; Cates, M. E. J. Phys.: Condens. Matter 1992, 4, 3719–3741. (18) Khatory, A.; Lequeux, F.; Kern, F.; Candau, S. J. Langmuir 1993, 9, 1456–1464. (19) Shikata, T.; Dahman, S. J.; Pearson, D. S. Langmuir 1994, 10, 3470– 3476. (20) Raghavan, S. R.; Kaler, E. W. Langmuir 2001, 17, 300–306. (21) Rehage, H.; Wunderlich, I.; Hoffman, H. Prog. Colloid Polym. Sci. 1986, 72, 51–59. (22) Wunderlich, I.; Hoffmann, H.; Rehage, H. Rheol. Acta 1987, 26, 532–542. (23) Hu, Y. T.; Wang, S. Q.; Jamieson, A. M. J. Rheol. 1993, 37, 531– 546. (24) Hartmann, V.; Cressely, R. Europhys. Lett. 1997, 40, 691–696. (25) Hu, Y. T.; Boltenhagen, P.; Pine, D. J. J. Rheol. 1998, 42, 1185– 1208. (26) Kim, W. J.; Yang, S. M. Langmuir 2000, 16, 4761–4765. (27) Vasudevan, M.; Shen, A.; Khomami, B.; Sureshkumar, R. J. Rheol. 2008, 52, 527–550. (28) Rehage, H.; Hoffmann, H. J. Phys. Chem. 1988, 92, 4712–4719. (29) Boek, E. S.; Jusufi, A.; Lowen, H.; Maitland, G. C. J. Phys.: Condens. Matter 2002, 14, 9413–9430. (30) Kim, W. J.; Yang, S. M. AdV. Mater. 2001, 13, 1191–1195. (31) Cates, M. E. Macromolecules 1987, 20, 2289–2296. (32) Cates, M. E.; Fielding, S. M. AdV. Phys. 2006, 55, 799–879. (33) Shelley, J. C.; Shelley, M. Y. Curr. Opin. Colloid Interface Sci. 2000, 5, 101–110. (34) Huang, C. C.; Xu, H.; Ryckaert, J. P. J. Chem. Phys. 2006, 125, 094901. (35) Huang, C. C.; Xu, H.; Ryckaert, J. P. Europhys. Lett. 2008, 81, 58002. (36) Hayter, J. B.; Penfold, J. Colloid Polym. Sci. 1983, 261, 1022– 1030. (37) Dorshow, R. B.; Bunton, C. A.; Nicoli, D. F. J. Phys. Chem. 1983, 87, 1409–1416. (38) Johnson, S. B.; Drummond, C. J.; Scales, P. J.; Nishimura, S. Colloids Surf., A 1995, 103, 195–206. (39) Manohar, C.; Rao, U. R. K.; Valaulikar, B. S.; Lyer, R. M. J. Chem. Soc., Chem. Commun. 1986, 379–381. (40) Shikata, T.; Hirata, H.; Kotaka, T. Langmuir 1988, 4, 354–359. (41) Shikata, T.; Hirata, H.; Kotaka, T. J. Phys. Chem. 1990, 94, 3702– 2706. (42) Hartmann, V.; Cressely, R. Colloids Surf. 1997, 121, 151–162. (43) Olsson, U.; Soderman, O.; Guering, P. J. Phys. Chem. 1986, 90, 5223–5232. (44) Ulmius, J.; Wennerstroem, H.; Johansson, L. B. A.; Lindblom, G.; Gravsholt, S. J. Phys. Chem. 1979, 83, 2232–2236. (45) Wennerstrom, H.; Lindman, B. Phys. Rep. 1979, 52, 1–86. (46) Kreke, P. J.; Magid, L. J.; Gee, J. C. Langmuir 1996, 12, 699– 705. (47) Mohanty, S.; Davis, H. T.; McCormick, A. V. Langmuir 2001, 17, 7160–7171. (48) Maillet, J. B.; Lachet, V.; Coveney, P. V. Phys. Chem. Chem. Phys. 1999, 1, 5277–5290. (49) Marrink, S. J.; Tieleman, D. P.; Mark, A. E. J. Phys. Chem. B 2000, 104, 12165–12173. (50) Yakovlev, D. S.; Boek, E. S. Langmuir 2007, 23, 6588–6597. (51) Piotrovskaya, E. M.; Vanin, A. A.; Smirnova, N. A. Mol. Phys. 2006, 104, 3645–3651. (52) Boek, E. S.; den Otter, W. K.; Briels, W. J.; Iakovlev, D. Philos. Trans. R. Soc. London, Ser. A 2004, 362, 1625–1638. (53) Padding, J. T.; Boek, E. S.; Briels, W. J. J. Phys.: Condens. Matt. 2005, 17, S3347-S3353. (54) Boek, E. S.; Padding, J. T.; Anderson, V. J.; Briels, W. J.; Crawshaw, J. P. J. Non-Newtonian Fluid Mech. 2007, 146, 11–21. (55) Padding, J. T.; Boek, E. S.; Briels, W. J. J. Chem. Phys. 2008, 129, 074903. (56) Schuler, L. D.; Daura, X.; Gunsteren, W. F. V. J. Comput. Chem. 2001, 22, 1205–1218. (57) Shang, B. Z.; Wang, Z.; Larson, R. G. J. Phys. Chem. B 2008, 112, 2888–2900. (58) Schuettelkopf, A. W.; van Aalten, D. M. F. Acta Crystallogr. 2004, 60, 1355–1363. (59) Pal, S.; Bagchi, B.; Balasubramanian, S. J. Phys. Chem. B 2005, 109, 12879–12890. (60) Jorge, M. Langmuir 2008, 24, 5714–5725. (61) Song, Y. H.; Guallar, V.; Baker, N. A. BioChem. 2005, 44, 13425– 13438.

13710

J. Phys. Chem. B, Vol. 113, No. 42, 2009

(62) Fisicaro, E.; Ghiozzi, A.; Pelizzetti, E.; Viscardi, G.; Quagliotto, P. L. J. Colloid Interface Sci. 1996, 184, 147–154. (63) Briels, W. J.; Mulder, P.; den Otter, W. K. J. Phys.: Condens. Matter 2004, 16, S3965-S3974. (64) Berendsen, H. J. C.; Postma, J. P. M.; Gunsteren, W. F.; Hermans, J. Interaction Models for Water in Relation to Protein Hydration. In Intermolecular Forces; Pullman, B., Ed.; Reidel: Dordrecht, The Netherlands, 1981; pp 331-342. (65) der Spoel, D. V.; Lindahl, E.; Hess, B.; Groenhof, G.; Mark, A. E.; Berendsen, H. J. C. J. Comput. Chem. 2005, 26, 1701–1718. (66) Berendsen, H. J. C.; Postma, J. P. M.; van Gunsteren, W. F.; DiNola, A.; Haak, J. R. J. Chem. Phys. 1981, 81, 3684–3690. (67) Hess, B.; Bekker, H.; Berendsen, H. J. C.; Fraaije, J. G. E. M. J. Comput. Chem. 1997, 18, 1463–1472. (68) Ryckaert, J. P.; Bellemans, A. Faraday Discuss. Chem. Soc. 1978, 66, 95–106.

Wang and Larson (69) Miyamoto, S.; Kollman, P. A. J. Comput. Chem. 1992, 13, 952– 962. (70) Essmann, U.; Perera, L.; Berkowitz, M. L.; Darden, T.; Lee, H.; Pedersen, L. G. J. Chem. Phys. 1995, 103, 8577–8593. (71) Sepulveda, L.; Cortes, J. J. Phys. Chem. 1985, 89, 5322–5324. (72) Berr, S.; Jones, R. R. M.; Johnson, J. S., Jr. J. Phys. Chem. 1992, 96, 5611–5614. (73) Sammalkorpi, M.; Karttunen, M.; Haataja, M. J. Am. Chem. Soc. 2008, 130, 17977–17980. (74) Ruckenstein, E.; Beunen, J. A. Langmuir 1988, 4, 77–90. (75) Cuccovia, I. M.; da Silva, I. N.; Chaimovich, H. Langmuir 1997, 13, 647–652. (76) Mysels, K. J.; Princen, L. H. J. Phys. Chem. 1959, 63, 1696–1700. (77) Imae, T. J. Phys. Chem. 1990, 94, 5953–5959.

JP901576E