Langmuir 2000, 16, 7493-7502
7493
Mesoscopic Simulation of Polymer-Surfactant Aggregation Robert D. Groot Unilever Research Vlaardingen, P.O. Box 114, 3130 AC Vlaardingen, The Netherlands Received January 4, 2000. In Final Form: June 8, 2000 The formation of a polymer-surfactant complex has been studied in bulk solution using dissipative particle dynamics simulation of a system containing polymer, surfactant, and water. Widom insertion is used to obtain the surfactant chemical potential. Two distinct modes of surfactant adsorption on the polymer have been observed, depending on the interaction mechanism. These are (a) continuous adsorption and (b) polymer-surfactant aggregation via discrete micelles. If the binding between polymer and surfactant is dominated by hydrophobic interaction (i.e., when the surfactant tails tend to adsorb on the polymer), continuous adsorption is found, leading to molecular bottlebrush or swollen cage conformations. If the binding between polymer and surfactant is dominated by the surfactant headgroups, micellar adsorption is found, leading to a necklace of micelles on a polymer backbone. For each interaction strength between polymer and surfactant headgroup, there is an interaction between polymer and surfactant tail where the mode of adsorption crosses over from micelle adsorption to molecular bottlebrush conformations. This crossover is comparable to a critical point in the condensation of gases to liquids. A schematic phase diagram indicating where the various adsorption modes occur is given. The polymer endpoint separation as well as its swelling exponent are found to pass through a minimum as the surfactant concentration is increased. When the bulk critical micelle concentration is reached, the polymer swells beyond its original size.
1. Introduction When surfactant is dissolved in a polymer solution, it often interacts with the polymers. The generally accepted picture is that complete micelles adsorb on the polymer,1-3 leading to a necklace of micelle pearls on a polymer backbone.4 The surfactant concentration where surfactant starts to adsorb on the polymer is the critical aggregation concentration (CAC), which is generally lower than the surfactant concentration, or critical micelle concentration (CMC), at which micelles are formed in absence of polymers. For a review of the recent literature on polyelectrolyte-surfactant interaction the reader is referred to Hansson and Lindman.5 Kwak edited another survey of the field, with overviews of the nature of the complex, phase behavior and applications.6 Apart from thermodynamic quantities such as CMC and CAC, another experimentally accessible quantity that characterizes the complex is the radius of gyration of a polymer in the presence of surfactant. For the poly(ethylene oxide) (PEO) and sodium dodecyl sulfate (SDS) system, it is found that on increasing SDS concentration the polymer initially reduces in size, but when the surfactant concentration is increased beyond a certain point the polymer swells again.7,8 A similar behavior is found for hydroxipropylcellulose with SDS and with cetylpyridinium chloride.9 Other direct experimental data on the structure of this complex has come forward from (1) Nagarajan, R. J. Chem. Phys. 1989, 90, 1980. (2) Ruckenstein, E.; Huber G.; Hoffmann, H. Langmuir 1987, 3, 382. (3) Goddard, E. D.; Ananthapadmanabhan, K. P. Interactions of Surfactants with Polymers and Proteins; CRC Press: London, 1993. (4) Cabane, B. J. Phys. Chem. 1977, 81, 1639. (5) Hansson, P.; Lindman, B. Curr. Opin. Colloid Interface Sci. 1996, 1, 604. (6) Kwak, J. C. T., Ed.; Polymer-Surfactant Systems; Marcel Dekker: New York, 1998. (7) Claesson, P. M.; Fielden, M. L.; Dedinaite, A.; Brown, W.; Fundin, J. J. Phys. Chem. B 1998, 102, 1270. (8) Mears, S. J.; Cosgrove, T.; Obey, T.; Thompson, L.; Howell, I. Langmuir 1998, 14, 4997.
SANS data on PEO/SDS mixtures by Chari et al.10 They suggest that the polymer resembles a swollen cage, rather than a necklace around SDS micelles. Fluorescence measurements on the same system indicate that the aggregation number of SDS is low at the onset of binding but increases with surfactant concentration where the aggregate forms an elongated rod.11 From rheology measurements the PEO coil is found to expand due to SDS adsorption.12 Direct microscopic observations of DNA in C16TAB solutions were made by Mel’nikov et al.13 They found a coil-globule transition taking place at the CAC and a small range of surfactant concentrations where coils and globules coexist. The polymers were found to be either completely in the coil state or completely in the globule state, but none were seen in an intermediate state. Adsorption of the cationic surfactants n-dodecyldimethylamine14 and alkyltrimethylammonium bromide15 on various nonionic polymers showed the aggregation number in the polymer-bound complex to reduce, and the binding cooperativity to increase when the electrostatic interaction is increased. To understand these effects the experimental trends have to be compared with the predictions of either a theory or simulation. One important contribution to the theory of micellization started with the work of Scheutjens et al. who calculated the CMC, as well as the micellar size and shape from self-consistent field theory.16 However, the generalization to the aggregation of surfactant in the (9) Shimabayashi, S.; Uno, T.; Oouchi Y.; Komatsu, E. Prog. Colloid Polym. Sci. 1997, 106, 136. (10) Chari, K.; Antalek, B.; Lin, M. Y.; Sintra, S. K. J. Chem. Phys. 1994, 100, 5294. (11) Vanstam, J.; Brown, W.; Fundin, J.; Almgren, M.; Lindblad, C. ACS Symp. Ser. 1993, No. 532, 194. (12) Brackman, J. C. Langmuir 1991, 7, 469. (13) Mel’nikov, S. M.; Sergeyev, V. G.; Yoshikawa, K. J. Am. Chem. Soc. 1995, 117, 2401, 9951. (14) Brackman, J. C.; Engberts, J. B. F. N. Langmuir 1992, 8, 424. (15) Hansson, P.; Almgren, M. J. Phys. Chem. 1996, 100, 9038. (16) Leermakers, F. A. M.; Scheutjens, J. M. H. M. J. Colloid Interface Sci. 1989, 136, 231.
10.1021/la000010d CCC: $19.00 © 2000 American Chemical Society Published on Web 08/25/2000
7494
Langmuir, Vol. 16, No. 19, 2000
presence of polymer is far from trivial; one contribution at this point has been made recently by Wallin and Linse.17 Another relatively recent theory of the interaction between polymers and micelles in dilute solutions was put forward by Nikas and Blankschtein,18 who use a thermodynamic approach similar to earlier work of Nagarajan1 and Ruckenstein et al.2 Their approach relies on the assumed shape of micelles and is based on previously derived semiempirical and theoretical relations to describe various contributions to the free energy. The approach taken by Sear19 is to calculate the partition function of a polymer between two neighboring adsorbed micelles, and from this to calculate the partition function of a polymer coil with a necklace of micelles. These approaches contain valuable elements for a theory of polymer-surfactant interaction, but they are not free from assumptions regarding the structure of the complex. The theory by Wallin and Linse17 and the thermodynamic approach18 allow for micelles of variable size, but they are assumed to be spherical. An alternative model has also been proposed, where surfactant molecules are assumed to adsorb on specific binding sites on the polymer. This is the one-dimensional Ising model, in this context also known as the Satake-Yang model, which can be solved exactly.20 The CAC in this model comes forward as a cooperative transition arising from nearest neighbor interaction. This physically models a bottlebrush conformation. For most experimentally studied polymer-surfactant systems it is accepted that the former picture applies, i.e., a polymer decorated by micelles, rather than a molecular bottlebrush or swollen cage. However, most experimental attention has been given to oppositely charged systems. As the nature of polymer-surfactant interaction in these systems is very similar, it is not surprising that the same model applies to all of these systems. Indeed, not all experimental systems confirm to this rule. For instance Zana mentions polysoaps and oppositely charged surfactants for which the surfactant aggregation number increases linearly with surfactant concentration.21 Another example was mentioned by Hines and Thomas, who studied the binding of various surfactants to sodium carboxymethyl cellulose, and found aggregation numbers varying between 1 and 16, depending on the surfactant type.22 The question thus arises when surfactant can be expected to adsorb on polymers in the form of discrete micelles, and when it can be expected to follow a continuous adsorption mode. Such a question can be answered by a simulation model in which no assumptions are made regarding the structure or size of the aggregate. A correct statistical mechanical approach is to calculate the chemical potential of surfactant from a simulation where polymer, surfactant and water are represented by freely moving molecules. This can be done with standard methods such as Widom insertion and Biased Sampling Monte Carlo techniques.23 From such simulations one finds the Gibbs free energy change of the system on adding one surfactant molecule. Using this method the critical micelle concentration (CMC), or in the presence of polymer, the CAC can be calculated straightforwardly. The purpose of this paper (17) Wallin, T.; Linse, P. Langmuir 1998, 14, 2940. (18) Nikas, Y. J.; Blankschtein, D. Langmuir 1994, 10, 3512. (19) Sear, R. P. J. Phys. Condens. Matter 1998, 10, 1677. (20) Shirahama, K. In Polymer-Surfactant Systems; Kwak, J. C. T., Ed.; Marcel Dekker: New York, 1998; Chapter 4. (21) Zana, R. In Polymer-Surfactant Systems; Kwak, J. C. T., Ed.; Marcel Dekker: New York, 1998; Chapter 10. (22) Hines, J. D.; Thomas, R. K., 1998, unpublished results. (23) Frenkel, D.; Smit, B. Understanding Molecular Simulation; Academic Press: London, 1996.
Groot
is to describe the application of this method to a polymersurfactant system. In the literature no report is made of simulations of this kind. When polymers and surfactants are simulated, usually one of the two following assumptions is made: the surfactants are either assumed to form micelles which are subsequently modeled by hard spheres, or the surfactant is assumed to bind to the polymer irreversibly. For instance, the interaction between polyelectrolytes and oppositely charged surfactant micelles has been simulated by Wallin and Linse.24 The modeling of explicit surfactant is not included in this study, and the spherical nature and size of the micelle is fixed. The other assumption has been adopted by, e.g., Balazs and Hu,25 who report a simulation study of the aggregation of polymers and surfactants on a 3D lattice. They are mainly concerned with the formation of loops and networks for polymers with discrete localized stickers. A similar model has been used by Chakrabarti and Toral,26 who studied cluster forming of polymersurfactant mixtures on a 2D lattice. A more recent study based on a continuum model was published by Saariaho et al.27 who model a polymer with many sidearms: molecular bottlebrushes. The endpoint distribution of bottlebrushes is a point of debate. Various theories exist, predicting a wide variety of scaling exponents, but what is clear is that the persistence length of the polymer increases with the length of the side chains and with their concentration, i.e., with the degree of substitution of the backbone.27 To predict when in a polymer-surfactant system such molecular bottlebrushes are formed, and when the surfactant adsorbs as micelles, a new simulation technique is employed here, where both the polymer and the surfactant molecules are represented by strings of soft spheres. For this model the chemical potential of surfactant in the presence of polymer can be obtained relatively easy using the Widom insertion method. In the present work a number of examples of polymer-surfactant interactions are described, from which we can deduce when we have micelle binding and when a continuous binding process is pertinent. 2. Simulation Model and Methods 2.1 Simulation Model. The length and time scales pertinent to polymer-surfactant aggregation are beyond the reach of a full atomistic simulation. However, in this study we are not concerned with atomistic details but rather with the structure of the whole complex. Therefore it is unnecessary to have a resolution to the angstrom scale: a resolution to the nanometer scale suffices. The simulation strategy is therefore to group atoms together into single “beads” and use these centers of mass as new simulation entities. This way we apply a real-space renormalization principle to simulate directly on the longer (mesoscopic) length scale that is relevant for the problem at hand. When atoms are grouped together, the effective interaction between the clusters becomes much softer than the original (Lennard-Jones) potential,28 because no atoms need to overlap if two centers of mass coincide. Consequently the equation of state does not depend too much on the precise shape of the interaction potential. Therefore we may choose a convenient functional form, (24) Wallin, T.; Linse, P. J. Chem. Phys. 1996, 100, 17873; 1997, 101, 5506. (25) Balazs, A. C.; Hu, J. Y. Langmuir 1989, 5, 1230. (26) Chakrabarti, A.; Toral, R. J. Chem. Phys. 1989, 91, 5687. (27) Saariaho, M.; Ikkala, O.; Szleifer, I.; Erukhimovich, I.; ten Brinke, G. J. Chem. Phys. 1997, 107, 3267. (28) Forrest, B. M.; Suter, U. W. J. Chem. Phys. 1995, 102, 7256.
Polymer-Surfactant Aggregation
Langmuir, Vol. 16, No. 19, 2000 7495
provided that we take the parameters such that the thermodynamics of the underlying system is reproduced correctly.29 For the present application the size difference between a polymer, surfactant molecules, and surfactant micelles has to be resolved, which poses an extra problem: a clash of length scales. To solve this problem we have to restrict ourselves to the simplest possible model, hence the surfactant molecules are represented by two beads tied together by a harmonic spring. The modeling strategy now has two aspects. First, surfactant molecules must be mapped onto this dumbbell model, and second the trends predicted by the model have to be checked to experimental data. If we would find that for whichever parameter set the dumbbell model cannot reproduce the experimental data, then we would have ruled out this simplest (dumbbell) model. In that case there would be no point in pursuing the first step. For this reason we concentrate on the second step here. The polymer is modeled by a string of beads; polymer and surfactant beads are connected by harmonic springs of spring constant C ) 4. The repulsive interaction potential between the beads is given by
uij(r) ) 1/2 aij(1 - r/Rc)2θ(1 - r/Rc)
(1)
where θ is the step function, Rc is the cutoff range for the potential, and i and j represent particle types. There is a close relationship between this soft repulsive sphere model and the well-known Flory-Huggins model for polymer interactions.29 This relationship allows us to translate the repulsion parameters aij into the more familiar FloryHuggins interaction parameters χij. If we choose the bead density equal to FRc3 ) 3, then a liquid with the compressibility of water is simulated when we take aww ) 25 kT (between water beads), and the solubilities of the other phases are matched by29 χij ) 0.306 (aij - 25). This model has previously been applied to polymers in the melt30,31 and in solution.32,33 These results show that even polymer chains as short as L ) 10 beads follow the correct endpoint distribution and are characterized by the correct scaling exponents. Some minor inaccuracies in the results of Kong et al. can be attributed to the fact that their simulations did not satisfy the fluctuation-dissipation theorem. The molecular weight dependence of the interfacial tension between polymer melts that comes forward from this simulation model compares better with experiments than mean-field lattice theory does.30 As mentioned, for this initial study we restrict ourselves to a minimal model, and use a polymer length of L ) 50 beads, following Wallin and Linse.24 The interaction with the solvent is chosen such that the polymer is at its θ point, i.e., it has a Gaussian endpoint distribution. A polymer of 50 segments is at its θ point at χc ) 0.65, hence according to the above translation rule the polymer-water repulsion parameter has to be chosen as29 apw ) 27.12. To model a two-bead surfactant that readily agglomerates in spherical micelles, the head-head repulsion was increased and the tail-tail repulsion was decreased relative to the water-water repulsion. Obviously, this is only a qualitative model of surfactant, but it does show which trends are to be expected in real systems. To model a range of polymer-surfactant interactions, various repulsions be(29) Groot, R. D.; Warren, P. B. J. Chem. Phys. 1997, 107, 4423. (30) Groot, R. D.; Madden, T. J. J. Chem. Phys. 1998, 108, 8713. (31) Groot, R. D.; Madden, T. J.; Tildesley, D. J. J. Chem. Phys. 1999, 110, 9739. (32) Kong, Y.; Manke, C. W.; Madden, W. G.; Schlijper, A. G. J. Chem. Phys. 1997, 107, 592. (33) Spenley, N. A. Europhys. Lett. 2000, 49, 534.
Table 1. Interaction Parameters Used in the Simulationa parameter set
χph
χpt
aph
apt
1 2 3 4 5 6 7
1.0 0.0 -3.0 -4.5 -4.5 -4.5 -4.5
-1.5 -1.5 -1.5 -1.5 0.0 2.0 18.0
28.3 25.0 15.2 10.3 10.3 10.3 10.3
20.0 20.0 20.0 20.0 25.0 31.5 83.8
a For each set the number of surfactant molecules present was varied from 0 to 50, unless stated otherwise. There are 57 runs in total.
tween the polymer beads and the surfactant tails and headgroups were studied. The following matrix of repulsion parameters was chosen
(2)
which captures at least qualitatively the relevant physics. Here (w, h, t, p) denote water, surfactant headgroup, surfactant tail, and polymer. From eq 2 it is clear that the headgroup has a strong affinity to water (χhw ) -3), which represents a strong hydration of the headgroup. Effectively this describes a property of ionic surfactant as this interaction strongly increases monomer solubility. Hence, although a real long-range electrostatic interaction is missing at least some of the characteristics are captured. It should be noted that the Flory-Huggins χ parameters used here are slightly different from convention, in that the interaction between particles of the same type can be characterized by a nonzero χ parameter. The χ parameters are calculated from29 χij ) 0.306 (aij - 25), rather than χij ∝ (aij - 1/2aii - 1/2ajj). The true parameters in the simulation are the repulsion parameters given in eq 2; the results are presented in χ parameters because of their historic use in polymer-surfactant chemistry. The parameters that have been varied are the repulsion parameters between polymer and surfactant headgroup, aph, and between polymer and surfactant tail, apt. The parameters of the modeled systems are summarized in Table 1. The parameter sets can be understood as follows: sets 1 and 2, attraction between polymer and surfactant tail, repulsion and strong repulsion between polymer and surfactant headgroup; sets 3 and 4, attraction between polymer and both head and tail; sets 5-7, attraction between polymer and surfactant head, repulsion between polymer and surfactant tail. The simulations comprised of one homopolymer (length L ) 50) in a box of size 10 × 10 × 10 Rc3, with Ns ) 0, 1, 2, 5, 10, 20, or 50 added surfactant molecules, unless stated otherwise. These molecules were immersed in a solvent of single beads adding up to a total of 3000 particles in each run. The time evolution of this system was calculated using the dissipative particle dynamics (DPD) method. This method is an improvement over ordinary molecular dynamics in that noise and friction are added in a pairwise additive way, such that the simulated system still represents a correct Gibbs-Boltzmann ensemble.29,34 These extra noise and friction forces lead to an exact inclusion of all hydrodynamic interactions.35 The coarsegraining of the model, the introduction of the pairwise (34) Espan˜ol, P.; Warren, P. B. Europhys. Lett. 1995, 30, 191. (35) Espan˜ol, P. Phys. Rev. E 1995, 52, 1734.
7496
Langmuir, Vol. 16, No. 19, 2000
Groot
friction, and the soft nature of the interaction potential all contribute to a very powerful simulation tool. Consequently, the accessible evolution times are much larger than in ordinary molecular dynamics. An estimate of the simulated time and length scales is given in the next section. 2.2 Analysis of Length and Time Scales. For the sake of clarity on the size of the time step and distance scale, an order of magnitude analysis is given here. The most important clock that determines how much time passes is the diffusion of micelles. In a simulation of 3000 DPD beads in total, with Ns ) 100 surfactant molecules (no polymer) in a volume V ) 1000 Rc3, four micelles are formed. Hence, the mean volume per micelle is πd3/6 ) 50/3000 × 1000 Rc3 ) 16.7 Rc3 from which a diameter d ≈ 3.2 ( 0.1 Rc follows, where Rc is the cutoff radius of the DPD interactions, the natural DPD unit of length. The error is obtained from the realization that the equilibrium number of micelles per 100 molecules must be between 3.5 and 4.5. Visual inspection of the simulated micelles corroborates this size. In real life the size of an SDS micelle obtained from surface force measurements7 is about 4 nm; hence, we have
Rc ≈ 1.25 ( 0.05 [nm]
(3)
According to Stokes-Einstein the diffusion constant of the micelles is (inserting the viscosity of water at room temperature)
Dmic )
4.3 × 10 kT ≈ 3πηd d
-19
[m2 s-1]
(4)
which works out at Dmic ≈ 10-10 m2 s-1. Although occasionally a molecule in the simulation is in the monomer state, the vast majority of the surfactant molecules is bound in micelles. Therefore the surfactant diffusion constant is equal to the micelle diffusion constant. In DPD units the surfactant diffusion constant is Ds ≈ 0.030 ( 0.003 Rc2 τ-1. Equating the diffusion constants in the simulation and in the physical system leads to
Dmic ≈ Ds ) 0.030 ( 0.003Rc2τ-1 ) 0.030 ( 0.003(d/d/)2τ-1 ) 4.3 × 10-19 d-1 [m2 s-1] (5) from which the time scale follows as
τ ) (7.0 ( 0.7) × 1016 d3 d-2 / ) 0.44 ( 0.08 [ns] (6) It follows that for the present polymer-surfactant simulation problem the elementary time steps taken (0.06τ) are some 25 ps. 2.3 Surfactant Chemical Potential. What we want is to describe a single polymer in an infinitely large surfactant solution. To calculate the bulk surfactant concentration that a simulated polymer-surfactant system is in equilibrium with, the first thing to establish is which quantities determine this equilibrium. This can be modeled by two different boxes: one containing polymer, surfactant, and water and the other containing only surfactant and water. A cartoon of the system that is considered is depicted in Figure 1. The total Gibbs free energy of the system comprising the sub-systems I and II is
G ) µpINpI + µsINsI + µwINwI + µsIINsII + µwIINwII (7)
Figure 1. Equilibrium between surfactant adsorbed at a polymer (I) and in bulk (II).
where I and II refer to boxes I and II and Np, Ns, and Nw are the number of polymer, surfactant, and water molecules present in the respective boxes. This thought experiment is an example of the Gibbs ensemble in which the total number of molecules in the two boxes is fixed, but molecules are allowed to move from box I to box II and vice versa.23 However, as an extra constraint, we impose that the total number of particles in each box is constant. The only allowed moves are swaps of a surfactant molecule from I to II and a simultaneous swap of an equal number of water beads from II to I and vice versa. The variation of the Gibbs free energy under these swaps is
δG ) (µsI - µsII)δNs - Ls(µwI - µwII)δNs
(8)
where Ls is the number of beads in a surfactant molecule. As in equilibrium δG ) 0, we find
µsI - LsµwI ) µsII - LsµwII ) µ
(9)
from which the proper chemical potential follows as
µ ) kT ln(Fs) - LskT ln(Fw) + ∆µs - Ls∆µw (10) Here Fs and Fw are the concentrations of surfactant and water molecules, respectively. Hence, to calculate the equilibrium concentration of surfactant in a system without polymer, we need to measure the excess chemical potentials of both water and surfactant. For this purpose the Widom insertion method is applied, which is particularly suitable for this model since the interaction potential is so soft. The practical application of this method to this model is detailed in Appendix A. As we only do virtual insertions to measure the water and surfactant chemical potentials, only one box is used in practice, and the chemical potentials in boxes I and II are calculated in separate runs. 3. Simulation Results 3.1 Polymer Endpoint Separation. The observable that is the easiest available from simulations is the polymer endpoint separation. Some typical results will be described here first, before the chemical potential simulations are discussed. It has been found that surfactants from all parameter sets studied show the same qualitative behavior as far as the polymer endpoint separation is concerned. For this reason we concentrate here on polymer-surfactant interactions of parameter set 2, given in Table 1. The simulations comprised of one homopolymer (length L ) 50) in a box of size 10 × 10 × 10, with Ns ) 0, 1, 2, 5, 10, 20, 33, 50, 73, 100, 200, 333, 500, and 800 added surfactant molecules. These were immersed in a solvent of single beads adding up to a total of 3000 particles
Polymer-Surfactant Aggregation
Langmuir, Vol. 16, No. 19, 2000 7497
micelles are discernible and the polymer coils from one micelle to another. This textbook state1-3 is alternated with a state where the polymer-surfactant complex forms a sausage where all surfactant molecules run across the polymer backbone collectively in a wavelike motion, see Figure 3b,c. This break-up of micelles is related to the strong attractive interaction between the polymer backbone and surfactant tails, as will be discussed in detail in the next section. To further analyze the system the endpoint distribution has been fitted to the scaling function37,38 Figure 2. Polymer endpoint separation as function of added surfactant in the simulation box. For small surfactant concentration all surfactant molecules are adsorbed on the polymer, hence Nagg ∝ φs for all but the last five points.
in each run. To equilibrate the systems and allow micelles to be formed, the systems were aged over 104 steps (0.25 µs), and then the polymer endpoint distribution function was averaged over the next 6 × 105 steps, or 16 µs. After 105 time steps (2.5 µs) the system with 800 surfactant molecules had formed a hexagonal phase, so this was discarded for further analysis. The results are shown in Figure 2. The endpoint separation has been scaled according to eq 3 to obtain a polymer size in nanometers, the error bars do not include the systematic error of this overall scale factor. This figure indicates a dramatic decrease in size of the polymer as the surfactant concentration increases, up to a certain point where precisely one micelle has formed at the polymer. From then on the polymer starts to swell again. At this point we mention that the radius of gyration of the polymer without surfactant is 3.56 ( 0.02 nm. If the solubility is slightly increased to χ ) 0.4, which is realistic for PEO, then RG ) 3.83 ( 0.02 nm is obtained. The experimental radius of gyration of PEO as reported by Cohen Stuart et al.36 follows the relation RG ≈ 0.0234 Mw0.565(0.018 where RG is in nanometers. Matching the polymer size to that of PEO leads to a corresponding polymer of molecular weight Mw ) 8000 ( 700. This number is given as a rough indication of the polymer size to which these simulations apply, but it is not to be understood that the results are restricted to PEO. Pictures of typical polymer conformations with 10 surfactant molecules added (less than one micelle) and 100 surfactant molecules added (more than one micelle) are shown in Figure 3. The conformations shown are some 100 and 300 ns old, respectively. In the Ns ) 10 system (Figure 3a) all surfactant molecules are already adsorbed on the polymer by 70 ns. What is observed in a movie of the Ns ) 100 simulation is that sometimes individual
ln(ψ(r)) ) a +
- 0.5 (1.026ν )ln(r) - br 1-ν
1/(1-ν)
(11)
where a and b are arbitrary fit parameters and ν is the swelling exponent. Examples of endpoint distribution functions are shown in Figure 4. On addition of surfactant the distribution first narrows (Ns ) 20) but for high surfactant concentration the polymer swells again. In Figure 5 the swelling exponent that is obtained this way is compared with the endpoint separation. The curves are qualitatively very similar. One can argue that for large amounts of surfactant the complex can be expected to behave as a stiff polymer, and therefore the simple scaling distribution as in eq 11 should not hold. However, the swelling exponent found this way should give a reasonable estimate of the swelling exponent on the relevant length scale, i.e., larger than a micelle. What this plot indicates is that a polymer that is initially marginally soluble (ν ) 0.5) goes through a coil-globule transition (ν < 0.4) when surfactant is added in a particular ratio. When yet more surfactant is added the polymer swells again, even more than a self-avoiding chain, ν ) 0.65. This should be contrasted to experimental observations. Chari et al.10 obtained the swelling exponent ν ) 0.65 for a saturated PEO/SDS system as we find here. Thompson et al.8 studied the thickness of an adsorbed layer of PEO and found an initial decrease, followed by a subsequent increase when the SDS concentration is increased, very similar to the present simulation results. Mel’nikov et al.13 report a more sudden coil-globule transition for DNA in a C16TAB solution and observe a size increase above 10-3 M surfactant concentration. Related to this, Shimabayashi et al.9 observe a minimum in the cloud point of hydroxypropylcellulose as function of the cationic surfactant cetylpyridinium chloride concentration, and Claesson et al.7 report a similar behavior for a cationic polymer interacting with SDS. A decrease
Figure 3. Conformations of polymer-surfactant complexes with 10 surfactant molecules present (a) and with 100 surfactant molecules present, for which two typical conformations are shown (b, c).
7498
Langmuir, Vol. 16, No. 19, 2000
Figure 4. Logarithm of endpoint distribution function for three amounts of surfactant.
Groot
Figure 6. Surfactant chemical potential without (b) and with ([) polymer. Most surfactants in these simulations are adsorbed in a single cluster of aggregation size Nagg ) 103Fs.
adsorbed on the polymer does not reflect the bulk surfactant concentration that the presently simulated systems are in equilibrium with. For a quantitative correspondence this is what we need to know. To this end we need to establish the chemical potential of the surfactant. This is the subject of the next section. 3.2 Surfactant Chemical Potential in PolymerSurfactant Complexes. 3.2.1 Hydrophobically Interacting Surfactants. In this section we concentrate on the chemical potential as defined in eq 10, again for surfactant of parameter set 2, given in Table 1. Some care must be taken in calculating the ideal part of the chemical potential when the number of surfactant molecules Ns is small. The ideal part for a simulation with Ns molecules, namely is23
µid ) kT ln((Ns + 1)/V) Figure 5. Swelling exponent (closed data points) compared with endpoint separation. For small surfactant concentration all surfactant molecules are adsorbed on the polymer, hence φs ∝ Nagg for all but the last five points.
in polymer size as seen in the present simulations is not the same as clouding and resolubilization (the collapsed polymer-surfactant complex is actually soluble for the present parameters) but the two phenomena are related for the following reason. When a polymer shrinks in size it follows that it sticks to the micelle. For certain polymer lengths, and surfactant and polymer concentrations, this may lead to micelle mediated polymer-polymer association and to the subsequent formation of an associative network. Simulations have shown that the formation of associative cross-links induces phase separation for sufficiently strong association.39 Therefore the present simulations are consistent with the observation of clouding, though much more work is necessary to establish the phase diagram for this model. At this point a qualitative correspondence between simulation and experiment has been established. However, the square-root dependence of the polymer size found as a function of surfactant concentration (or rather as a function of aggregation number) is quite unexpected. A rationale for this behavior is not easily found, but it should be realized that the number of surfactant molecules (36) Cohen Stuart, M. A.; Waaijen, F. H. W. H.; Cosgrove, T.; Vincent, B.; Crowley, T. L. Macromolecules 1984, 17, 1825. (37) Groot, R. D.; Bot, A.; Agterof, W. G. M. J. Chem. Phys. 1996, 104, 9202. (38) de Cloizeaux, J.; Jannink, G. Polymers in Solution, Their Modelling and Structure; Clarendon: Oxford, 1990. (39) Groot, R. D.; Agterof, W. G. M. J. Chem. Phys. 1994, 100, 1649.
(12)
because the Widom insertion method calculates the chemical potential of a system of Ns + 1 molecules as a perturbation from an Ns molecule system (see Appendix A). Therefore the concentration to be used in a ∆µ - ln(Fs) curve should be Fs ) (Ns + 1)/V rather than Fs ) Ns/V. This is the standard that has been adopted for the curves in Figure 6. In principle the difference between two chemical potentials can be measured by changing one molecule into another, but in the present case that would give unreliable statistics. When micelles are formed, there will be very few water molecules in the hydrophobic region (the micelle core), and this is the region that gives the largest contribution to the surfactant chemical potential on virtual insertion of a molecule. For this reason the water and surfactant chemical potentials have been measured separately in the same simulated system and are subtracted from each other afterward. Insertion of surfactants in a bath of water without polymer or surfactant is not very accurate using this method. Accuracy can be gained by changing water into surfactant, but this has not been pursued here. Figure 6 shows the chemical potential of surfactant in the presence of polymer ([) and without polymers present ([). The simulation results for insertion of surfactant without any previous surfactant present are indicated by the infinite dilution curves, as the excess chemical potentials measured in this way are in fact the excess chemical potentials at infinite dilution. We first concentrate on the system without polymer. When we draw a horizontal line at the minimum of the chemical potential (when a single micelle is present), we find the CMC of the pure surfactant at the intercept with the infinite dilution
Polymer-Surfactant Aggregation
Figure 7. Polymer endpoint separation as function of equilibrium bulk surfactant concentration. This is the surfactant concentration that a polymer-surfactant aggregate is in equilibrium with.
limit,16 as indicated in Figure 6. The CMC in this example is about 3 × 10-4. The unit of concentration used is Rc-3; using eq 3 this can be converted to moles/l by: 1 Rc-3 ≈ 0.85 M. Below the CMC, and at not too high concentrations above it, the surfactant molecules are clustered in a single micelle. So what this curve represents is the relative stability of surfactant micelles of different aggregation number, µ(Nagg) ≈ µ(103Fs), showing a minimum at a preferred micelle size. Similar results have been reported by Mackie et al.40 for a lattice model. The reason the curve bends upward at higher aggregation numbers must be attributed to the difference in head-head and tail-tail repulsion. This difference results in a spontaneaous curvature of the surfactant-water interface toward the tail side, and consequently to a free energy increase for large (i.e., flat) micelles. This in turn leads to a preference of small spherical micelles over elongated rods. For the system with polymer present, the chemical potential is a continuously rising function for 0 < Ns < 50. The initial slope 1/n ≈ 1/4 suggests that aggregates of 4 molecules bind to the polymer rather than complete micelles. Consequently, when we intercept the Ns ) 0 chemical potential with the infinite dilution curve without polymer, we find the critical (bulk) surfactant concentration at which surfactant starts to adsorb on the polymer: the CAC. In this example the CAC is about 1 × 10-4, see Figure 6. In the same way, for every polymer-surfactant system we can determine the equilibrium bulk concentration, by determining the horizontal intercept of the actual chemical potential with the infinite dilution limit. When this is combined with the polymer endpoint separation results (Figure 2), we find the curve shown in Figure 7. When compared to the curve in Figure 2, it is obvious that the square root dependence on the surfactant concentration (which is roughly the aggregation number, up to a factor V ) 103) is only pertinent when a small NVT ensemble is considered, but not in an open ensemble. In Figure 7 data points above the CMC are not included, because there the chemical potential is nearly independent of surfactant concentration, and therefore the corresponding bulk concentration is undetermined. To obtain this equilibrium concentration another simulation method must be used, viz. the Gibbs ensemble.41 In this method the two boxes shown in Figure 1 are simulated explicitly, (40) Mackie, A. D.; Panagiotopoulos, A. Z.; Szleifer, I. Langmuir 1997, 13, 5022. (41) Willemsen, S. M.; Vlugt, T. J. H.; Hoefsloot, H. C. J.; Smit, B. J. Comput. Phys. 1998, 147, 507.
Langmuir, Vol. 16, No. 19, 2000 7499
Figure 8. Surfactant chemical potential for systems 1, 3, and 6 (see Table 1) and for surfactant without any polymer present. As most surfactants are aggregated into one cluster, the concentration axis roughly denotes the clustersize via Nagg ) 103Fs.
i.e., the swaps of molecules between the boxes are real swaps rather than just virtual insertions. 3.2.2 General Polymer-Surfactant Interaction. In Figure 8 results are shown for some of the parameter sets indicated in Table 1. For parameter set 1 the same slope of µ against ln(Fs) is reproduced as for set 2, hence again aggregates of some four surfactants bind to the polymer. When the headgroup repulsion is changed from a net repulsion into a net attraction, the curve shifts down (parameter set 3) but still has the same initial slope at small surfactant concentration. However, when the interaction of the polymer with the tail is turned repulsive (χpt is increased from -1.5 to +2.0) a dramatic change is observed: for parameter set 6 the chemical potential first decreases with surfactant concentration, passes through a minimum, and then increases again. This behavior signals micelle adsorption. When for parameter set 6 the surfactant concentration in the bulk is slowly increased, the adsorbed amount will follow the infinite dilution limit in the presence of polymer, up to the point where on average 0.3 molecules are adsorbed. Here the chemical potential equals the minimum chemical potential at which adsorbed micelles can exist. At this point the adsorbed amount jumps up to the value corresponding to one adsorbed micelle, see Figure 9. Such behavior is predicted to occur when a strong interaction between surfactant headgroup and polymer is pertinent, for instance for the binding of a cationic surfactant to an anionic polymer. The loop structure seen for parameter set 6 will not appear in an open system, i.e., in a large system with many polymers present, but it tells us that adsorption on each polymer is quantized by multiples of discrete micelles. It also implies that a range of surfactant concentrations must exist where extended polymer conformations coexist with collapsed ones. This is precisely what is found by Mel’nikov et al.13 for the interaction of DNA (a strong polyanionic) with CTAB (a cationic surfactant). This result also implies that when the total surfactant concentration is continuously increased in a system with many polymers present, the bulk surfactant concentration will remain constant at the value corresponding to the CAC, until a micelle is adsorbed on every polymer in the system. At the CAC, the system contains a mixture of polymers without any surfactant molecules adsorbed and polymers with a complete micelle adsorbed. These micelles are smaller than the micelles that are formed in the bulk.
7500
Langmuir, Vol. 16, No. 19, 2000
Groot
Figure 9. Amount of surfactant adsorbed to the polymer as function of the chemical potential, for systems 1, 3 and 6 (see Table 1).
Figure 10. Chemical potential of surfactant that adsorbs on a polymer primarily via its tail (χpt ) -1.5) for various headgroup repulsions (χph ) 1, 0, -3, and -4.5) from top to bottom. As most surfactants are aggregated into one cluster, the concentration axis roughly denotes the clustersize via Nagg ) 103Fs.
This can be seen already in Figure 8: the micelle size corresponds to the number of surfactant molecules in the simulation where the chemical potential is at its minimum. For the polymer-free system this is Nagg ) 30.6 ( 0.6 molecules whereas for parameter set 6 this is Nagg ) 7.9 ( 0.8 molecules. The absolute numbers are lower than for, e.g., SDS because of the large head-head repulsion chosen, but the qualitative influence of the polymer on the micelle size is apparent. To highlight the effects of polymer-surfactant association via the surfactant headgroup and via the surfactant tail, the interactions were varied systematically. In Figure 10 the polymer-surfactant tail interaction is kept constant at a net attractive value, χpt ) -1.5, whereas the polymersurfactant headgroup interaction is varied from χph ) 1 (top curve) to χph ) -4.5 (lower curve). The fit curves have the functional form
µ ) µCAC + a ln(n) + bn
(13)
where n ) Ns + 1 ) 103Fs is the number of (real or virtual) surfactant molecules in the system. The constant a has the interpretation of an inverse cluster size, and the constant b is the second virial coefficient of the adsorbed surfactant. The numerical value of a is independent of
Figure 11. Second virial coefficient of surfactant adsorbed on polymer, see eq 13, correlated with the chemical potential difference of CMC and CAC.
Figure 12. Chemical potential of a surfactant that interacts with a polymer primarily via its headgroup (χph ) -4.5). The polymer-tail interaction decreases from top to bottom (χpt ) 18, 2, 0, and -1.5). As most surfactants are aggregated into one cluster, the concentration axis roughly denotes the clustersize via Nagg ) 103Fs.
the polymer-headgroup interaction; for the four systems it is 0.23, 0.22, 0.24, and 0.27, respectively, their spread being well within the uncertainty of the fit parameters ((0.06) and apparently equal to 1/4. All curves converge by the chemical potential of the bulk CMC (see Figure 10). Parameter b correlates linearly with µCAC, or alternatively, since the chemical potential at the CMC is µCMC ) -3.19, it can be corrrelated with µCMC - µCAC, see Figure 11. This indicates that the surfactant fills up the available adsorption area on the polymer continuously, and the chemical potential simply has to increase to the CMC value when there is no more room left on the polymer. A completely different picture emerges when the surfactant interacts primarily with the polymer via headgroup interactions. This would in practice be the case for an anionic polymer such as sodium (carboxymethyl)cellulose interacting with a cationic surfactant like C14TAB.5,24 In Figure 12 the χ parameter between polymer and surfactant headgroup is kept constant at χph ) -4.5, while the interaction with the surfactant tails is varied as χpt ) 18, 2, 0, -1.5. These results indicate a distinct minimum in the chemical potential, the position of which depends on the interaction of polymer with the tail. For a very strong repulsion, the micelles hardly bind to the polymer, so the tail repulsion overrules the headgroup attraction. For a moderate tail repulsion the micelles do
Polymer-Surfactant Aggregation
Figure 13. Schematic phase diagram showing the parameters of the performed polymer-surfactant simulations and showing which type of complex was subsequently found. The parameters are χpt, the interaction energy between polymer and surfactant tail, and χph, the interaction energy between polymer and surfactant headgroup.
bind to the polymer via the headgroup, but the size of the micelles is distinctly smaller than what is found in the bulk. The size of the micelles appears to decrease continuously as the polymer-tail repulsion is decreased. The minima in the chemical potential for the systems in Figure 12 occur at Nagg ) 25 ( 3, Nagg ) 7.9 ( 0.8, Nagg ) 3.8 ( 0.4, and Nagg ) 1.3 ( 0.4, respectively. For the last curve the micelle size has decreased to one molecule, and we are back to a continuous adsorption mode. It should be mentioned that for the parameter sets shown in Figure 12 the chemical potential cannot be fitted with the functional form given in eq 13. The behavior of this binding mode is qualitatively different from the previous type of thermodynamic behavior. To map out roughly when each of the binding modes is pertinent, a number of short runs has been performed throughout the (χpt, χph) parameter space, for 1 polymer and 30 surfactant molecules. The result is shown in the “phase diagram” in Figure 13. This diagram should be considered as a qualitative picture, capturing the relevant physical phenomena. To predict such a diagram for a particular surfactant requires a careful tuning of the parameters that have been kept fixed here, viz. the headgroupheadgroup and the headgroup-tail interactions. 4. Discussion and Conclusions Polymer-surfactant aggregation has been studied using the DPD simulation technique. The polymer was represented by a string of 50 soft spheres, and the surfactant was modeled by dimers of two soft spheres, immersed in a solvent of soft spheres. This polymer roughly corresponds to a real polymer of molecular weight in the order of Mw ≈ 8000. The interaction between polymer and surfactant has been varied systematically. This has led to the following observations: (1) Polymers shrink in size when small amounts of binding surfactants are added. The size passes through a minimum as the surfactant concentration is increased, and by the time the bulk CMC is reached the polymer swells beyond its original size. (2) There are two distinct modes of surfactant adsorption on polymers, depending on the interaction mechanism. These are (a) continuous adsorption of surfactant on the polymer and (b) surfactant adsorption by discrete micelles.
Langmuir, Vol. 16, No. 19, 2000 7501
(3) If the interaction between polymer and surfactant is dominated by the surfactant tails (i.e., when the tails tend to adsorb on the polymer) then mechanism a is pertinent. The tail interaction dominates when it is net attractive, even if there is also a net attraction via the headgroups. (4) If the binding between polymer and surfactant is dominated by the surfactant headgroups, mechanism b is pertinent. Hence headgroup interaction leads to micellar adsorption. (5) For each headgroup attraction there is a tail attraction where the mode of adsorption crosses over from micelle adsorption to molecular bottlebrush conformations. This crossover is comparable to a critical point in the condensation of gases to liquids. A schematic phase diagram indicating the various adsorption modes is shown in Figure 13. (6) For headgroup dominated polymer-surfactant interaction (e.g., anionic-cationic), the micelle size decreases as the polymer contains more hydrophobic parts, e.g., in polysoaps. An experimentally observable measure of polymersurfactant interaction is the endpoint separation of an isolated polymer in the presence of added surfactant. On increasing surfactant concentration the polymer is found to reduce in size. The reason the polymer decreases in size is observed in the simulation: the polymer coils around the micelle when only one micelle is formed. When the surfactant concentration is increased beyond this point the polymer size increases again. The reason for this increase is that now two or more micelles are adsorbed, so that the polymer does not have to fold completely around one micelle in order to maximize the number of polymersurfactant contacts. From these simulations we learn that there are two different classes of polymer-surfactant interaction with respect to the structures that are formed. An important question is whether the assumption of spherical micelle adsorption on polymers in a thermodynamic theory1,2 or in a simulation24 is justified. It now appears that this is justified when the interaction takes place via the headgroups of the surfactant, which is the case for an anionic polyelectrolyte interacting with a cationic surfactant. However, is not necessarily correct to assume that the micelle size is independent of the presence of polyelectrolyte. If hydrophobic interaction between polymer and surfactant is the dominant binding mechanism, then the assumption of spherical micelle adsorption on the polymer is incorrect. In that case, molecular bottlebrushes or adsorbed micelles alternating with the bottlebrush conformation appear. The reason most polymer-surfactant systems reported in the literature aggregate via discrete micelles is that predominantly systems have been studied that interact via the surfactant headgroup. The qualitative agreement between the available experimental data and the presently studied dumbbell surfactant model warrants a further study to map chemically specific surfactant systems onto the dumbbell model. Appendix A. The Widom Insertion Method The most straightforward way to obtain the excess chemical potential ∆µ of a molecule from a simulation is to use the Widom insertion method.23 An alternative method to calculate the chemical potential is to use the overlapping distribution method,41 but due to the softness of the interaction potential this is not necessary for the present model. In the Widom insertion method the ratio
7502
Langmuir, Vol. 16, No. 19, 2000
Groot
of the partition functions containing N + 1 and N molecules is calculated. This ratio is the excess chemical potential, which for a dimer follows exactly from
∫
∆µ ) -kT ln〈e-u(r1)/kT e-u(r2)/kTd3r2〉 - ideal chain (A1) where ideal chain is the logarithm of the configuration integral of an ideal chain. For an ideal chain, ∆µ should vanish exactly. This defines a state of reference. For flexible chains we can choose any ideal reference chain because the partition function of the reference chain is subtracted again. An obvious choice is to take a chain that is connected by such a spring force that the beadbead distribution function is equal to (or as close as possible to) the distribution of the chain that is to be studied. Let the interaction potential between neighboring beads of the reference chain be denoted by uref(r), and let the reference distribution be Pref(r) ) exp(-uref(r)/kT). The excess chemical potential of the dimer is then
∆µ ) -kT ln〈e
∫e
-u(r1)/kT
-uref(r2-r1)/kT-[u(r2)-uref(r2-r1)]/kT 3
d r2〉
- ideal chain
∫
) -kT ln〈e-u(r1)/kT e-[u(r2)-u
ref(r -r )]/kT 2 1
×
∫
Pref(r2 - r1)d3r2/ Pref(r2 - r1)d3r2〉 (A2) ) -kT ln〈e-u(r1)/kT〈e-[u(r2)-u
ref(r -r )]/kT 2 1
〉ref〉
where the last average is taken for points r2 chosen from the reference distribution, under the condition that the bead to which it is connected is located at r1. This average is also known as the Rosenbluth weight for point 2. The excess chemical potential of a chain longer than two beads is found analogously, by taking the average of the product
of all Rosenbluth weights from bead 2 up to bead L, multiplied by the Boltzmann factor for bead 1. After a Rosenbluth weight is determined via Monte Carlo averaging, one point from the equilibrium distribution is chosen (not the reference distribution), which then serves as a fixed point for the determination of the next Rosenbluth weight. This way chains are “grown” into the ensemble according to the correct equilibrium distribution. Because the chains we try to insert are biased toward the correct distribution, the method is known as biased sampling Monte Carlo. The chosen reference segment distribution function has the following form
Pref(r) ) exp(-1/2Cr2/kT - 1/2 aref(1 - r)2θ(1 - r)) (A3) where C ) 4, θ is the step function, and aref is chosen according to
aref ) 0.75(aj,j+1 - 18)
(A4)
where aj,j+1 is the repulsion parameter between next beads of the chain that is inserted. This choice is prudent if the repulsion parameter between the water molecules is aww ) 25. For repulsion parameters in the range 20 < aj,j+1 < 70, the reference distribution is a very accurate representation of the real bead-bead distribution as obtained from the simulation. To check the method, the excess pressure was plotted as function of the excess chemical potential for melts of chain length L ) 1-3, for segment densities 0.5 e FL e 6, for repulsion parameter a ) 25. A fit curve through the points was used to take the derivative dp/dµ, which according to the Gibbs-Duhem relation should equal the density. The deviation was found to be of the order 0.3%. Since the Gibbs-Duhem equation is satisfied, we conclude that the chemical potential is sampled correctly. LA000010D