Langmuir 2004, 20, 5123-5126
5123
Simulation of Phase Equilibria in Lamellar Surfactant Systems Martin Turesson,* Jan Forsman, Torbjo¨rn A° kesson, and Bo Jo¨nsson Theoretical Chemistry, University of Lund, P.O. Box 124, S-221 00 Lund, Sweden Received January 15, 2004. In Final Form: April 7, 2004 The coexistence of two lamellar liquid crystalline phases has been investigated by means of Monte Carlo simulations. The surfaces of the negatively charged bilayers formed by the surfactant molecules are modeled as planar infinite walls with a uniform surface charge density. Water is treated as a dielectric continuum, and only electrostatic interactions are considered. The counterions are mono- and divalent point ions, and their ratio is allowed to vary. Monovalent counterions lead to a repulsive osmotic pressure at all separations, while an attractive region exists when the counterions are divalent. In the latter case, one would expect a phase separation to take place, although it is not observed experimentally due to the limited stability of the lamellar phase at high water content. In a system with mixed counterions, however, the osmotic pressure exhibits a van der Waals loop under such conditions that two phases can coexist. A phase diagram is constructed, and the agreement with experimental data is excellent.
Introduction When amphiphilic molecules are mixed with a solvent (e.g., water), they spontaneously self-organize, at some critical concentrations, to form spherical micelles, cylindrical micelles, planar bilayers, bicontinuous structures, etc.1 A system consisting of lamellar liquid crystalline phases formed by the anionic amphiphile bis(2-ethylhexyl)sulfosuccinate (AOT) has been investigated by Khan et al.2,3 They found that NaAOT forms a lamellar phase if the water content is between 25 and 80 wt %. Substituting sodium with calcium or magnesium does not effect the lower limit, but the ability to incorporate water diminishes from 80 to 40 wt %. The ternary system NaAOT-water-Ca(AOT)2 was also investigated,4 and coexistence between two lamellar phases with different water contents was discovered. At the same time, Guldbrand et al.5 performed Monte Carlo simulations of systems with parameters chosen to mimic the AOT system. A nonmonotonic behavior of the pressure curve was obtained for systems with an intermediate electrostatic coupling, indicating the possibility of a phase separation into two lamellar phases. The same problem has also been studied by Moreira and Netz using both Monte Carlo simulations and a field theoretic approach.6,7 The origin of water depletion in the Ca(AOT)2-water system and the nonmonotonic behavior of the simulated force curve was attributed to electrostatic correlations between the ions. For weakly coupled systems, reasonably correct osmotic pressures can be predicted from the PoissonBoltzmann (PB) equation, although the mean field treatment overestimates the repulsion in the system.8 For more highly coupled systems with divalent counterions or higher surface charge densities, the PB equation is qualitatively (1) Evans, D. F.; Wennerstro¨m, H. The Colloidal DomainsWhere Physics, Chemistry, Biology and Technology Meet; VCH Publishers: New York, 1994. (2) Khan, A.; Fontell, K.; Lindman, B. J. Colloid Interface Sci. 1984, 101, 193. (3) Khan, A.; Fontell, K.; Lindman, B. Colloids Surf. 1984, 11, 401. (4) Khan, A.; Jo¨nsson, B.; Wennerstro¨m, H. J. Phys. Chem. 1985, 89, 5180. (5) Guldbrand, L.; Jo¨nsson, B.; Wennerstro¨m, H.; Linse, P. J. Chem. Phys. 1984, 80, 2221. (6) Moreira, A. G.; Netz, R. R. Phys. Rev. Lett. 2001, 87, 8301. (7) Moreira, A. G.; Netz, R. R. Eur. Phys. J. 2002, E8, 33.
wrong, due to the neglect of ion-ion correlations. Increasing the counterion valency or the surface charge, as well as decreasing the dielectric constant of the solvent, leads to stronger correlation effects. In a surfactant system, it is usually difficult to tune the electrostatic coupling by varying the surface charge density or the dielectric permittivity of the solvent and the counterion valency can only be changed in discrete steps. A mixture of monoand divalent counterions, on the other hand, is an alternative where the “average electrostatic coupling” can be continuously varied. In this report, we use Monte Carlo (MC) simulations in a variant of the grand canonical ensemble to calculate the phase diagram of a mixture of two surfactants, one with monovalent counterions and the other with divalent counterions. The focus is on the lamellar liquid crystalline part of the phase diagram, where two phases of different water contents can exist in equilibrium. A similar type of phase coexistence has also been observed in the micellar phase region.9 Model and Simulations The ternary NaAOT-water-Ca(AOT)2 system is modeled by two infinite parallel surfaces with a uniform surface charge density, σ; see Figure 1. The neutralizing counterions are allowed to move freely in the region between the surfaces. The simulations were performed at a temperature of 298 K, with a surface charge density of one unit charge per 65 Å2. In total, there were 300 charges distributed between mono- and divalent ions. By doubling the size of the system and the simulation length, the convergence with respect to the ensemble averages was checked, and no noticeable change was observed. A primitive model is adopted, in which the solvent is treated as a continuum, with a dielectric permittivity equal to that of bulk water, r ) 78. The counterions, with valencies of zi ) 1 or 2, are treated as point charges. The energy of interaction, u(rij), between two particles, i and j, separated (8) Jo¨nsson, B.; Linse, P.; Åkesson, T.; Wennerstro¨m, H. In Surfactant in Solution; Mittal, K. L., Lindman, B., Eds.; Plenum Publishing Corporation: New York, 1984. (9) Svensson, A.; Piculell, L.; Karlsson, L.; Cabane, B.; Jo¨nsson, B. J. Phys. Chem. B 2003, 107, 8119.
10.1021/la049858i CCC: $27.50 © 2004 American Chemical Society Published on Web 05/11/2004
5124
Langmuir, Vol. 20, No. 12, 2004
Turesson et al.
Figure 2. Counterion distribution profiles, as a function of distance, z, from one wall, for monovalent ions (long-dashed curve) and divalent ions (solid curve). β∆µ ) 4. Figure 1. Schematic picture of a lamellar phase with AOT as the building block of the bilayers. The parallel lines, separated by a distance of h, define the charged walls in the simulation with a surface charge density, σ.
by a distance of rij is calculated according to Coulomb’s law.
u(rij) )
zizje2 4π0rrij
(1)
0 and e are the permittivity of vacuum and the elementary charge, respectively. The counterions also interact with the charged surfaces within the MC box, and these interactions are treated as described elsewhere.10 Electrostatic interactions are long ranged, and a tail correction due to Coulombic interactions ranging outside the MC box is included, using a technique described by Valleau et al.11 Repulsive undulation forces from the nonrigid bilayers12 and steric repulsions from protruding headgroups13 are not taken into account. At separations shorter than 2030 Å, oscillatory forces due to solvent packing can appear,14 but they are not included in this study. To observe a phase separation, a mixture of mono- and divalent counterions is necessary according to experiments.4 This means that the chemical potential of the two types of counterions will be the same in the two phases, although their relative compositions will not. In the simulations, we maintain chemical equilibrium by a modified grand canonical simulation technique, where two monovalent particles are substituted with one divalent particle and vice versa.15 This technique implies a varying composition of counterions when the separation between the surfaces is changed. The parameters are chosen to represent the AOT system. The difference in chemical potential, ∆µ ≡ µ2+ - 2µ+, between mono- and divalent ions was fixed at all separations, which together with electroneutrality ensures that the chemical potentials for the two salt species remain constant. One-third of the trial configurations were generated by ordinary Metropolis displacements,16 and two-thirds, by equal parts of trial (10) Jo¨nsson, B.; Wennerstro¨m, H.; Halle, B. J. Phys. Chem. 1980, 84, 2179. (11) Valleau, J. P.; Ivkov, R.; Torrie, G. M. J. Chem. Phys. 1991, 95, 520. (12) Helfrich, W. Z. Z. Naturforsch., A: Phys. Sci. 1978, 33, 305. (13) Israelachvili, J. N.; Wennerstro¨m, H. Langmuir 1990, 6, 873. (14) Tang, Z.; Scriven, L. E.; Davis, H. T. J. Chem. Phys. 1994, 100, 4527. (15) Allen, M. P.; Tildesley, D. J. Computer Simulation of Liquids; Oxford University Press: Oxford, U.K., 1989. (16) Metropolis, N. A.; Rosenbluth, A. W.; Rosenbluth, M. N.; Teller, A.; Teller, E. J. Chem. Phys. 1953, 21, 1087.
Figure 3. Curves that show how the fraction of monovalent ions, fmon ) Cmon/(Cmon + Cdi), in the slit increases with separation at constant chemical potentials of the counterions. β∆µ ) 5. The area/headgroup for the different curves is indicated in the figure.
insertions of two monovalent ions and deletion of one divalent ion or vice versa. The osmotic pressure was calculated across the midplane of the slit:17
∑i Ci(mp) + pes
Posm ) kBT
(2)
where Ci(mp) is the counterion concentration of species i at the mid-plane, pes is the electrostatic force per unit area acting across the mid-plane,11 and kBT is the thermal energy. Results From the simulations, the distribution of ions in the slit is obtained. As can be seen in Figure 2, the divalent ions are on average more attracted to the walls compared to the case of the monovalent ions, leading to a lower mid-plane concentration of divalent ions. The overall pressure can therefore be negative due to the ion-ion correlation attraction; see eq 2. The input parameter ∆µ can be used to regulate the ratio of mono- and divalent counterions at some fixed distance. The system is more strongly coupled at short separations, and as a consequence, the fraction of divalent ions, fdi ) Cdi/(Cmon + Cdi), increases as the separation is reduced, as demonstrated in Figure 3. Cmon and Cdi are the average concentrations of monovalent and divalent counterions in the slit, respectively. Similarly, fdi increases with increasing surface charge density at a fixed separation and ∆µ. The reduction of the number of monovalent ions at short separations is partially a trivial effect; that is, a fluid that can form dimers at increased concentration will do so. In this case, the effect is more pronounced (17) Wennerstro¨m, H.; Jo¨nsson, B.; Linse, P. J. Chem. Phys. 1982, 76, 4665.
Phase Equilibria in Lamellar Surfactant Systems
Figure 4. Pressures for the limiting cases when fdi ) 1 (longdashed curve) and fdi ) 0 (solid curve). The short-dashed curve shows a system with β∆µ ) 4.
Figure 5. Osmotic pressure as a function of separation for varying surface areas/headgroups. The chemical potential difference between mono- and divalent ions is constant, β∆µ ) 5, for each curve; see Figure 3.
Figure 6. Enlarged curve from Figure 5; the area/headgroup is 81 Å2 and β∆µ ) 5. The horizontal line is the result from a Maxwell construction and is chosen to make areas 1 and 2 equal. It intersects the curve at two separations, h1 and h2, marked with circles, which determines the water content of the two coexisting lamellar phases.
because divalent ions are energetically favored in narrow slits. This is due to the increased correlation of the divalent ions at short separations, where the density of counterions is high. Simulated pressure curves for systems with fdi ) 1 and fdi ) 0 are seen in Figure 4 together with an intermediate pressure curve corresponding to β∆µ ) 4. With a mixture of mono- and divalent ions, one obtains pressure curves which sometimes display a nonmonotonic behavior with a van der Waals loop. Another set of simulated pressure curves are shown in Figure 5. The attraction increases with coupling strength, which in this case is modulated by the surface charge. Figure 5 shows that, for sufficiently coupled systems, the force between the surfaces can become attractive. An enlargement of one of the curves, Figure 6, demonstrates the van der Waals loop from which equilibrium distances can be calculated. These can be found via a Maxwell equalarea construction. Once the equilibrium distances, h1 and h2, and the composition of counterions at these separations
Langmuir, Vol. 20, No. 12, 2004 5125
Figure 7. Two-phase regions for three values of β∆µ; from the bottom to the top, β∆µ ) 5, 6, and 11, respectively. In the same direction, the number of divalent ions increases. The shortdashed horizontal line indicates equilibrium distances calculated by a Maxwell construction.
Figure 8. Comparison between experimental (b) and simulated (+) data with an area/headgroup of 65 Å2. The asterisks (/) show data where 10 Å has been added to the simulated equilibrium distances.
are determined, a phase diagram can be constructed and compared with experimental data. In Figure 7, the two-phase regions calculated from the Maxwell constructions are shown with the area/headgroup and the chemical potential difference, β∆µ, as variables. Thus, a low area/headgroup, that is, a high surface charge density, and a large β∆µ corresponds to a strongly coupled system. From the bottom to the top in Figure 7, the chemical potential and consequently the number of divalent ions in the system increases. If fdi is increased and we still want to maintain approximately the same water content in the two phases, we have to reduce the surface charge density. On the other hand, if we reduce the surface charge density at a constant β∆µ, we will eventually reach a critical point. It is difficult to accurately establish the critical point, which therefore is only schematically drawn in the figure (long-dashed lines in Figure 7). To construct a phase diagram from the simulated data, necessary experimental parameters, determined by Fontell,18 were used. The thickness of the bilayer was set to 20 Å. The densities of the water and bilayer regions were estimated as 1.0 and 1.2 g/cm3, respectively. The headgroup area of the fully dissociated AOT molecule was set to 65 Å2. One problem associated with such quantitative comparisons is that experimentally measured bilayer separations are rather poorly defined and cannot be unambiguously mapped onto our model system. In Figure 8, we see that, if 10 Å is added to each equilibrium separation, a much better agreement with the experimental data is obtained. Given the simplicity of our “hard wall” model system and the ambiguous way in which the separation (18) Fontell, K. J. Colloid Interface Sci. 1972, 44, 318.
5126
Langmuir, Vol. 20, No. 12, 2004
is defined experimentally, the agreement between theory and experiment is excellent. In ref 4, a theoretical model is described where a repulsive force from the Poisson-Boltzmann equation is combined with an attractive van der Waals force. A phase diagram for a mixed surfactant system is calculated with this model. Their calculated phase diagram looks very similar to our simulated one, but the Hamaker constant needed in order to obtain agreement with experiment is ∼20 times larger than those for typical hydrocarbonwater systems. This indicates that the repulsive pressure predicted from the PB approximation is much too high and that the mean field theory does not provide a plausible mechanism for the phase separations. Note also that the deficiency of the PB equation not only comes from the complete neglect of the correlation term, pes, in eq 2 but also from a too high mid-plane concentration, kBTCi(mp). It is actually the latter that is numerically most significant.21 Conclusions
Turesson et al.
primitive model system, show the existence of a two-phase region; that is, two lamellar phases with different water contents are found to be in equilibrium. Thus, the additional attractive force not accounted for in the standard Derjaguin-Landau-Verwey-Overbeek (DLVO) theory19,20 can play an important role. Despite the simplicity of the model, with only Coulombic interactions between the ions, the agreement between experiment and theory is surprisingly good. The mechanism for the additional attraction can be understood in terms of ion-ion correlations. These act in two ways: (i) they tend to accumulate the counterions close to the charged bilayer, thereby reducing the entropic repulsion, and (ii) the correlations of ions across the midplane of the aqueous layer give rise to a direct attractive interaction. In numerical terms, the former effect dominates and for really strongly coupled systems the entropic term can be several orders smaller than that predicted by the PB equation.
Phase diagrams for lamellar liquid crystalline systems, calculated from grand canonical simulations of a restricted
LA049858I
(19) Verwey, E. J. W.; Overbeek, J. Th. G. Theory of the Stability of Lyophobic Colloids; Elsevier Publishing Company, Inc.: Amsterdam, The Netherlands, 1948. (20) Derjaguin, B. V.; Landau, L. Acta Phys. Chem. 1941, 14, 633.
(21) Jo¨nsson, B.; Wennerstro¨m, H. When ion-ion corelations are important in charged colloidal systems. Electrostatic effects in soft matter and biophysics; Kluwer Academic Publishers: Dordrecht, The Netherlands, 2001.
