From Homogeneous Dispersion to MicellesA Molecular Dynamics

CNRS-ENSIC-INPL, 1 rue Grandville, BP 451, 54001 Nancy, France, and. Institut Universitaire de France, Maison des Universite´s,. 103 Boulevard Saint-...
0 downloads 0 Views 242KB Size
Langmuir 2005, 21, 5223-5229

5223

From Homogeneous Dispersion to MicellessA Molecular Dynamics Simulation on the Compromise of the Hydrophilic and Hydrophobic Effects of Sodium Dodecyl Sulfate in Aqueous Solution Jian Gao,†,‡ Wei Ge,*,† Guohua Hu,§,| and Jinghai Li*,† Multiphase Reaction Laboratory, Institute of Process Engineering, Chinese Academy of Sciences, Beijing 100080, People’s Republic of China, Graduate School of the Chinese Academy of Sciences, Beijing 100039, People’s Republic of China, Laboratory of Chemical Engineering Sciences, CNRS-ENSIC-INPL, 1 rue Grandville, BP 451, 54001 Nancy, France, and Institut Universitaire de France, Maison des Universite´ s, 103 Boulevard Saint-Michel, 75005 Paris Received November 23, 2004. In Final Form: February 1, 2005 The structural and functional diversity of surfactant systems has attracted simulation works in atomistic, coarse grain, and mesoscopic models (Bandyopadhyay, S.; et al. Langmuir 2000, 16, 942; Senapati, S.; et al. J. Phys. Chem. B 2003, 107, 12906; Maiti, P. K.; et al. Langmuir 2002, 18, 1908; Srinivas, G.; et al. J. Phys. Chem. B 2004, 108, 8153; Groot, R. D.; et al. J. Chem. Phys. 1999, 110, 9739; Rekvig, L.; et al. Langmuir 2003, 19, 8195). However, atomistic models have suffered from their tremendous computational cost and are, so far, not able to simulate the structural behaviors in sufficient spatio-temporal scales (Shelley, J. C.; Shelley, M. Y. Curr. Opin. Colloid Interface Sci. 2000, 5, 101). The other two approaches are not microscopic enough to describe the configurations of the surfactants that determine their behaviors (Shelley and Shelley). In this study, we propose to simplify atomistic models based on the observation that the compromise of the hydrophilic and hydrophobic effects (Li, J.; Kwauk, M. Chem. Eng. Sci. 2003, 58, 521-535) and molecular structures of surfactants are the dominant factors shaping their structures in the systems. With this simplification, we are able to simulate with moderate computing cost the whole process of micelle formation from an initially uniform dispersion of sodium dodecyl sulfate (SDS) in aqueous solution. The resulting micelle structures are different from those predicted by atomistic simulations that started with a predefined micelle configuration at the same surfactant concentrations. However, if we use their initial micelle configuration, micelle structures the same as theirs are obtained. Analyses show that our results are more realistic and that the results of the atomistic simulations suffer from artificial initial conditions. Therefore, our model may serve as a reasonable simplification of atomistic models in terms of the general structure of micelles.

1. Introduction Above the critical micelle concentration (CMC), surfactant molecules in aqueous solution aggregate to form a wide variety of assemblies, such as spherical or rodlike micelles, or more complex morphologies under different conditions.8 As a result of the amphiphilic nature of surfactant molecules, the aggregation process is sponta* To whom correspondence should be addressed. Tel.: +86-108261 6050. Fax: +86-10-6255 8065. E-mail: [email protected] (W.G.); [email protected] (J.L.). † Institute of Process Engineering, Chinese Academy of Sciences. ‡ Graduate School of the Chinese Academy of Sciences. § CNRS-ENSIC-INPL. | Institut Universitaire de France. (1) Bandyopadhyay, S.; Tarek, M.; Lynch, M. L.; Klein, M. L. Langmuir 2000, 16, 942-946. (2) Senapati, S.; Berkowitz, M. L. J. Phys. Chem. B 2003, 107, 1290612916. (3) Maiti, P. K.; Lansac, Y.; Glaser, M. A.; Clark, N. A. Langmuir 2002, 18, 1908-1918. (4) Srinivas, G.; Shelley, J. C.; Nielsen, S. O.; Discher, D. E.; Klein, M. L. J. Phys. Chem. B 2004, 108, 8153-8160. (5) Groot, R. D.; Madden, T. J.; Tildesley, D. J. J. Chem. Phys. 1999, 110, 9739-9749. (6) Rekvig, L.; Kranenburg, M.; Vreede, J.; Hafskjold, B.; Smit, B. Langmuir 2003, 19, 8195-8205. (7) Shelley, J. C.; Shelley, M. Y. Curr. Opin. Colloid Interface Sci. 2000, 5, 101-110. (8) Gelbart, W. M.; Ben-Shaul, A.; Roux, D. Micelles, Membranes, Microemulsions, and Monolayers; Springer-Verlag: New York, 1994.

neous with a reduction of free energy. As the phenomena appear on tiny spatio-temporal scales (in nanometers and nanoseconds), direct experimental investigations have been fairly difficult, if not impossible.9 Microscopic simulations may serve as a powerful tool for the explorations in this area with their flexibility in isolating specific acting factors from the complexity with idealized models and exhaustive data analysis in a truly nonintrusive manner. Microscopic simulations could be different in their levels of description. In the most fundamental atomistic models, 9 all atoms of the system are explicitly represented, and the Lennard-Jones (LJ) potential is often used to estimate the potential energy of the ensemble. These models are capable of obtaining detailed information on micelle configuration and morphologies. However, because of their tremendous computational cost, they usually have to start with a single micelle of a certain size and shape in the solution; that is, only the short-time dynamics of relaxation to equilibrium can be studied. As a matter of fact, it is not even known whether the resulting structures correspond to equilibrium conditions or not.10-12 (9) Marrink, S. J.; Tieleman, D. P.; Mark, A. E. J. Phys. Chem. B 2000, 104, 12165-12173. (10) Tieleman, D. P.; van der Spoel, D.; Berendsen, H. J. C. J. Phys. Chem. B 2000, 104, 6380-6388. (11) Bogusz, S.; venable, R. M.; Pastor, R. W. J. Phys. Chem. B 2000, 104, 5462-5470. (12) Bast, T.; Hentschke, R. J. Phys. Chem. 1996, 100, 12162-12171.

10.1021/la047121n CCC: $30.25 © 2005 American Chemical Society Published on Web 04/21/2005

5224

Langmuir, Vol. 21, No. 11, 2005

Gao et al.

Scheme 1. Simplified Molecular Structure of SDS

Table 1. LJ Parameters Used for the Simulations  (kcal/mol) σ (Å)

water

CHn

O

S

0.155 3.166

0.118 3.905

0.200 3.150

0.250 3.550

Table 2. Parameters c, s, and d Used for the Simulations

models13

Coarse grain and dissipative particle dynamics14 models offer much larger scales for simulating spontaneous micelle formation. However, in these models several atoms or a section of chained radicals is represented by a single computational particle. It is, thus, impossible to establish relationships between molecular surfactant structure and its amphiphilic properties. Comparison of different models was made by Rajagopalan15 and Shelley and Shelley.7 In this paper, we focus on one of the most commonly used surfactants, the anionic detergent sodium dodecyl sulfate (SDS), and try to develop a hybrid of these models that can simulate both the formation process and the structure details of the micelles formed. 2. Simulation Methodology Some authors have pointed out that many surfactant systems follow similar patterns.16 This is a general motivation and justification for our simplification of atomistic models. It mainly consists of the following aspects: Simplification of the Molecular Structure (See Scheme 1). The headgroup includes one sulfur atom and four oxygen atoms, which is explicitly represented by five particles. For the alkane chain, however, each CHn radical in SDS is modeled as a single particle Cx. Nevertheless, the framework of the molecule is well preserved, especially for factors affecting the hydrophile-lipophile balance of the amphiphilic molecule, such as the relative position between the chained CHn particles, and the structure of the hydrophilic SO4 radical. Each water molecule is also simulated as a single isotropic particle. Simplification of the Interactions between the Model Particles. We concentrate on the dominant interactions shaping the hydrophilic and hydrophobic characteristics of the surfactant, which have been identified as the chemical bonds supporting the framework of the molecule and the interactions between the surfactant and the water molecules. For water molecules, the extended simple point charge (SPC/E) potential model of Berendsen et al.17 has a single LJ site for the oxygen atom that carries a charge of -0.8476e, and two hydrogen atoms are also simply represented by point charges, each carrying a charge of +0.4238e. The σ (size parameter) and  (energy parameter) of the hydrogen atoms in the model are 0. So, an oxygen atom can represent a water molecule, if one does not consider electrostatics. In our model, a water molecule is represented by a single particle with LJ potential.18 The LJ parameters of water molecules are the same as those of the oxygen atom of the SPC/E potential model, but the Coulombic charge on the oxygen atom is ignored. (13) Shelley, J. C.; Shelley, M. Y.; Reeder, R. C.; Bandyopadhyay, S.; Klein, M. L. J. Phys. Chem. B 2001, 105, 4464-4470. (14) Whittle, M.; Dickinson, E. J. Colloid Interface Sci. 2001, 242, 106-109. (15) Rajagopalan, R. Curr. Opin. Colloid Interface Sci. 2001, 6, 357365. (16) Esselink, K.; Hilbers, P. A. J.; van Os, N. M.; Smit, B.; Karaborni, S. Colloids Surf., A 1994, 91, 155-167. (17) Berendsen, H. J. C.; Grigera, J. R.; Straatsma, T. P. J. Chem. Phys. 1987, 91, 6269-6271. (18) Smit, B.; Hilbers, P. A. J.; Esselink, K.; Rupert, L. A. M.; van Os, N. M.; Schlijper, A. G. Nature 1990, 348, 624-625.

O

c

s

d

0.9

0.8

1

S

c

s

d

0.8

0.8

1

O-S

c

s

d

1.1

0.9

1

Table 3. Bond Lengths and Angles of SDS Molecules bond

bond length (Å)

bond

bond angle (deg)

C-C C1-O4 O4-S S-O1-3

1.53 1.42 1.58 1.46

C-C-C C-C-O C-O-S O(ester)-S-O O-S-O

111.0 109.5 112.6 102.6 115.4

For the SDS molecule, the electrostatic forces are considered implicitly by adjusting the potential between the ionic SO4 radicals upon neglecting the long-range tails. The simplification should be reasonable for the simulations of primary clustering structures such as micelles where surfactant molecules are closely packed. In molecular dynamics (MD) simulations, the electrostatic interaction between two particles is, in general, calculated by qiqj/rij, where qi represents the Coulombic charge on the atom i. That is,

U(r bij) )

[( ) ( ) ]

qiqj σij + 4ij rij rij

12

-

σij rij

6

(1)

In our LJ potential function, intermolecular and intramolecular interactions (excluding bonded sites) are modeled using a cut and shifted LJ potential:

bij) ) Uij(r

{

4dij 0

[(

σij

s(rij -r0)+r0

) ( 12

-c

σij

s(rij -r0)+r0

)] 6

cut - Uij(rcut ij ) rij e rij

rij > rcut ij

(2)

where rij is the distance separating two particles and r0 is the value of rij when Uij is 0. The values of those parameters are chosen so that SDS is best approximated.19 In particular, we believe that electrostatic interactions also contribute to the hydrophobic and hydrophilic effects,3 and, therefore, the parameters c and s are adjusted to mimic the electrostatic interactions of O-O, S-S, and O-S, qualitatively. The corresponding LJ parameters are compiled in Tables 1 and 2. The bond angle bending potential is described by a harmonic function of the form U(θ) ) 1/2K(cos θ - cos θ0)2, where K is a constant and θ0 is the balance bond angle,20 and the bond bending force constant is K ) 520 kJ/mol.21 The bond length is kept by the SHAKE algorithm.22 Table 3 shows the values of the bond lengths23 and angles in a SDS molecule. The geometry adopted for the headgroup agrees well with X-ray dif(19) Schweighofer, K. J.; Essmann, U.; Berkowitz, M. J. Phys. Chem. B 1997, 101, 3793-3799. (20) Salaniwal, S.; Cui, S. T.; Cochran, H. D.; Cummings, P. T. Langmuir 2001, 17, 1773-1783. (21) Brown, D.; Clarke, J. H. R. J. Chem. Phys. 1994, 100, 16841692. (22) Rapaport, D. C. The Art of Molecular Dynamics Simulation; Cambridge University Press: Cambridge, 1995. (23) Shelley, J.; Watanabe, K.; Klein, M. L. Int. J. Quantum Chem., Quantum Biol. Symp. 1990, 17, 103-117.

MD Simulation of SDS in Aqueous Solution

Langmuir, Vol. 21, No. 11, 2005 5225

Table 4. Key Parameters for Three Previous Simulations Shelley et al.23 SDS molecule number water molecule number concentration temperature time software Imax/Imin radius of gyration, Rg

42 1901 1.1 M 298 K 182 ps 1.13

MacKerell31 60 4398 298 K 120 ps CHARMM 1.02 16.02 Å

Bruce et al.28 60 7579 0.4 M 300 K 5 ns AMBER 1.05 16.2 Å

fraction data.24 Furthermore, as in many coarse grain models, the torsional forces between the CHn radicals are neglected. Two different cutoff radii are used in the simulations. 1/6 The first one, rcut ij ) 2 σij, results in a purely repulsive potential, which is used to calculate water-CHn, CHn-O, and CHn-S interactions. The second one is rcut ij ) 2.5σij, that is, the potential has a short-range attractive part in addition, which is used to calculate the interactions of water-water, water-O, water-S, and O-S. The use of such a modified LJ potential is a significant simplification and as a consequence will hardly allow exact quantitative studies, but rather it will qualitatively mimic the true properties of surfactant systems. To test our water model, we performed a MD simulation of pure water. For comparison, our simulation conditions are set to be those of ref 25. The system includes 216 water molecules, temperature T ) 34.3 °C, and density F ) 1.06 g/cm3. The self-diffusion coefficient measured in our simulation was D ) 2.55 ( 0.33 × 10-5 cm2/s, whereas in ref 25, D ) 4.2 × 10-5 cm2/s, and the experimental value cited therein was D ) 2.85 ( 0.15 × 10-5 cm2/s. In view of the rapid variation of the experimental D with temperature, both our simulation values and those of ref 25 should be considered favorable. 3. Results and Discussion We performed MD simulations for two systems: one includes 60 SDS molecules and 2730 water molecules, and the other has 240 SDS molecules and 10 920 water molecules. The simulations were performed in the NVT ensemble (the number of particles, the volume, and the temperature are constant), and the conventional periodic boundary conditions in MD simulations are applied. The equations of motion are integrated via the Verlet algorithm26,27 with a time step of 3.3 fs. The constraint methods are used to keep a constant temperature22 of 298 K. The code used in the simulations was programmed and parallelized by ourselves. For the system of 60 SDS molecules, the simulation for 2 ns took about 37 h on 4 PIII 1.13G CPUs. For the system of 240 SDS molecules, the CPU time is 152 h. Three “large-scale” investigations into the behavior of the ionic surfactant SDS using atomistic models have been reported, as summarized in Table 4. That of Bruce et al.28 is, to our knowledge, the largest scale simulation of its kind so far. At the start of this simulation, the micelle is so constructed that the innermost methyl group in each of the hydrophobic tails is placed at the apex of a ball with (24) Jarvis, J. A. J. Acta Crystallogr. 1953, 6, 327-330. (25) Rahman, A.; Stillinger, F. H. J. Chem. Phys. 1971, 55, 33363359. (26) Allen, M. P.; Tildesley, D. J. Computer Simulation of Liquids; Clarendon Press: Oxford, 1990. (27) Rychaert, J. P.; Ciccotti, G.; Berendsen, H. J. C. J. Comput. Phys. 1977, 23, 327-341. (28) Bruce, C. D.; Berkowitz, M. L.; Perera, L.; Forbes, M. D. E. J. Phys. Chem. B 2002, 106, 3788-3793.

Figure 1. Breakdown of the potential energy for the system with 60 SDS molecules.

a radius of 3.5 Å28 and the remainder of the monomer extending outward. The number of SDS molecules necessary for forming an aggregate was reported to be in the range of 60-70 at room temperature.29 Thus, it is considered appropriate to choose 60 SDS molecules for forming a micelle. This number is both small enough to be computationally feasible and large enough to mimic all of the physical properties of a micelle. All those studies reported flexible hydrocarbon tails, stable micelle structure, and low water penetration into the micelle. They also compared simulated results with experiments. 3.1. Micelle Shape. We first simulated a single SDS micelle and compared it with their results. However, to simulate the micelle formation process and to simplify the initialization process, we started with the SDS molecules distributed uniformly in an array. The final result was an elliptical micelle with the ratios of the principle moments of inertia (I1/I3, I2/I3, I1/I2) being 1.52 ( 0.02, 1.29 ( 0.03, and 1.18 ( 0.05, respectively. This is in contrast to the values reported in Table 4 which are much closer to unity. This notable difference can be ascribed either to the simplifications in our model or to the difference in the initial conditions. We noticed in particular that the SDS concentrations in the previous simulations were considerably higher than its CMC, which is only about 0.008 M at room temperature. According to Gelbart et al.,8 the shape of micelles should be rodlike if the SDS concentrations are high enough, for example, one to two orders of magnitude above its CMC. We also compare the potential energy between water and the headgroups of the SDS molecules and that between water and the alkanes. As shown in Figure 1, the latter is considerably higher, demonstrating that the amphiphilic nature of the SDS molecules has been well presevered. Therefore, the results given by our simulation approach seem to be more reasonable than those of the previous simulations in terms of the micelle shapes at SDS concentrations well above its CMC. However, the validity of the previous simulations can be justified for studying the statistical properties of existing spherical micelles formed with SDS concentrations just above the CMC. In that case, only the surrounding water is still influential to the SDS molecules, at least for a short period of time compared with the lifetime of the micelle. In other words, the SDS concentration is no longer a critical parameter (29) Chen, J. M.; Su, T. M.; Mou, C. Y. J. Phys. Chem. 1986, 90, 2418-2421.

5226

Langmuir, Vol. 21, No. 11, 2005

Figure 2. Variation of the ratios of the principle moments of inertia of the SDS micelle simulated.

Gao et al.

Figure 3. Spherical micelle of 60 SDS molecules in aqueous solution with 2730 water molecules simulated. Coloring scheme: CHn in red, O or S in blue. Water molecules are not shown for clarity.

Scheme 2. Relaxation of a Predefined Micelle to a Spherical Initial Configuration

for the dynamics of an existent micelle, and the micelle structure can be reasonably studied by a sample volume cut from the solution with very few adjacent water molecules. Although experimentally impossible, the stability of the micelle solution suggests that we can start with a state sufficiently close to the target micelle, and a reasonable simulation should lead to a stable spherical micelle with correct statistical properties, as has been done in the simulations in Table 4. To compare the micelle structure with other reports, we “forged” the elliptical structure we have obtained to a spherical initial conformation. We moved every SDS molecule radially to where their headgroups have equal distance to the micelle center of mass (MC). To avoid overlapping with the water, the computational domain was also expanded slightly and the water molecules were mapped linearly to their new position. Then, the system was relaxed and the computational domain was gradually shrunk to the original one (as shown in Scheme 2). This process was repeated several times with decreasing amplitudes. Finally, we obtained a spherical micelle with (Ii/Ij)max ) 1.04. Thereafter, the simulation lasted 2 ns and, after a slight rebond at first, the geometry of the micelle is maintained at I1/I3 ) 1.09 ( 0.02, I2/I3 ) 1.05 ( 0.01, and I1/I2 ) 1.04 ( 0.03. Figure 2 shows the variation of Ii/Ij in this period of time, and Figure 3 is the snapshot of the micelle at the end. In Figure 4, the total potential energy of the elliptical micelle is -0.401 ( 0.005 kcal/mol in the last 1 ns, and that of the spherical micelle is -0.362 ( 0.004 kcal/mol in the last 1 ns. The potential energy of the elliptical micelle is much lower than that of the

Figure 4. Variation of total potential energies in systems forming elliptical and spherical micelles, respectively.

spherical one, showing that the former is more stable. In what follows, we will demonstrate that the more stable elliptical micelle structure is also more reasonable in other statistical aspects. 3.2. Micelle Size. One of the primary characteristics of micelle structure is its size. Our simulated radius of gyration is 15.5 ( 0.1 Å. This value is in good agreement with the usually referenced experimental value of 15.4 Å for a lithium dodecyl sulfate (LDS) micelle.30 In comparison, the value reported by MacKerrell31 is 16.02 Å. For the positions of the headgroups, as designated by the sulfur atom, our root-mean-square distance from the MC is 18.8 ( 0.2 Å on average, whereas that of MacKerrell31 is 19.70 Å. This is to be compared with the experimental value of Bendedouch et al.30 which is 18.9 Å for LDS using small angle neutron scattering. We think that the size of the LDS micelle should be almost the same as the SDS micelle, because the counterions are dissociated in aqueous solution; that is, the micelles do not include any Na+ or Li+, and the rest of a LDS molecule is identical to that of a SDS molecule. The effect of the cations in the solution on the size of the micelle is not crucial. In summary, the overall micelle geometry predicted by our simulation approach is in satisfactory agreement with experimental and comparable to atomistic simulation results. 3.3. Density Distribution of the Micelle Group. The micelle structure is now analyzed in terms of the density (30) Bendedouch, D.; Chen, S.-H.; Koehler, W. C. J. Phys. Chem. 1983, 87, 153-159. (31) MacKerell, J. A. D. J. Phys. Chem. 1995, 99, 1846-1855.

MD Simulation of SDS in Aqueous Solution

Figure 5. Density distribution of the group-to-MC distance: CHn (0), S (b), and water (2).

Figure 6. Probability distribution of the group-to-MC distance: S (9), C1 ([), C6 (2), and C12 (O).

of selected groups as a function of the distance from the MC. As shown in Figure 5, the carbon density decreases rapidly between 15 and 20 Å, while the sulfur density decreases between 21 and 25 Å. These density profiles are in good agreement with the experimental value of 16.7 Å for the alkyl radius and of 22.3 Å for the total radius from X-ray scattering.32 Comparing the carbon density profiles with the water one, we find that the interior of the micelle is completely void of solvent (water), which is again in good agreement with experiments30,33,34 and atomistic simulations23,28,31 in Table 4. 3.4. Probability Distribution of the Micelle Group. We studied more detailed structure of the micelle through the probability distribution of the C1, C6, C12, and S with respect to the micelle MC. Figure 6 shows the profiles of the four groups in the micelle. The probability distributions of the CHn in the hydrocarbon tails tend to broaden as the distance from the headgroup increases. For the C1 group the distribution is similar to that of the sulfur atom in the headgroup. An increase in the width of distribution occurs on the C6 group, while the C12 distribution even has a tail in the region of the sulfur distribution. This trend was (32) Itri, R.; Amaral, L. Q. J. Phys. Chem. 1991, 95, 423-427. (33) Jones, R. R. M.; Maldanado, R.; Szajdzinska-Pietek, E.; Kevan, L. J. Phys. Chem. 1986, 90, 1126-1129. (34) McManus, H. J. D.; Kang, Y. S.; Kevan, L. J. Phys. Chem. 1992, 96, 5622-5628.

Langmuir, Vol. 21, No. 11, 2005 5227

Figure 7. Hydration numbers of the S and CHn groups in the hydrocarbon tails, which are counted for all water molecules within 3.5 Å of a group.31

also observed in the studies of Shelley et al.23 and MacKerrell.31 In general, the hydrophobic tails of the surfactant molecules are curled rather than straight and not perfectly pointing to the MC of the micelle.28 3.5. Micelle-Water Interactions. The interactions between the micelle and the water can be characterized in Figure 5 and and also in Figure 6. To analyze the interactions between the alkyl chains and the water, the hydration numbers of the individual group are determined. As shown in Figure 7, it decreases rapidly from C1 to C3. This is similar to the report of Shelley et al.23 and MacKerell.31 A small increase occurs with the C12 methyl group, which is consistent with the broad probability distribution of the C12 group in Figure 6 and the report of MacKerell,31 supporting previous theoretical descriptions of micelle structure.35 3.6. Spontaneous Formation of SDS Micelles. The foregoing discussion suggests that our model yields reasonable micelle structures very close to those given by atomistic models under similar conditions. This justifies its use for exploring the spontaneous micelle formation process in these systems. The second system contains 10 920 water molecules and 240 SDS molecules. The concentration of SDS is the same as that of Shelley et al.23 A major difference between their work and ours is that in the latter the micelle is not built in advance. Instead, the SDS molecules are evenly distributed in a uniform matrix. The variation of the potential energy in the simulation was monitored and is shown in Figure 8. We can find that the energy drops sharply at the beginning and then fluctuates about a constant mean value. To test whether the system had reached the equilibrium, we re-initialized the velocities of every particle of the system under the same temperature at t ) 2 ns and then simulated the system for another 2 ns. No visible changes in the conformation as well as in the potential energy have been observed. Two rodlike micelles are formed and stabilized at last, and the aggregation numbers are 119 and 102, respectively. These two values are much higher than that at the CMC. This is also in accordance with classical experimental results36 but is certainly not accessible by previous simulations without enough SDS molecules. A snapshot of the micelles is shown in Figure 9. We compared our simulation results (35) Dill, K. A.; Koppel, D. E.; Cantor, R. S.; Dill, J. D.; Bendedouch, D.; Chen, S.-H. Nature 1984, 309, 42-45. (36) Oh, S. G.; Shah, D. O. J. Am. Oil Chem. Soc. 1993, 70, 673-678.

5228

Langmuir, Vol. 21, No. 11, 2005

Figure 8. Variation of the potential energy in the system with 240 SDS molecules.

Gao et al.

Figure 11. Probability distribution of the group-to-MC distance: S (9), C1 ([), C6 (2), and C12 (O).

Figure 9. Two rodlike micelles of SDS. The system consists of 10 920 water molecules and 240 SDS molecules. The coloring scheme is the same as for Figure 3. Figure 12. Hydration numbers of the S and CHn groups.

Figure 10. Density distribution of the group-to-MC distance: CHn (0), S (b), and water (2).

with the phase diagram of SDS in ref 37 and found that our system located correctly in the micellar region, showing, at least, a qualitative agreement with the experiments therein. Figure 10 shows the density distribution of group CHn, S, and water in a micelle in Figure 9. It is similar to that of a spherical micelle. There are no water molecules inside (37) Ke´kicheff, P.; Grabielle-Madelmont, C.; Ollivon, M. J. Colloid Interface Sci. 1989, 131, 113-132.

the micelle, and hydrocarbon chains form the core of the micelle. The cross-sectional radius of the rodlike micelle is about 14 Å, less than that of the spherical micelle. The density of the hydrocarbon chain is larger than that of the spherical micelle. So the interior of the rodlike micelle is more tightly packed than that of the spherical one. Figure 11 shows the probability distribution of groups C1, C6, C12, and S for the rodlike micelle. It is again similar to that for the spherical one. The group C12 has the broadest distribution in the interior of the micelle and extends to the range of group S. The groups C1 and S have broader distributions inside the rodlike micelle than inside the spherical one. The interactions between the rodlike micelle and the water can be analyzed from the group hydration numbers. Figure 12 shows the hydration numbers of different groups. They are similar to those of the spherical micelle. However, the hydration number of each group in the rodlike micelle is slightly larger than that of the spherical one. This shows that the rodlike micelle is more stable than the spherical one, because a higher hydration number corresponds to more contact between the hydrophilic head of the SDS and water, leading to a more significant decrease in free energy of the system. 4. Conclusion We have simplified atomistic models of SDS by maintaining the framework of the molecular structure and

MD Simulation of SDS in Aqueous Solution

atomistic interactions influential to the amphiphilic characters of the molecule while smoothing out other details. The resulting model has produced reasonable statistical properties of predefined spherical micelles in aqueous solutions, which are almost the same as those of detailed atomistic models. Simulation of spontaneous micelle formation in these systems can now be carried out with moderate computation and with little loss in essential structural information. Results show that the shape of SDS micelles under concentrations much higher than its CMC is not spherical but rodlike, which is in agreement

Langmuir, Vol. 21, No. 11, 2005 5229

with experimental observations but has never been reproduced at an atomistic level. Acknowledgment. We thank Prof. Huizhou Liu and Prof. Guanghui Ma (Institute of Process Engineering, Chinese Academy of Sciences) for many illuminating discussions. This work is supported by the National Natural Science Foundation of China under Grants 20336040 and 20221603 and the Chinese Academy of Sciences under its Outstanding Overseas Research Team Project. LA047121N