Simulated Phase Behavior of Model Surfactant Solutions - Langmuir

The simulated systems comprise up to 100 000 sites with a maximum of 12 sites per surfactant molecule. We examine the effect of different conditions ...
0 downloads 0 Views 397KB Size
Download PDF
1626

Langmuir 2000, 16, 1626-1633

Simulated Phase Behavior of Model Surfactant Solutions B. Fodi and R. Hentschke* Max-Planck-Institut fu¨ r Polymerforschung, Postfach 3148, 55021 Mainz, Germany Received July 1, 1999. In Final Form: October 6, 1999 We performed molecular dynamics simulations in an attempt to model the phase behavior of a nonionic surfactant. We employed a coarse-grained model, parametrized to yield the phase behavior of the corresponding physical system, instead of using a realistic model with its severe spatial and temporal limitations. Within that model our surfactants are composed of chains of Lennard-Jones sites connected by harmonic springs, and the solvent is modeled as a Lennard-Jones fluid. The simulated systems comprise up to 100 000 sites with a maximum of 12 sites per surfactant molecule. We examine the effect of different conditions (temperature and concentration) and different model parameters (head/tail size and ratio, strength of interaction between the different sites). Our simulations show that some of the real system’s characteristics can be mirrored by the coarse-grained model, whereas some aspects cannot.

1. Introduction The amphiphilic character of surfactant molecules in aqueous solutions or in oil/water mixtures leads to a variety of self-assembled aggregates such as spherical or rodlike micelles, bilayers, or more complex morphologies. The tendency of amphiphiles to self-organize so that the tails are shielded from the water finds its fullest expression in the lyotropic phases which occur when the concentration of amphiphiles is large. Examples include the cubic phase (a cubic packing of normal or inverted micelles), the hexagonal phase (a close-packed array of long cylindrical micelles), and the lamellar phase (which consists of sheets of amphiphiles). The considerable industrial importance of these systems motivates research activities on both the experimental1 and the theoretical2,3 side (a general overview can be found in ref 4). However, our theoretical understanding of real surfactant solutionssespecially at higher concentrationssis still fairly limited. Different theoretical approaches utilizing both analytical and simulation techniques are used to examine surfactant structures. The appropriate method depends on the time and length scales underlying the physical property of interest. One can divide the models into three groups,5 namely microscopic models, time-dependent Ginzburg-Landau equations, and membrane theories. Most investigations have focused on the low concentration regime in the vicinity of the critical micelle concentration (cmc). In this context numerous studies have been devoted to predicting micelle properties from the molecular composition. Molecular modeling methods, molecular dynamics (MD) and Monte Carlo (MC) simulations, as well as self-consistent field calculations have been applied.6-41 Recent review articles by Larson,42 and Karaborni and Smit3 nicely summarize the application of * Permanent address for correspondence: Bergische Universita¨tGesamthochschule, FB 8 Physik, Postfach 100127, D-42097 Wuppertal, Germany. (1) Fletcher, P. Curr. Opin. Colloid Interface Sci. 1996, 1, 101-106. (2) Gelbart, W. M.; Ben-Shaul, A. Micelles, Membranes, Microemulsions and Mololayers; Springer-Verlag: New York, 1994. (3) Karaborni, S.; Smit, B. Curr. Opin. Colloid Interface Sci. 1996, 1, 411-415. (4) Jo¨nsson, B.; Lindman, B.; Holmberg, K.; Kronberg, B. Surfactants and Polymers in Aqueous Solution; Wiley: New York, 1998. (5) Gompper, G.; Schick, M. In Phase Transitions and Critical Phenomena; Domb, C., Lebowitz, J. L., Eds.; Academic Press: London, 1994; Vol. 16. (6) Tanford, C. The hydrophobic effect. Formation of micelles and biological membranes; Wiley: New York, 1980.

computer simulation to the field of self-assembly, not only at the low concentration regime. (7) Nagarajan, R.; Ruckenstein, E. Adv. Colloid Interface Sci. 1986, 26, 205. (8) Puvvada, S.; Blankschtein, D. J. Chem. Phys. 1990, 92, 3710. (9) Marcelija, S. S. Biochim. Biophys. Acta 1974, 367, 165. (10) Dill, K. A. J. Phys. Chem. 1982, 86, 1498. (11) Gruen, D. W. R. J. Phys. Chem. 1985, 89, 146. (12) Ben-Shaul, A.; Szleifer, I.; Gelbart, W. M. J. Chem. Phys. 1985, 83, 3597. (13) Ben-Shaul, A.; Szleifer, I.; Gelbart, W. M. J. Chem. Phys. 1986, 85, 5345. (14) Leermakers, F. A. M.; Scheutjens, J. M. H. M. J. Chem. Phys. 1989, 93, 7417. (15) Leermakers, F. A. M.; Scheutjens, J. M. H. M. J. Colloid Interface Sci. 1990, 136, 231. (16) Leermakers, F. A. M.; van der Schoot, P. P. A. M.; Scheutjens, J. M. H. M.; Lyklema, J. Surfactants in Solution; Plenum: New York, 1990; Vol. 7. (17) Bo¨hmer, M. R.; Koopal, L. K. Langmuir 1990, 6, 1478. (18) Bo¨hmer, M. R.; Lyklema, J.; Koopal, L. K. J. Phys. Chem 1991, 89, 9569. (19) Linse, P. J. Phys. Chem. 1993, 97, 13896. (20) Linse, P. Macromolecules 1993, 26, 4437. (21) Linse, P. Macromolecules 1994, 27, 2685. (22) Linse, P. Macromolecules 1994, 27, 6404. (23) Jo¨nsson, B.; Eldholm, O.; Teleman, O. J. Chem. Phys 1986, 85, 2259. (24) K. Watanabe.; Ferrario, M.; Klein, M. L. J. Phys. Chem. 1988, 92, 819. (25) Wendoloski, J. J.; Kimatian, S. J.; Schutt, C. E.; Salemme, F. R. Science 1989, 243, 636. (26) Shelley, J. C.; Watanabe, K.; Klein, M. L. Int. J. Quantum Chem. Quantum Biol. Symp. 1990, 17, 103. (27) Shelley, J. C.; Sprik, M.; Klein, M. L. Langmuir 1993, 9, 916. (28) Laaksonen, L.; Rosenholm, J. B. Chem. Phys. Lett. 1993, 216, 429. (29) Bo¨cker, J.; Brickmann, J.; Bopp, P. J. J. Phys. Chem. 1994, 98, 712. (30) Smit, B.; Hilbers, P. A. J.; Esselink, K.; Rupert, L. A. M.; van Os, N. M. Nature 1990, 348, 624. (31) Smit, B.; Hilbers, P. A. J.; Esselink, K.; Rupert, L. A. M.; van Os, N. M. J. Phys. Chem. 1991, 95, 6361. (32) Leermakers, F. A. M.; Lyklema, J. Colloids Surf. 1992, 67, 239. (33) Haan, S. W.; Pratt, L. R. Chem. Phys. Lett. 1981, 79, 436. (34) Owenson, B.; Pratt, L. R. J. Phys. Chem. 1984, 88, 2905. (35) Pratt, L. R.; Owenson, B.; Sun, S. Adv. Colloid Interface Sci. 1986, 26, 69. (36) Larson, R. G. J. Chem. Phys. 1992, 96, 7904-7918. (37) Rodrigues, K.; Mattice, W. L. J. Chem. Phys. 1991, 95, 5341. (38) Wang, Y.; Mattice, W. L.; Napper, D. H. Langmuir 1993, 9, 66. (39) Haliloglu, T.; Mattice, W. L. Chem. Eng. Res. 1994, 49, 2851. (40) Rector, D. R.; van Swol, F.; Henderson, J. R. Mol. Phys. 1994, 82, 1009-1031. (41) Adriani, P.; Wang, Y.; Mattice, W. L. J. Chem. Phys. 1994, 100, 7718. (42) Larson, R. G. Curr. Opin. Colloid Interface Sci. 1997, 2, 361364.

10.1021/la990862h CCC: $19.00 © 2000 American Chemical Society Published on Web 12/17/1999

Phase Behavior of Model Surfactant Solutions

Langmuir, Vol. 16, No. 4, 2000 1627

In this work we confine ourselves to molecular simulations of binary surfactant/water mixtures. However, instead of using a realistic model that enables the observation of processes only on the nanosecond time scale, we employ a coarse-grained model in an attempt to access time scales more relevant to the study of surfactant phase behavior. Within our model surfactants are composed of chains of head and tail Lennard-Jones (LJ) sites connected by harmonic springs. In addition, the solvent is modeled by a single LJ site. Similar models that contain only the basic features of an amphiphilic system were recently used in MC simulations by Larson36,43-45 and MD calculations by Smit et al.30,31,46 The intention of our MD study is the investigation of “off lattice” self-organization at higher concentrations. Furthermore we want to explore whether the present simple model reproduces general features of the experimentally observed behavior. The “real” systems which are underlying our model are poly(oxyethylene) surfactants. This group of nonionic surfactants is widely used as emulsifying agents and detergents and still is of great scientific interest.47-55 Like ionic systems, poly(oxyethylene) surfactants form micelles at small concentrations as well as complex mesophases at higher concentrations. But, because this is only a crude model which we expect to reproduce only the general features of amphiphilic behavior, other single-tailed surfactants such as anionic, cationic, or zwitterionic as well as diblock-copolymer could also be devised. We examine the effect of different thermodynamic conditions, especially the concentration dependence and the influence of the model parameters, i.e., the head/tail size, their corresponding ratio, and the strength of the interactions between the different sites. The resulting morphologies are discussed in terms of local curvature measured via a simple algorithm, which appears to be a useful tool for analyzing the relation between the different model architectures and the observed structures. 2. Simulation Methodology In our model we limit ourselves to the basic features of a (linear) amphiphilic molecule with its two building blocks: hydrophilic units (head sites, H) and hydrophobic units (tail sites, T). This means that a model surfactant consists of a certain number of H and T sites, joined via harmonic potentials. Furthermore, a third LJ-type (W) is used to represent water molecules. The intermolecular and intramolecular interactions (excluding bonded sites) are modeled using a cut and shifted LJ potential,

Uij(r bij) )

{

4 0

[( ) ( ) ] σij rij

12

-

σij rij

6

cut - Uij(rcut ij ), rij e rij

(1)

rij > rcut ij

where rij is the separation between two particles, and  and σ are the well depth and the size parameter of the LJ potential, (43) Larson, R. G. J. Chem. Phys. 1989, 91, 2479. (44) Larson, R. G. Chem. Eng. Sci. 1994, 49, 2833-2850. (45) Larson, R. G. J. Phys. II 1996, 6, 1441-1463. (46) Smit, B. Phys. Rev. A 1988, 37, 3431. (47) Laughlin, R. G. Aqueous Phase Behaviour of Surfactants; Academic Press: New York, 1994. (48) Clerd, M.; Levelut, A. M.; Sadoc, J. F. J. Phys. II 1991, 1, 12631276. (49) Kahlweit, M. Ber. Bunsen-Ges. Phys. Chem. 1994, 98, 490-497. (50) Schubert, K.-V.; Strey, R.; Kline, S. R.; Kaler, E. W. J. Chem. Phys. 1994, 101, 5343-5355. (51) Funari, S. S.; Holmes, M. C. J. Phys. Chem. 1994, 98, 30153023. (52) Funari, S. S.; Rapp, G. J. Phys. Chem. B 1997, 101, 732-739. (53) Sallen, L.; Oswald, P.; Ge`minard, J. C. J. Phys. II 1995, 5, 937961. (54) Strey, S. Ber. Bunsen-Ges. Phys. Chem. 1996, 100, 182-189. (55) Sallen, L.; Oswald, P.; Sotta, P. J. Phys. II 1997, 7, 107-138.

Figure 1. Correspondence between the detailed molecular and the coarse-grained model used in this work. A head site represents an oxyethylene unit, whereas a tail site represents 3 CH2 units. The size parameters of the head and tail sites reproduce approximately the lengths and diameters of a certain “real” surfactant. On the right side an H4T4 surfactant is shown corresponding to the “real” C12EO4 molecule. The bottom row shows the representation of water by a single LJ site. respectively. Two different cutoff radii are used in the simulations. 1/6 The first is rcut ij ) 2 σij, which results in a purely repulsive interaction. The second is rcut ij ) 2σij, i.e., there is a short-range attractive part in addition to the repulsive part of the potential. Similar models using this type of potential in a related context can be found in refs 40, 56-59. The harmonic potential for bonded sites i and j is given by Uijb(rij) ) 1/2K(rij - r0)2, where K is a constant. The use of such a cut and shifted LJ potential with a small cut off radius represents a significant simplification and as a consequence will hardly allow exact quantitative studies, but rather will qualitatively mimic the true behavior of surfactant solutions. To enhance “reality” we chose the size parameter σ to reproduce approximately the lengths and diameters of a certain “real” surfactant, i.e., poly(oxyethylene), CmEOn. Here the oxyethylene units make up the hydrophilic part and the alkyl chain forms the hydrophobic part. In our simulations a simple head site (H) is chosen to represent one oxyethylene unit whereas a tail site (T) represents 3 CH2 groups (see Figure 1). The corresponding LJ parameters are compiled in Table 1. The LJ values for the water-water interactions are chosen to reasonably reproduce diffusion and density. The σ values for the interaction between different types of sites are calculated using the Lorentz combining rule (σij ) (σi + σj)/2). The constant K for bound particles is set to 837 kJ/mol Å2, and r0 is 3.4 Å in all simulations. We employ conventional MD simulations, numerically solving bi ) -∇bıU(r b1, ..., b rN), for all Newton’s equations of motion, (d2/dt2)r particles i ) 1, ..., N. The atomic equations of motion are integrated via a half-step Verlet algorithm60 with a time step of 2.5 fs. The most time-consuming part in the force evaluation at each MD step is the determination whether two atoms are within the cutoff distance. To reduce the resulting computational effort, which scales such as N2, where N is the number of atoms, to O(N), we employ the cell method.61,62 The cell method divides the simulation box into a lattice of M3 cells with M ) L/lcell, where L is the box length (for the case of a cubic simulation box) and lcell is the cell length with lcell > rcut. Each particle is assigned to its cell, an (56) Gottberg, F. K.; T. A. Hatton, K. A. S. J. Chem. Phys. 1997, 106, 9850-9857. (57) Siepmann, J. I.; Karaborni, S.; Smit, B. Nature 1993, 365, 330. (58) Smit, B.; Esselink, K.; Hilbers, P. A. J.; van Os, N. M.; Rupert, L. A. M.; Szleifer, I. Langmuir 1993, 9, 9-11. (59) Palmer, B. J.; Liu, J. Langmuir 1996, 12, 746-753. (60) Allen, M. P.; Tildesley, D. J. Computer Simulation of Liquids; Clarendon Press: Oxford, 1990. (61) Quentrec, B.; Brot, C. J. Comput. Phys. 1975, 13, 430-432. (62) Hockney, R. W.; Eastwood, J. W. Computer simulation using particles; McGraw-Hill: New York, 1981.

1628

Langmuir, Vol. 16, No. 4, 2000

Fodi and Hentschke

Table 1. LJ Parameter Sets Used in this Work P1

water-water water-head water-tail head-head head-tail tail-tail

(W-W) (W-H) (W-T) (H-H) (H-T) (T-T)

σ [Å]

rcut/σ

 [kJ/mol]

2.77 3.885 3.885 5.0 5.0 5.0

2 2 21/6 2 21/6 2

3.35 3.35 3.35 3.35 3.35 3.35

∆σ [Å] P2

P3

P4

P5

P6

(W-T) (H-H) (H-T) (T-T) (W-T) (H-H) (H-T) (T-T) (W-H) (W-T) (H-H) (H-T) (T-T) (W-T) (H-H) (H-T) (T-T) (W-H) (W-T) (H-H) (H-T) (T-T)

rcut/σ

21/6

∆ [kJ/mol] -1.67 -1.67 -1.67 -1.67 -1.67 -1.67 -1.67 -1.67

-0.75 -1.67 -1.67 -1.67 -1.67 -1.67 -1.67 -1.67 -1.67

-1.5 -0.75 -0.75 -0.75 -1.5 -0.75 -1.5 -0.75

21/6

-1.67 -1.67 -1.67 -1.67

a The σ values for the interactions between different particles (for example W-H) are evaluated using the Lorentz combining rule. Parameter sets P2-P6 are variations of P1 showing the nonzero differences relative to the P1 values, with exception of rCut/σ, for which the absolute values are shown.

operation of O(N), and the cutoff criterium is checked for particles in same and in neighboring cells only. The cell list is updated at fixed intervals (here between 20 and 40 time steps). The simulations are performed under NPT conditions (the number of particles, the pressure, and the temperature are constant) using the weak coupling method according to Berendsen et al.63 (temperature relaxation time 0.3 ps; pressure relaxation time 0.4 ps). The number of molecules varies between 30 000 and 100 000 for the different systems. The reference pressure is 1 bar. If not indicated otherwise, the simulations are carried out at T ) 300 K.

3. Results Influence of Concentration. Temperature-solute mole fraction (T-x) phase diagrams and the effect of chain length on the phase behavior have been intensely studied experimentally for a series of poly(oxyethylene) surfactants by Sjo¨blom and Mitchell and co-workers.64,65 Our first aim is to determine, whether we observe different aggregation behavior with increasing surfactant concentration. In these simulations we use parameter set P1 (cf. Table 1). Water-tail and head-tail interactions are purely repulsive, i.e., the corresponding cutoff radius is rcutij ) 21/6σij. All other interactions (water-water, water-head, head-head, and tail-tail) include an attractive part, i.e., rcutij ) 2 σij. First the influence of the surfactant mole fraction on the phase behavior is studied. The model surfactant is composed of four head and four tail sites, H4T4, which corresponds to the “real” C12EO4 system (cf. (63) Berendsen, H. J. C.; Postma, J. P. M.; van Gunsteren, W. F.; DiNola, A.; Haak, J. R. J. Chem. Phys. 1984, 81, 3684-3690. (64) Sjo¨blom, J.; Steinius, P.; Danielsson, I. In Nonionic surfactants; Schick, M. J., Ed.; Marcel Dekker: New York, 1987. (65) Mitchell, D. J.; Tiddy, G. J.; Warning, L.; Bostock, T.; McDonald, M. P. J. Chem. Soc., Faraday Trans. 1 1983, 79, 975.

Figure 1). We consider four weight fractions. Table 2 gives an overview of the simulated systems. Representative snapshots of the four simulated mole fractions are shown in Figure 2. All systems are started from a “randomized” isotropic configuration (cf. the inset), corresponding to a state far away from equilibrium. The upper left snapshot (11.6 wt %) shows the formation of a spherical micelle next to a rodlike one (consisting of 55 and 144 monomers, respectively). Even for very long simulation times (1 500 000 MD steps), we did not observe fusion of the two aggregates. At 20 wt % we obtain a single rodlike aggregate. The lower left snapshot shows the 66 wt % system. Here we observe the formation of a bicontinuous structure, where the hydrophilic heads form a pipelike network shielding the hydrophobic tails from the water. At 95 wt % the system favors a slightly perturbed lamellar morphology. Increasing the temperature by about 40 K did not change the structure of either the 66 or the 95 wt % system significantly. In addition to identifying the resulting structures via direct visual examination, we also obtain structural information from different site-site correlation functions. Selected examples are shown in the Figures 3 and 4. Note that the first head atom in the surfactant is labeled H1, the second is H2 andsproceeding along the surfactants the last atom is labeled T4 (see also Figure 1). When comparing the different weight fractions one must keep in mind that they are normalized to the corresponding “ideal gas state”, i.e., the state at which the particles are distributed uniformly over the entire volume. Consequently, we calculate the coordination numbers (CN) by counting all neighboring sites up to the first minimum of g(r) (first neighbor shell). For the H1-H1 correlation, the first minimum is at about 7.0 Å (lower panel of Figure 3). As expected, there is a much closer packing of the heads in the lamellar phase (CNH1-H1 ) 6.3) in comparison to the other three (66 wt %: CNH1-H1 ) 2.3; 20 wt %: CNH1-H1 ) 0.4; 11.6 wt %: CNH1-H1 ) 0.4). This is reasonable, because in a spherical micelle the heads are surrounded by water in addition to other head atoms driving the surfactant sites apart. This observation is confirmed by the pronounced peak at about 8.6 Å in the H1-H1 correlation of the 11.6 wt % system. This peak correlates with the head-water-head distance. To monitor the strength of solvation of the heads as a function of the mole fraction, the water coordination numbers of H1-H4, T1, and T4 are calculated. Again the first solvation shell (up to the first minimum at about 5.6 Å in the X-W radial correlation) is considered. Figure 5 shows that the results are complementary to the above CNH1-H1, i.e., the water content between heads is large for 11.6 and 20.0 wt %. It is very low at 95 wt %. The figure also illustrates the strong decrease of the water coordination going from H1 to H4, with almost no water found in the neighborhood of the tail end (T4). Figure 6 monitors the length of the end-to-end vector b ree as function of the different mole fractions. The increase of b ree with increasing surfactant concentration is a packing effect and is consistent with the findings of Figure 5. The thickness of the double layer can also be deduced from the H1-H1 correlation, i.e., the broad peak at about 5.3 nm. The T1-T1 correlation (upper panel of Figure 4) shows a pronounced peak in all four systems, with coordination numbers CNT1-T1 ) 7.1 for the 95 wt % system, CNT1-T1 ) 5.9 for the 66 wt %, and CNT1-T1 ) 5.3 for the 20 wt % and also the 11.6 wt % system. This is because in all three structures the terminal tail sites are in close proximity to each other, with the closest packing in the lamellar-like phase. The inset in the upper panel

Phase Behavior of Model Surfactant Solutions

Langmuir, Vol. 16, No. 4, 2000 1629

Table 2. Listing of All Simulation Snapshots Shown in this Work system

weight fraction

mole fraction

H4T4 H4T4 H4T4 H4T4 H2T6 H6T4 H8T4

11.6 % 20.0 % 66.0 % 95.0 % 66.0 % 66.0 % 66.0 %

0.007 0.013 0.094 0.5 0.094 0.077 0.064

parameter set (fig.) P1 (2) P1 (2) P1 (2) P1 (2) P2 (9) P2 (9) P2 (9)

P2 (8) P6 (9)

P3 (8)

N

Nsurf

30000 16600 30000 30000 35000 40000 28538

200 200 1700 3500 1986 1813 1080

a The masses of the corresponding C EO surfactant and of water are used to calculate the weight fractions. N denotes the total number m n of sites (T, H, and W), Nsurf is the number of surfactant molecules. PX refers to the parameter sets in Table 1. The numbers in brackets are the corresponding figures.

Figure 3. Radial pair correlation functions g(r). Here rW-H1 is the distance between water and first head atoms (upper panel), and rH1-H1 is the distance between first head atoms (lower panel). Figure 2. Snapshots of four different H4T4/water mixtures. For 11.6 and 20 wt % the water is omitted. Otherwise the shading corresponds to the model depicted in Figure 1. The inset shows the initial random configuration for the 95 wt % case.

of Figure 4 shows the pronounced ordering in the layers in the lamellar-like system. The lower panel in Figure 4 shows the pair correlation function for H2-T3. For the two highly concentrated systems two pronounced steps can be distinguished: the first corresponds to the intramolecular H2-T3 distance, the second to the distance between H2 and T3 in neighboring molecules orientated tail-to-tail to each other. These first results (micelles at low mole fractions, lamellar morphology at high mole fraction) of our model surfactant system are in accord with the general features of the phase behavior of our target surfactant system CmEOn, m ≈ 12, n ≈ 4. Furthermore, the stability of the morphologies against moderate temperature changes (see above) is also observed experimentally. Figure 7 shows the experimental phase diagrams of C12EO4 and C12EO6, taken from the study of Mitchell and co-workers.65 One main difference between our (limited) results thus far and the experimental phase diagram of C12EO4 64,65 is that this apparently lacks the bicontinuous phase found in the simulation. In the experiment the lamellar phase, LR, extends almost over the entire concentration regime on the surfactant-rich side of the phase diagram. Contrary to our findings, recent MC simulations of Larson,45 show with increasing surfactant a hexagonal phase followed by

Figure 4. Radial pair correlation functions, g(r), for the terminal tail sites (upper panel) and between the second head atom and the third tail atom (lower panel).

the gyroid and subsequently the lamellar phase. Experimentally, with increasing head length (EO5, EO6), LR is suppressed in favor of a hexagonal phase, accompanied by a narrow bicontinuous regime. We try to capture this tendency with our model in the next section. Parameter Dependence of the Phase Behavior. In this section we examine the influence of several model parameters on the phase behavior. Larson, in a number

1630

Langmuir, Vol. 16, No. 4, 2000

Fodi and Hentschke

Figure 5. Coordination numbers CNXW of the water solvation for different sites in the H4T4 system vs inverse surfactant weight fraction. H1 (circles), H2 (squares), H3 (up triangles), H4 (down triangles), T1 (crosses), T4 (stars).

Figure 6. Length of the end-to-end vector b ree in the H4T4 system vs inverse weight fraction.

of studies36,42,44,45 examined the phase behavior of singletail surfactants and its dependence on the chain length and the head-to-tail ratio. In contrast to Larson, who used MC lattice simulations, we also exploit the additional freedom of “off lattice” simulations which allows us to reproduce the size proportions of surfactant and water by using different surfactant size parameters (see below). Here we focus on the concentration regime, where we observe the bicontinuous structure. From the experimental point of view this region is particularly interesting because the lamellar, hexagonal, cubic, and bicontinuous phases can appear in this concentration regime, depending on the molecular architecture of the monomer. Guided by the results obtained with the parameter set P1 (see Table 1), we change the well depth for selected site-site interactions. More specifically, the water-water and the headwater parameters remain unaltered, whereas the W-T, H-H, H-T, and the T-T parameters are reduced, i.e., all water interactions aresin relationsstrengthened, and so the solvent’s role on the phase behavior becomes

Figure 7. Experimental phase diagrams of the C12EO4 (top) and the C12EO6 system (bottom) taken from ref 65. L1 ) micellar solution, LR ) lamellar phase, H1 ) hexagonal phase, V1 ) normal bicontinuous cubic phase, W ) aqueous solutions, L2 ) liquid surfactant containing dissolved water, L3 ) sponge phase, S ) solid surfactant.

more pronounced. This choice of parameters will be referred to as parameter set P2. The upper snapshot in Figure 8 demonstrates the drastic change in the morphology (cf. Figure 2, lower left). Instead of a bicontinuous structure a highly ordered lamellar phase is found. As mentioned above, this is what is observed experimentally for the “real” system (cf. Figure 7). As an aside we examined the effect of the parameter set P2 on the 11.6 wt % system. Here, three aggregates are found of similar size, which are less compact (not spherical, more “spongy”). Of course, the small number of assembled aggregates does not allow any quantitative analysis. The lower snapshot in Figure 8 shows again, for the 66 wt % case, what happens if, in addition to the altered well depths, the head-head interaction is made purely repulsive (P3). As before, a lamellar structure with a slightly curved interface (H1-H1 coordination number is diminished) is obtained. As mentioned above, one of the advantages of using “off lattice” MD simulation is that one can change not only the number of head and tail sites but also their sizes. To check the influence of size, we studied the H4T4 system at 66 wt % reducing σHX and σTX

Phase Behavior of Model Surfactant Solutions

Langmuir, Vol. 16, No. 4, 2000 1631

Figure 9. Different HmTn surfactants at 66 wt % using parameter sets P2 and P6 (reduced head size and only the repulsive part of the head-head potential).

Figure 8. Snapshots of the H4T4 system at 66 wt %, obtained with the indicated parameter sets.

(X ) W, T, H) (cf. parameter sets P4 and P5 in Table 1). For the smaller head sites the lamellar structure is preserved (even if, in addition, the head-head interaction is made purely repulsive), whereas for the reduced tail sites, a bicontinuous structure replaces the lamellar structure. To ensure that the morphologies in our simulations are not solely a consequence of finite size frustrations, we repeated selected calculations for different box sizes (e.g., 1.5 times the original box length) with up to 100 000 particles. Our findings are that the overall structures are largely independent of the box size except for some additional distortion at the smaller box sizes. Note that the large systems require more than 6 weeks of computing time on a DEC alpha workstation with 500 MHz for one point in the phase diagram. To study the influence of the head-to-tail ratio on the assembly process and how this is related to the experimental phase behavior, we consider three additional model

Figure 10. Schematic drawing of a defect in a hexagonal phase with bridges connecting neighboring columns. The drawing is taken from ref 72.

systems, H8T4, H6T4, and H2T6. Again we restrict ourselves to the 66 wt % case using parameter set P2 (with one exception). The resulting morphologies of our calculations are displayed in Figure 9. All three P2 systems (H8T4, H6T4, and H2T6) form bicontinuous structures (in contrast to the H4T4 system in Figure 8 (top)). To find possible higher symmetries in these structures we examined slices and projections of each system. One could possibly expect a distorted hexagonal columnar phase with bridge-defects, as illustrated in Figure 10, or even ordered cubic bicontinuous structures such as the double diamond or the gyroid-type phase. These structures normally occur between the lamellar and the hexagonal phase and often

1632

Langmuir, Vol. 16, No. 4, 2000

Fodi and Hentschke

are regarded as the structures by which the lamellar phase changes into the hexagonal phase. In addition to these phases, noncubic so-called intermediate phases with a wide variety of possible structures such as disrupted lamellar or rhomboedrical networks have been reported experimentally.51,52,66 For example, C16EO6 and C22EO6 (roughly corresponding to H6T5 and H6T7) are known to form such intermediate phases. However, comparing our structure with the aid of projections and slices to the ordered cubic phases leads to no definite conclusion. In ref 45, Larson describes the particulars of the simulation of the intermediate phases, where he outlines that the appearance of the gyroid and the rhomboedricallike intermediate phases are sensitive to the box dimensions and are very slow to equilibrate. Note that in order to avoid frustration effects in the H8T4 case, instead of using cubic box dimensions, we choose fixed box-axes with a ) b‚tan30° ) c (a, b, and c denote the x, y, and z dimensions of the box) corresponding to the expected hexagonal symmetry. If we examine the H6T4 system, an interesting feature is shown in the lower right snapshot of Figure 9. This simulation, in contrast to the one above, is carried out with parameter set P6, i.e., the head-head repulsion is purely repulsive and the size of the head sites is smaller. This change results in a lamellar morphology. This is reasonable and consistent with the results discussed above. Reducing the head size leads to a monomer architecture comparable to the H4T4 system, which also self-assembles to yield a lamellar structure under otherwise identical conditions. To distinguish between the different bicontinuous structures, we take into account the local curvature as follows. For each H1 site, i, a list of neighboring H1 sites, j, is computed. These sites fulfill the condition rij < rcrit, i.e., all j sites within a distance rcrit (“first” neighbors) are included. Here, rcrit ) 2.5 × 21/6 × σHH. Furthermore, all H1 sites, k, which fulfill rjk < rcrit are also added (“second” neighbors). To avoid the inclusion of H1 sites j or k which do not belong to the same “pipe” or “layer” as H1 site i, an additional criterium must be fulfilled. We check that at least one of the tail sites of every j-surfactant neighbor lies within rcrit of one of the tail sites of i. Analogously, we check that at least one tail site of every k surfactant lies within rcrit of the tail sites of surfactant j. Having computed the neighbor list for all H1 sites, i, we calculate the quantity

S Bi )

(

1

Nineig

)

∑ br li

l)1 Nneig i

(2)

for all i which have more than two first neighbors. Nneigi is the number of first and second neighbor sites of i and the l -summation includes all the neighbors in the list. To characterize the local curvature, we compute the projection Bi onto the end-to-end vector b riee ) b riTX - b riH1, where Siz of S TX stands for the terminal tail site (cf. the sketch in Figure 11). For a locally convex water-surfactant interface, Siz > 0, whereas for a concave water-surfactant interface, Siz < 0. A locally flat interface (including a saddle point) results in Siz ≈ 0. Notice that the rather large number of neighbors smoothes out the local roughness, e.g., the staggering due to interdigitation of neighboring surfactant molecules. The result of this analysis is shown in Figure 11, where the relative frequency, f, of Siz is plotted for the H4T4, H2T6, H6T4, and H8T4 systems. Except for H4T4, where (66) Funari, S. S.; Holmes, M. C. J. Phys. Chem. 1992, 96, 110291038.

Figure 11. The relative frequency, f, of Siz (see eq 2) for four different surfactant systems at 66 wt %. Dotted line: H4T4 (P1); solid line: H2T6 (P2), dashed line: H6T4 (P2), crosses: H8T4 (P2). The left inset illustrates the calculation of Siz. The right inset shows the result for H4T4 (P2). Table 3. The mean projection, 〈Siz〉, and Its Root-Mean-Square Deviation ∆Sz (See Eq 2) for Selected Systems system

parameter

〈Siz〉 [Å ]

∆Sz [Å ]

H2T6 (66 wt %) H4T4 (66 wt %) H6T4 (66 wt %) H8T4 (66 wt %) H4T4 (66 wt %) H4T4 (11.6 wt %)

P2 P2 P2 P2 P1 P1

0.59 0.30 1.65 2.49 1.86 3.25

2.30 2.02 3.10 2.46 2.23 2.43

parameter set P1 is used, these results are obtained using P2. To clarify our analysis, an inset is added in Figure 11, where we plot Siz for the lamella morphology obtained for H4T4 using P2. In addition to Figure 11, we calculate the mean projection, 〈Siz〉i, and the root-mean-square deviation (〈S2iz〉 - 〈Siz〉2)1/2 for a number of our systems. The results are compiled in Table 3. For H4T4 (P2) the value of 〈Siz〉 in Table 3 indeed is close to zero as expected even though the distribution is rather broad. For the bicontinuous structures, one can see that all curves are shifted toward positive values. This means that in all four systems the surfactants prefer to reside on a curved interface with the headgroup on the (predominantly) convex side. Even in the H2T6 case there is no change toward an inverse structure as is often observable when the surfactant tail is long compared to the head.45 But in comparison to the other systems, Siz for H2T6 is significantly diminished. The effect of the lengthening of the surfactant’s head can be explained by curvature arguments, too.67-70 In the H6T4 and the H8T4 cases, the heads prefer curved structures, where there is more space for the heads. So, in contrast to the H4T4 system, which undergoes the ordering transition to the lamellar structure when the P2 parameters are used, the H6T4 and H8T4 systems do not. When the (67) Israelachvili, J.; Mitchell, D. J.; Ninham, B. W. J. Chem. Soc., Faraday Trans. 1 1976, 72, 1525-1568. (68) Helfrich, W. J. Phys. I 1987, 48, 291-295. (69) Sadoc, J. F.; Charvolin, J. J. Phys. I 1986, 47, 683-691. (70) Hyde, S. T. Colloq. Phys. 1990, 51, C7-209-C7-228.

Phase Behavior of Model Surfactant Solutions

headgroup is longer, higher surfactant concentrations are necessary to induce the shift to morphologies with less curvature. 4. Conclusion We performed MD simulation studies of the phase behavior of surfactant solutions using a coarse-grained model. Our aim was to gain insight in how far the MD technique in combination with such a model can be used to observe self-organization of amphiphilic systems at higher concentration. Even though we only employed a simple cut and shifted LJ potential with two different types of LJ sites, the model showed some of the characteristic features of real surfactants’ phase behavior. We found that within our simplified model the surfactants self-organize into very different micellar phases. This is consistent with the conclusions of other authors,71 who found that the influence of the water structure including hydrogen bonding is, at least for the organization of certain structures, not the driving force. (71) Karaborni, S.; van Os, N. M.; Esselink, K.; Hilbers, P. A. J. Langmuir 1993, 9, 1175-1178.

Langmuir, Vol. 16, No. 4, 2000 1633

We investigated the concentration dependence of our model system at four mole fractions. for a certain symmetric (H4T4) model surfactant. For one of the weight fractions (66 wt %) we explored the effect of varying the model parameters on the morphology. This was done bearing in mind that small differences in physicochemical details can lead to significant changes in the phase behavior. First we reduced the well depth parameter  for some of the site-site interactions. In addition we examined the effect of the head (tail) size and the head-to-tail ratio. At these highly concentrated conditions (66 wt %) we found lamellar and bicontinuous structures. Our algorithm for analyzing the local curvature appeared to be a sensitive tool for the distinction of these phases. The fact that we cannot find hexagonal ordering in the concentration regime where it is seen experimentally may be explained by the simplicity of the model or by the still relatively small system size causing distortion of the morphology. LA990862H (72) Sallen, L.; Sotta, P.; Oswald, P. J. Phys. Chem. B 1997, 101, 4875-4881.