Interfacial properties and phase transitions in ternary symmetric homopolymer–copolymer blends: A dissipative particle dynamics study

Polymer 54 (2013) 2146e2157

Contents lists available at SciVerse ScienceDirect

Polymer journal homepage: www.elsevier.com/locate/polymer

Interfacial properties and phase transitions in ternary symmetric homopolymerecopolymer blends: A dissipative particle dynamics study Zhiqiang Bai a, b, c, Hongxia Guo a, b, c, * a

Beijing National Laboratory for Molecular Sciences, Institute of Chemistry, Chinese Academy of Sciences, Beijing 100190, China Joint Laboratory of Polymer Sciences and Materials, Institute of Chemistry, Chinese Academy of Sciences, Beijing 100190, China c State Key Laboratory of Polymer Physics and Chemistry, Institute of Chemistry, Chinese Academy of Sciences, Beijing 100190, China b

a r t i c l e i n f o

a b s t r a c t

Article history: Received 26 October 2012 Received in revised form 7 January 2013 Accepted 9 February 2013 Available online 15 February 2013

We use dissipative particle dynamics simulations to study the interfacial properties and their relevance to the phase transitions in ternary symmetric blends. For the blend in the lamellar (LAM) state, by analyzing the distribution of inter-bilayer distances as well as the position fluctuations and correlations of bilayers, we found that the “Discrete Harmonic” theory can describe the fluctuations of the LAM phase in ternary blends in a satisfactory way when the LAM phase is not too close to the phase transition to the bicontinuous microemulsion (BmE). In particular, as the LAM stack is swollen increasingly with the addition of homopolymers, the bilayers in the LAM phase experience a crossover from the single coherent fluctuation mode to the coexisting of coherent and incoherent fluctuation modes wherein the incoherent free membrane regime starts at larger length scales and becomes more pronounced. Meanwhile, the in-plane correlation length increases, the bending modulus and compressibility modulus decrease, indicating that as the more homopolymers are added, fluctuations between adjacent bilayers become less coherent, bilayers become more flexible, and the highly swollen LAM system becomes more susceptible to thermal fluctuations. Therefore, the phase transition from LAM to BmE is expected to occur when the bending modulus turns out to be within KBT or when the persistent length is no larger than the lamellar spacing. Also the two criteria yield the same predictions for the transition point. While for the macrophase-separated (2P) phase, we found that Helfrich model gives a good description of the fluctuation spectrum of the copolymer monolayer at long wavelengths if the 2P systems are not too close to the phase transition to BmE. With the addition of copolymers the interfacial tension reduces but the bending modulus increases. Moreover, the phase transition of 2P to BmE coincides with the saturation of the interface between homopolymer phases, i.e. attaining a vanishing interfacial tension. In a word, our results justify the notion that there exists a direct correlation between peculiar changes in interfacial properties and the phase transitions from 2P or LAM to the BmE state. Ó 2013 Elsevier Ltd. All rights reserved.

Keywords: Symmetric ternary blend Interfacial property Phase transition

1. Introduction Blending polymers provides plentiful opportunities to create new materials with tailored properties [1e3]. However, in most cases, the miscibility between different polymers is rather low, polymeric blends are therefore unstable and tend to macrophase separate into their own domains, which deteriorates the mechanical properties of polymeric materials. Generally, polymeric blends are deemed as an assembly of the interfaces, and thus understanding and controlling the properties of these internal interfaces

* Corresponding author. Beijing National Laboratory for Molecular Sciences, Institute of Chemistry, Chinese Academy of Sciences, Beijing 100190, China. E-mail address: [email protected] (H. Guo). 0032-3861/$ e see front matter Ó 2013 Elsevier Ltd. All rights reserved. http://dx.doi.org/10.1016/j.polymer.2013.02.011

is curial for tailoring the material properties for practical applications. It is well established that adding copolymers, e.g. block copolymers and random copolymers, is an effective approach to stabilize the immiscible polymers, as copolymers are preferentially absorbed at the interfaces, reduce the interfacial tension and hence the driving force toward the macrophase separation [4e10]. Furthermore, copolymers can be utilized to direct the morphologies of polymeric blends and improve their properties. By the addition of substantial amounts of copolymers, the macrophase separation of polymer blends is prevented and new types of disordered or ordered microphase-separated structures, such as microemulsions and a variety of ordered mesophases are formed [11e14]. Microemulsion can take on two main morphologies: droplet and bicontinuous. The bicontinuous microemulsion (denoted as BmE), in which the blend is stabilized by the thermal

Z. Bai, H. Guo / Polymer 54 (2013) 2146e2157

fluctuations of its internal interface, contains an intriguing morphology composed of an extensive amount of undulating interfaces with near-zero spontaneous curvature [11,15]. Because the two homopolymers in BmE form continuous and randomly intertwined domains in the entire system, BmEs are receiving ever-increasing attention and are extensively utilized for preparing composites with improved mechanical or transport properties [16e22], such as interpenetrated elastomers, tunable porous membranes [19] or as templates for nanoporous materials [17,18]. It is well accepted that the formation of BmE depends crucially on the properties of its internal interface [23e29]. Additionally, adding sufficiently higher amounts of copolymers has a propensity to create a wide range of copolymer-rich ordered mesophases with many internal interfaces [12,13]. In such a case, knowing the interfacial properties is particularly important for understanding and consequently for tailoring the properties of these polymeric materials. Recently, the promise of tuning the morphology of the immiscible polymers by using copolymers has driven extensive studies into the thermodynamics of these ternary blends. As the simplest case, a symmetric ternary system, which contains symmetric AeB block copolymers and equal amounts of two homopolymers A and B of same lengths, is often used for fundamental phase behavior studies [11,15,30e33]. In a typical phase diagram, the blend is in a disordered state at high temperature. While, at low temperature, a macrophase-separated phase (denoted by 2P) and a lamellar phase (denoted by LAM) could form on the homopolymerrich side and on the copolymer-rich side, respectively. Between the 2P and LAM phases, a bicontinuous microemulsion phase could form in appropriate conditions. In comparison with the large research effort into the phase behavior of A/B/A-B ternary blends, the interfacial properties of such blends, especially for the interfaces that are part of a mesophase (i.e. lamellar phase in symmetric ternary blend) has received comparatively less attention. In principle, the interfacial tension g and the elastic constants (i.e. the bending modulus k and the interlayer compressibility modulus B) are key parameters for the interfacial properties. g is a measure of the energy cost of increasing the area of an interface by one unit [34], while k describes the energy cost of a flat interface toward bending [35] and B is directly related to interlayer interactions (hydration force and van der Waals forces) in lamellar phase [36,37]. For example, the BmE is characterized by a near-zero interfacial tension of g ¼ 0 and a low bending modulus of k  KBT [26]. Actually, these features not only allow a BmE to accommodate an extensive amount of internal interfaces but also permit thermal fluctuations to stabilize it [11,26,38e40]. Physically, it has been hypothesized [11] that the phase transition of 2P to BmE corresponds to the saturation of the A/B interface with A-B diblock copolymers and a vanishing interfacial tension, while the phase transition from LAM phase to BmE coincides with the flexibility of an isolated copolymer-laden interface exceeds a threshold value (i.e. k attaining a value that is comparable to KBT). Unfortunately, relatively little work has been done to systematically investigate this conjuncture. Particularly, despite the ability of block copolymers to reduce the interfacial tension between the coexisting homopolymer phases is well documented [5e10], their effect on the flexibility of the interface that separate the homopolymer phases is not completely understood [23e29]. Also, as just mentioned, little attention has been paid to the interfacial properties of the LAM phase in ternary blends. In light of the great interesting of ternary blends in both natural and engineering sciences, it is instrumental to make systematic investigations on the interfacial properties in both 2P and LAM phases and their relevance to the phase transitions both from 2P to BmE and from LAM to BmE.

2147

With the rapid development in computational power, computer simulation has become a promising tool to study the thermodynamic, dynamic, structural and interfacial properties of the polymers. However, the long-time nature and metastability involved in these ternary polymeric materials limits the straightforward application of atomic-level molecular dynamics (MD) simulations in studying their phase behavior and interfacial properties, since realistic atomic-based MD simulations are too expensive for the investigation of these interesting phenomena that lie within the mesoscopic spatio-temporal scale. As an alternative to the MD simulation, dissipative particle dynamics (DPD) [41,42], which is a relatively new coarse-grained mesoscopic simulation method, offers a computationally more efficient way of investigating the phase behavior and the interfacial properties of multi-component polymeric systems. Unlike MD simulations, DPD involves the use of soft repulsive interactions and a momentum conserving thermostat, and thus allows one to access large systems on a long time scale at a reasonable computational cost. Currently, DPD has been extensively applied to study the self-assembly process [43e47] or phase behavior and transition of various complex fluids [48e59], the interfacial properties of a monolayer [35,60] or a bilayer [61], and the fusion and fission of membranes or vesicles [62,63], etc. Given the instinct strength of the DPD method and the fundamental problem we considered, we adopt the DPD simulation technique and study the interfacial properties of ternary blends until approaching the phase transitions from 2P to BmE and from LAM to BmE and try to shed more light on the correlation between these phase transitions and peculiar changes in interfacial properties. To elucidate the basic mechanism of the behavior and interfacial character of such systems, a DPD simulation with simplified or generic coarse-grained (CG) models is expected to be sufficient to provide useful information. On the one hand, these simple models involve much fewer degrees of freedom and softer interactions than the atomistic models, and thus allow for the simulations at considerably large distance and time scales. On the other hand, these simple models only keep the most salient properties of the molecules of investigation (i.e. excluded volume, connectivity and the tendency to un-mix), therefore the correlation between changes in interfacial properties and the proximity to phase transitions from 2P to BmE and from LAM to BmE we draw from simulations is held for whole classes of ternary polymer systems. Actually, such generic CG models have already been used in DPD simulations of the phase separation of polymer blends and the microphase separation of block copolymers and led to many important results [48,49,51,55e58]. Like in MD simulations, the evaluation of the interfacial tension g in DPD simulations is relatively straightforward. One can either calculate g from the difference between the normal and tangential pressure tensor components across the interfaces that separate the macrophases or that are part of a mesophase [64], or estimate g from the spectral intensities of undulating modes if the interfaces in these systems are not highly curved [35,60,65]. Actually, the use of the undulation spectrum implies that these interfaces display obvious long-wavelength fluctuations and thereby their free energies can be well represented by the Helfrich’s curvature model in terms of interfacial tension and bending rigidity [66]. As a note of cause, this method can only be applied when the interfacial fluctuations are gentle and hence cannot be used to explore the larger fluctuations beyond thermally accessible amplitude or fluctuations of a 2P phase too close to the phase transition to the BmE. Specifically, several studies found that the interfacial tension obtained via the fluctuation spectrum is in good agreement with the value computed from the pressure tensor [35,60]. Meanwhile, the bending modulus k for the mono- or bi-layer system can be also straightforward extracted from the undulation spectrum by

2148

Z. Bai, H. Guo / Polymer 54 (2013) 2146e2157

monitoring the height fluctuations in the amphiphilic (surfactant or copolymer) monolayer between coexisting phases [25,35,38,60,67] or in the middle sheet of a single membrane bilayer [65,68]. Although this method as a means to extract k has now become commonplace, understanding of the effect of molecular structure and composition of amphiphilies on the k of an amphiphilic monolayer or a lipid bilayer is far from complete [25,35,38,60]. For example, there are competing opinions regarding the dependence of k on the amphiphilic interface coverage [25,35,38,60]. Therefore, a further systematic study of the variation of interfacial properties of an amphiphilic monolayer in the 2P phase until the amphiphilic density at the interface approaches to that of zero interfacial tension will help to resolve this issue. While, for the multilayers such as those appearing in an ordered lamellar morphology, it is exceedingly difficult to compute k in much the same way as just mentioned for single layer systems. As the lamellar phase exhibits a stack of bilayers, it is a difficult task to identify specific layers, follow those layers in time, and measure their fluctuation amplitude. Additionally, the equilibrium lamellar phase is tensionless and exhibits the meso-scale thermal undulations due to the flexibility, which give rise to an entropic repulsion between neighboring monolayers and thus to a swelling of lamellar phase. However, the undulations are restricted due to the presence of other membranes, so that additional factor for keeping the lamellar spacing to its equilibrium value such as the compressibility modulus B has to be taken into account. Indeed, up to now little work has been done about the undulations and thickness fluctuations in the lamellar phase of ternary blends. Recently, Loison et al. [69] have pointed out that the discrete harmonic model, a mesoscopic theory, which is used to measure the bending modulus Kc and the compressibility B for stacked membranes, can describe the fluctuations in lamellar stacks of a binary amphiphile-solvent mixture in a satisfactory way. Moreover, by combining the phenomenological parameter Kc/B with the analyses of other quantities, they showed it possible to determine the elastic constants of Kc and B separately [69]. In this paper, we extend this type of study to investigate the thermal fluctuations of the LAM phase formed in the symmetric ternary blend and extract their elastic constants and furthermore explore the physics of the transition from LAM to BmE. The paper is organized as follows: In the next section, we recapitulate the main points of the elastic theories for measurements of interfacial properties in single layer systems and in the LAM phase. In Section 3, we describe the simulation method, model, and simulation details. After that, we present and discuss our results in Section 4. Here, we first discuss briefly the phase behavior of our model. Then by analyzing the fluctuations of bilayers in the lamellar phase and of amphiphilic monolayers in the 2P phase in details, we study the behavior of the interfacial tension g and the elastic constants upon approaching the phase transitions from LAM to BmE and from 2P to BmE and try to shed more light on the topic of the correlation of the peculiar changes in interfacial properties to these phase transitions. Finally, we close with a summary of our findings in Section 5. 2. Theoretical backgrounds 2.1. Elasticity of an amphiphilic monolayer between coexisting phases On a mesoscopic length scale the properties of interfaces between coexisting phases are well described by the Helfrich theory [66,70], which regards the amphiphilic monolayers in the 2P phase as single smoothly undulating surfaces with free energies of

Z f ¼

i h k dA g þ ðc1 þ c2  2c0 Þ2 þ kc1 c2 2

(1)

where c1 and c2 are the principal curvatures. The four coefficients, g, c0, k and k, are denoted as the interfacial tension, spontaneous curvature, bending modulus, and Gaussian modulus, respectively. In our study, c0 vanishes because amphiphilic monolayers in symmetric ternary systems are symmetric with respect to the two macrophases. Additionally, in Eq. (1) we neglect the term involving the Gaussian modulus, since in our study the topology does not change and kdoes not affect the energy fluctuations in a monolayer [60]. We further suppose that the interfacial fluctuations are gentle, then Eq. (1) can be rewritten in terms of a height fluctuation function h(x,y), which measures the displacement of the interface from its mean position [66],

Z f ¼ A

!2 # "   g vh vh 2 k v2 h v2 h dxdy þ þ þ 2 vx vy 2 vx2 vy2

(2)

Note that the free energy is expanded to the second order h. Via a Fourier transform, h(x,y) can be written in terms of wavelength dependent undulation modes and

~f ¼ gq2 hðqÞ ~ 2 þ k q4 hðqÞ ~ 2 2 2

(3)

~ where q ¼ (qx,qy) is the lateral wave-vector and h(q) is the Fourier component of h(x,y). Apparently, the last term illustrates that in the Helfrich model the bending free energy incurred by every undulation mode is directly proportional to the bending modulus k and the square of the fluctuation amplitude. Then, according to the equipartition theorem, we can relate an average energy per mode to the spectral intensities of each undulation mode, such that [35,65]

D  2 E 1 K T ~  ¼ B gq2 þ kq4 h q  A

(4)

where KB is Boltzmann’s constant, T is the temperature, and A is the area of the monolayer. Hence, from the simulation, we can obtain 2 ~ Aiby monitoring the interfacial structure factor SðqÞ ¼ hjhðqÞj h(x,y), and then estimate the interfacial tension g and the bending modulus k by fitting the resulting undulation spectrum at small values of the wave vector q to Eq. (4) [35,60,65]. Recently, this method has been used to compute the bending modulus or the interfacial tension of many typical interfaces, such as surfactant monolayers [35,38,48,60,67], copolymer monolayers [25], biological membranes [65,71], and its feasibility has been well demonstrated. 2.2. Elasticity of the lamellar phase Generally, the Helfrich’s curvature model of Eq. (1) is often used as the starting point for more complex elastic theories of a single bilayer and a lamellar stack by incorporating additional factors such as the membrane thickness and interlayer coupling [72]. As the classical mesoscopic theory for fluctuations in stacked membranes, the “Discrete Harmonic” model (DH) [73,74], describes the lamellae as a discrete set of two-dimensional fluctuating layers, stacked in the z direction with an average interlayer distance of d and extended continuously in the (x,y) plane. The position of the nth layer can be parameterized by a unique height function Zn(x,y) and the fluctuation of the height of the nth layer from its mean position is characterized by the local displacement hn(x,y) ¼ Zn(x,y)  nd. The lamellar phase at equilibrium is tensionless (i.e. the interfacial tension of stacked membranes formed in binary [69] and quaternary [74] mixtures almost vanishes) and the contribution of the Gaussian curvature is a constant for such a stack with fixed

Z. Bai, H. Guo / Polymer 54 (2013) 2146e2157

topology [75]. For the simplest case, wherein only the interactions between adjacent layers are taken into consideration, the free energy of a lamellar stack with zero spontaneous curvature is [69,76]

f ¼

N1 X n¼0

!2 ) ( Kc v2 hn v2 hn B 2 þ þ  h Þ ðh n nþ1 2 vx2 2 vy2

Z

dxdy A

(5)

where N is layer number. The two terms on the right describe the energy cost due to local layer deformations and local deviations from the average interlayer distance. Thus, the elasticity of the stacked membranes is characterized by two elastic constants, bending modulus Kc and compressibility modulus B. Accordingly, the in-plane correlation length x for the thermal fluctuations of the lamellar phase is given as x ¼ (Kc/B)1/4. Conforming to Loison et al. the analysis of the fluctuations of hn(x,y) in a lamellar stack of bilayers, is performed in both Fourier and real spaces [69]. As addressed by Loison et al. the Fourier modes decouple in Eq. (5) and the equipartition principle applies [69]. Hence two types of Fourier transformations, including a discrete transformation in the z direction and continuous transformation in the x and y direction, are taken as below

hðqt ; qz Þ ¼ Z hn ðqt Þ ¼

X n

hn ðqt Þeiqz nd

drhn ðrÞ eIqt r ¼

A

(6) 1 X hðqt ; qz Þeþiqz nd N qz

dhn ðx; yÞ2 ¼_

dhn ðrÞ2

NLx Ly KB T 2B½1  cosðqz dÞ þ Kc q4t

(8)

X   1 X iq nd   2 1 N1 sn qt ¼_ 2 e z < h qt ; qz  > ¼ N j¼0 N qz < hj ðqt Þ$hnþj ðqt Þ* >

(9)

Apparently, the quantity s0 is the auto-correlation fluctuation spectra describing correlations within a single layer, whereas sn(n > 0) is the cross-correlation fluctuation spectra characterizing correlations between bilayers. When N / N (i.e. an infinitely thick stack of layers) and inserting Eq. (8), we obtain Ref. [69] N/N Lx Ly KB T Kc q4t



4 1=2 1þ H



H 1 pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi n N/N sn ðqt Þ ¼ s0 ðqt Þ 1 þ  HðH þ 4Þ 2 2

1   1 NX hnþj ðx; yÞ  hj ð0; 0Þ2 N j¼0

(10)

(11)

(12)

Again in conformity with Loison et al. this quantity can be calculated from back transforming the sn(qt) into the real space (x,y) [36,69,74,76],

D

Here hi denotes thermal averages. Because of the large statistical error of simulation results forhjhðqt ; qz Þ2 ji, extracting the elastic constants Kc and B by a direct comparison with Eq. (8) is impossible [69]. As an alternative method to determine these k elastic constants, the trans-bilayer structure factor is considered, which describes correlations between bilayer positions in different bilayers and is defined as [69].

s0 ðqt Þ ¼

where H is a dimensionless parameter H ¼ (xqt)4. Additionally, Eq. (11) indicates that the ratio sn/s0 depends only on H. From the statistics of the bilayer positions, the trans-bilayer structure factors and their ratio sn/s0 can be calculated. By comparing these results with the predictions of Eqs. (10) and (11), the in-plane correlation length x (that is, the phenomenological parameter Kc/B) can be derived. Furthermore, two regimes are expected with a crossover at qc w x1 in the curves of s0 and sn/s0 vs qt [69]. If qt is much larger than qc, fluctuations of different bilayers are incoherent, the fluctuation spectra s0 is proportional to q4 t and the ratio sn/s0 of the cross correlation between different bilayers decay exponentially like Hn with the averaged inter-bilayer spacing nd. Whereas, in the large-wavelength (qt << qc) regime, fluctuations of different bilayers are coherent, s0 is proportional to q2 t , sn/s0 tends toward unity in the infinite wavelength limit (i.e. H ¼ 0). In real space, the heighteheight correlation function is defined as [69,74,76]

(7)

where q is the wave vector, qz and qt are the z component and the projection into the (x,y)-plane of q, respectively. In simulations with finite box dimensions L, the components of wave vector take only discrete values qa ¼ ka(2p)/La, where a is the Cartesian indices, kt¼ (k2x þ k2y )1/2, kz is given by the number of bilayers N. According to the equipartition theorem, one obtains the average amplitudes of fluctuations [77]

D E jhðqt ; qz Þj2 ¼

2149

E

¼

2h1 q21

ZN 0

ds

1  J0

 pffiffiffiffiffihpffiffiffiffiffiffiffiffiffiffiffiffiffiffi i2n r 2s 1 þ s2  s x pffiffiffiffiffiffiffiffiffiffiffiffiffiffi s 1 þ s2

(13)

Here r ¼ (x2 þ y2)1/2, J0 is the zeroth-order Bessel function, q1 is the position of the first diffraction peak (q1 ¼ 2p/d), and h1 is the Caillé parameter,

h1 ¼

4pKB T pffiffiffiffiffiffiffiffi 8d2 BKc

(14)

Note that all heighteheight correlations are proportional to h1 and h1 has often been used to characterize the width of diffraction peaks in studying X-ray scattering spectra [73]. Similar fits in the real space allow us to extract the Caillé parameter h1. By combining this quantity with the resulting phenomenological parameter Kc/B in the Fourier space, the bending modulus Kc and compressibility modulus B can be determined. We should point out that theories of this type can well describe the fluctuations of the lamellar phase in binary amphiphilee solvent mixtures, wherein the amphiphiles are arranged in stacks of bilayers interspersed by solvent layers and possess a relatively high degree of order [69,72,78]. In symmetric ternary blends, the symmetrically swollen lamellar phase is formed, which consists of roughly equidistant copolymer rich bilayers, separated from each other by homopolymer-rich domains [79]. Similarly, this symmetrically swollen bilayer stack exhibits the liquid-like behavior along the lamellar plane and the position ordering in the direction perpendicular to the lamellar plane. In particular, for the system where the LAM phase is stable and not too close to the transition from LAM to BmE, the thermal fluctuations are restricted due to the presence of other bilayers and the LAM phase is only swollen to a certain degree. Thus, it is conceivable that the above DH elastic theories should work for the symmetrically swollen lamellar stacks of bilayers formed in symmetric ternary blends. Also, the method applied here to obtain the values of the elastic constants Kc and B is expected to work for measurements in this symmetrically swollen lamellar system as well.

2150

Z. Bai, H. Guo / Polymer 54 (2013) 2146e2157

3. Simulation method and model Dissipative particle dynamics (DPD) is a particle-based mesoscopic simulation technique which was first introduced by Hoogerbrugge and Koelman [42] in 1992 and improved by Espanol and Warren [80] in 1995. Each DPD bead represents a group of atoms or molecules. For easy numerical handling, the cutoff distance, bead mass, and temperature were set as the reduced units, that is, mffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ¼ kBT ¼ 1. Then the reduced time unit was defined as rc ¼ p s ¼ mrc =kB T . The dynamics of the DPD bead is governed by Newton’s equations of motion.

vri vv ¼ v i ; mi i ¼ f i vt vt

(15)

where ri, mi, and vi are the position, mass, and velocity of the ith bead, respectively. The force fi acting on the ith particle is a sum of three pairwise forces: a conservative force FC, a dissipative force FD, and a random force FR.

fi ¼

X R ðFCij þ FD ij þ Fij Þ

(16)

isj

where this sum runs over all neighbor particles within a cutoff radius rc. Note that the pairwise characteristic of all these forces ensures momentum conservation so that hydrodynamics is preserved. The conservative force is a soft repulsion [41],

( FCij ¼

  r aij 1  rij =rc b   ij 0 rij  rc

  rij < rc

(17)

where aijis the repulsion parameter which presents the maximum repulsion strength between particle i and j. rij ¼ rirj, rij ¼ jrij j, and b r ij ¼ rij =jrij j. The other two forces act together as a thermostat to the systems [41]

   D b b FD ij ¼ 4w rij r ij $vij r ij

(18)

  r ij FRij ¼ swR rij xij Dt 1=2 b

(19)

Here s and 4 are noise amplitude and dissipative parameter, respectively. vij ¼ vivj. wD(rij) and wR(rij) are distance dependent weight functions for dissipative force and random force, respectively. xij is a random number with zero mean and unit variance. In order to generate a correct equilibrium GibbseBoltzmann distribution, the dissipative and random forces have to satisfy the following relations [41,80]

h  i2   wD rij ¼ wR rij ;

s2 ¼ 24KB T

(20)

where KB is the Boltzmann constant and T is temperature. Following Groot and Warren [41], the two weight functions are chosen to be in the same form as the conservative force for simplicity,

h  i2   w rij ¼ wR rij ¼ D

(

2   rij < rc 1  rij =rc   0 rij  rc

(21)

In the present work, our main objective is to obtain a fundamental understanding of the interfacial properties and their relevance to the phase transitions in symmetric ternary blends. We consider a simple symmetric ternary system, which is composed of one type of short symmetric block copolymers and equal amounts of two short immiscible homopolymers of same

lengths, A2/B2/A4B4. Pioneering experiments by Bates et al. [11,15,31,81] have shown that the general features of phase diagrams experimentally observed for high molecular weight (MW) ternary blends are essentially preserved for systems of lower MWs from a few hundred to a few thousand. This thermodynamic correspondence between high- and low-MW systems indicates that low-MW systems, which equilibrate quickly, could serve as valuable model systems for studying the thermodynamic and dynamic processes of ternary polymer blends [11,15,31,79,81]. In addition, a theoretical work on symmetric binary polymer blends reported that the short homopolymer models such as A2 or B2 are polymeric-like with respect to thermodynamics, even if the Gaussian chain character is not preserved in them [82]. Moreover, on the basis of extensive lattice Monte Carlo simulations, Larson [83] found that the phase diagram obtained from short A4B4 diblock copolymers in monomeric A1 and B1 mixtures agrees well with that predicted by SCF theory for higher MW polymeric A/B/AB melts. Thus, despite the chains we used are not very long, our work on low-MW ternary blends has value in itself, since in light of the universal feature of the rich phase behavior in ternary polymer blends, many qualitative features of ternary polymer blends can be gained from the study of the low-MW ternary systems wherein the kinetic limitations in the high-MW blends are avoided and the system can readily equilibrate. Note that in addition to the conservative, dissipative and random forces, a harmonic spring force of FSij ¼ Crij is introduced to link the adjacent particles on a polymer backbone together, where, C is the spring constant and C ¼ 4.0 in this study. In the simulation, the system density was fixed at 3.0, thus the repulsion parameters between identical types of beads were set at the usual value of aAA ¼ aBB ¼ 25 according to the reference [41]. To introduce a bias that induces phase separation and to produce a representative phase behavior at intermediate incompatibility, we restricted the repulsive interaction parameter for the unlike particle pairs to a fixed value of aAB ¼ 50 and focused on the A4B4/A2/B2 system with various homopolymer volume fractions FH. It is worthwhile to note that with this set of parameter space, we can study the interfacial behavior of LAM and 2P phases as a function of FH and explore its relevance to the phase transitions from 2P and LAM to BmE. In this study, we used the velocity-Verlet algorithm to integrate the equations of motion with a time step of Dt ¼ 0.04s, and the noise amplitude was fixed at s ¼ 3. All the simulations were performed in a NVT ensemble with periodical boundary conditions applied in all three dimensions. To identify the phases and phase boundaries of symmetric A4B4/A2/B2 ternary blends as a function of FH at fixed interaction parameters, a series of DPD simulations were carried out in a box with a size of Lx  Ly  Lz(h30rc). At the beginning, these simulations were started by equilibrating the systems in the athermal limit, i.e. by setting aij ¼ 25 for all DPD beads. Then, those totally disordered configurations were quenched to the desired aij value and equilibrated again. The equilibration is considered to be reached when those quantities such as the total energy, the interfacial tension, the global properties of the chains and molecular aggregates, as well as the density profiles do not change with the simulation time. Additionally, to reduce the statistical errors all the calculated results in this paper were obtained by averaging over 100 configurations after the simulated systems were well equilibrated. To check the influences of box size, we performed additional simulations in a larger simulation box of Lx ¼ Ly ¼ Lz ¼ 40rc. For such a system size we do not find finite size corrections to the whole shape of phase diagram. Thus, even though we have not performed the same simulation in a much larger system, we are sure that our results are not an artifact of finite-size effect, but of general relevance for real ternary blend systems. Measurement of

Z. Bai, H. Guo / Polymer 54 (2013) 2146e2157

the interfacial tension and elastic constants from the undulation spectrum requires simulations on a relatively large system, as this method is only applicable in the ‘‘long-wavelength’’ limit. By examining the effect of box size, we found that for the two employed system sizes of Lx ¼ Ly ¼ Lz ¼ 40rc and Lx ¼ Ly ¼ Lz ¼ 30rc the structure factors derived from monolayer fluctuations nearly fall onto a single curve and can be fitted in the long wavelength limit to the theoretical predictions of Eq. (4). Also, both systems yield matching estimates for elastic modulus and the calculated interfacial tension does not show pronounced system dependence (the detailed description can be found in Supporting Information). Recent studies indicated that the bending modulus is essentially independent of the system size if employed system size is much larger than the interfacial thickness [68,84]. All these suggest that Lx ¼ Ly ¼ Lz ¼ 30rc is a reasonable box size and therefore all the monolayer fluctuation analysis for the 2P phase were performed in systems of this size. In addition, we found that for the LAM phase not too close to the transition point from LAM to BmE when simulation boxes containing as many as eight lamellae the simulation results can be well compared to the above mentioned continuum theory for stacked membranes. Note that in practice the lamellar phase at equilibrium is tensionless and thus producing equilibrium lamellar structures with the bilayers parallel to the one face of the simulation cell is crucial for further fluctuation analysis. Following an efficient procedure for generating equilibrium lamellae in diblock copolymer systems [85], the equilibrium lamellar phase for ternary system A4B4/A2/B2 at a given repulsion parameter can be easily generated in the following three steps. Firstly a starting lamellar configuration with the equilibrium layer spacing, which is determined by measuring the variation of the diagonal pressure tensor elements and the order parameter, is “built by hand”. Then it is equilibrated at a large repulsion parameter of aij ¼ 75, equivalent to a low temperature, deep in the lamellar phase. Afterward, the repulsion parameter is decreased to the desired value of aij ¼ 50 equivalent to heating to the desired temperature and equilibrated for a sufficient time. Finally the equilibrium flat lamellar systems in ternary blends are generated. During the simulation, we monitor the structure factor to trace the morphology of the resulting phases. The structure factor is defined as [86]:

SðqÞ ¼

( *" X X 1 q

L3

    exp iq$rj f 4j

rj

#2 +), X

1

(22)

q

where L is the length of simulation box. f(4j) is the Fourier component of the concentration fluctuation, f ð4j ðrÞÞ ¼ 4jB  4jA  h4jB  4jA i, and hi denotes thermal statistical average. As mentioned in the Introduction, a thermodynamically stable BmE phase could be formed in symmetric ternary polymer blends [11,15,81,87]. Teubner and Strey [88] have derived an expression of S(q) for the BmE morphology on the basis of phenomenal Landau theory, which is so-called T-S model:

  1 S q w a2 þ c1 q2 þ c2 q4

(23)

where all coefficients are system dependent and satisfy a2 > 0, c1 < 0, c2 > 0, 4a2c2c21 > 0, 1 < c1/(4a2c2)1/2 < 0. The above model well describes the observed broad peak and the distinct q4 decaying at the high q range in the SANS (Small Angle Neutron Scattering) scattering spectra of BmE structure [89]. Consequently, we use Eq. (23) to identify whether BmE-like structure is formed in symmetric ternary blends. In addition, the resulting coefficients are often employed to calculate the amphiphilicity factor

2151

fa ¼ c1/(4a2c2)1/2, which measures the stability of the microemulsion and has a magnitude of 1 w 0 [81]. 4. Results and discussion In this section, we firstly characterize the different phases and estimate the phase boundaries of ternary blend A4B4/A2/B2 as a function of total volume fraction of homopolymers FH at fixed temperature. Note that in the vicinity of the phase transitions from 2P to BmE and from BmE to LAM the incremented step of DFH is set to 2.5% in order to improve the precision in the location of phase transition points Then, we analyze the thermal fluctuations in the lamellar phase and in the 2P phase, and discuss the behavior of the interfacial tension g and the elastic constants upon approaching the phase transitions from LAM to BmE and from 2P to BmE. 4.1. Identification of phases and phase regions Various phases are identified by snapshots, density profiles and structure factors, as typically illustrated in Fig. 1. In the lamellar (LAM) phase, copolymers form periodic layers with a long-range order and homopolymers are distributed in the corresponding lamellar layers, as typically shown in Fig. 1(a). In the phaseseparated (2P) state, the copolymers form two discrete monolayers at the interfaces between homopolymer phases, due to the periodic boundary conditions [90], as shown in Fig. 1(b). For the BmE phase, homopolymers A and B form co-continuous interwoven networks and are divided by undulating copolymer monolayers. The snapshot in Fig. 1(c) indicates that these monolayers are disorderedly oriented without any long-range order. When the structure factors in the BmE phase are fitted by the T-S equation, the fitting coefficients fully satisfy the criteria of “a > 0, b < 0, c > 0 and 4acb2 > 0”, as indicated in Fig. 1(c). Furthermore, the amphiphilicity factor fa we measured is approximately equal to 0.89, which also satisfies the criteria of “1 < fa < 0”, indicating that the BmE structure formed in our simulation is stable. By carefully characterizing the phases formed in A4B4/A2/B2 blends with different FH at fixed aAB ¼ 50 and aAA ¼ aBB ¼ 25, we can determine the phase regions, as listed in Table 1. In accordance with the experimental [11,14,15,81,87] and theoretical results [30e33], a transition sequence of LAM to BmE to 2P phase with the increase of VH is seen. 4.2. Elastic properties of LAM phase 4.2.1. Fluctuation analysis of LAM system with VH ¼ 40% Here, we choose an LAM system with VH ¼ 40% as an example to explain the whole process of fluctuation analysis for the LAM phase. Based on the variation of diagonal pressure tensor elements and order parameter with bilayer spacing under the quasistatic uniaxial compression and dilation strain Refs. [64], the equilibrium bilayer spacing for the tensionless flat lamellar phase at VH ¼ 40% was determined as d ¼ 8.32rc. For comparisons with the behavior described by DH theory [69,74] for an infinitely thick stack of bilayer, two typical equilibrium flat lamellar configurations with Lx ¼ Ly ¼ 30rc and different lengths of Lz ¼ Nd and N ¼ 8 and 4 along the z direction were generated [85]. For the fluctuation analysis, we parameterized each bilayer n by a set of the local height Zn(x,y) and thickness Tn(x,y) functions. In practice, we just considered discrete values of x and y, i.e. x ¼ nxLx/Nx, y ¼ nyLy/Ny, and here we assigned Nx ¼ Ny ¼ 30. For each position (x,y), the height Zn(x,y) and thickness Tn(x,y) of the nth bilayer were defined as the mean and the difference of two height z(x,y) values of the upper and lower copolymer monolayers. Here, this z value at each point (x,y) for each copolymer monolayer was given by the average of the z-coordinates of the joint

2152

Z. Bai, H. Guo / Polymer 54 (2013) 2146e2157

Fig. 1. Proper physical characterization and corresponding snapshots for typical phases selected from A2/B2/A4B4 blends in our simulations. Density profiles of (a) an LAM system with FH ¼ 40% and (b) a 2P system with FH ¼ 92.5%. The blocks A, blocks B, homopolymers A and B, and also their corresponding density profiles are shown in green, red, pink and cyan colors, respectively. The abscissa represents the normalized length along the layer normal. (c) The structure factor of the BmE phase formed from the system with FH ¼ 77.5% and box length of 40rc. The red smooth line corresponds to the fitting curve of T-S model and fitting parameters are a2 ¼ 0.00412 > 0, c1 ¼ 0.08232 < 0, c2 ¼ 0.51243 > 0, 4a2c2c21 ¼ 1.67Ee3 > 0 (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.).

points of copolymer molecules that have ((xkx)2þ(yky)2)1/2 < R. It should keep in mind that the size of R is varied for different systems in order to ensure that at least one diblock copolymer molecule is assigned to each grid, and we typically chose R ¼ 1e1.5rc in our calculations. Then the local deviation of the height of the nth bilayer from its average value was presented ashn ¼ Zn ðx; yÞ  Z n , where the mean position Z n was determined separately in each csonfiguP ration (Z n ¼ x;y Zn ðx; yÞ=ðNx Ny Þ). First, we analyze the bilayer fluctuations in real space. Fig. 2(a) shows the distribution of inter-bilayer distances P n hjZn ðx; yÞ  Z0 ðx; yÞji for the lamellar stack of 8 bilayers. The peaks are periodic, implying the smectic ordering of the bilayers along the layer normal. Moreover, the nth peak corresponds to the distance between bilayers, which are separated from each other by n layer(s) of homopolymer A-rich and B-rich domains. From the width of the first peak we determined the value of the Caillé parameter, h1 ¼ 0.47. By fitting each peak to a Gaussian function, the mean and the variance are derived, as shown in Fig. 2(b) and (c). Obviously, the mean inter-bilayer distance is proportional to n, hjZn ðx; yÞ  Z0 ðx; yÞji ¼ nd with d ¼ 8.32rc. The variances correspond to the heighteheight fluctuations and are hence compared with the prediction of Eq. (13) using h1 ¼ 0.47 in Fig. 2(c). Similar to the results on the lamellar phase in an amphiphileesolvent mixture, the DH theory describes our simulation data relatively well for small n (n < 4) and at large n some deviations between theory and simulation data are observed, which are presumably stemmed from the finite-size effects [69]. Next, the position fluctuations and correlations of bilayers are investigated in Fourier space. Fig. 3(a) shows the trans-bilayer structure factor sn(qt) for n ¼ 0, 1, 2, as a function of qt (qt¼(q2x þ q2y )1/2). It should be noted that the continuum theory is no longer valid below molecular length scales, thus we should carefully choose the qt scope for analyzing our simulation data. Generally, the qtmin is selected to be equal to 2p/L ¼ 0.20944rc1. While for the choice of qtmax, we should keep in mind that 2p/qtmax must be larger than the molecular length in any cases. In our system, the molecular length scale is about 3.2rc (the root of mean-squared end-to-end distance of A4B4 molecule), thus the Table 1 The position of phases formed in A4B4/A2/B2 blends as a function of FH at fixed aAB ¼ 50 and aAA ¼ aBB ¼ 25.

FH

0 w 75%

77.5 w 82.5%

85 w 100%

Phase

LAM

BmE

2P

value of qtmax should smaller than 2rc1. Fig. 3(b) shows the ratios of s1,2/s0 as a function of qt, which are fitted by Eq. (11) with only one fit parameter x. Remarkably, the theoretical predictions fit our data on an LAM stack of 8 bilayers very well. Additionally, two fits give the same values of in-plane correlation length within the errors x z 2.05  0.06rc. Finally, from the two fitting parameters of Caillé parameter h1 and in-plane correlation lengths x, the bending modulus Kc and compressibility modulus B for a symmetrically swollen lamellar phase in A4B4/A2/B2 blends with VH ¼ 40% are determined to be 0.842KBT and 0.044 KBT$rc4, respectively. To illustrate possible finite size effects, we also performed the position fluctuation analysis for an LAM stack of 4 bilayers in real space and in Fourier space, respectively. The corresponding results are presented in Fig. 2(b), (c) and Fig. 3(c), respectively. As illustrated in Fig. 2(b) and (c), the bilayer spacing d does not depend significantly on the system size, but the variances are reduced in the small systems, which is presumably due to the finite-sizeeffects [69]. As a result, the Caillé parameter is underestimated in the small systems (4 bilayers). Furthermore, a careful analysis of the ratios of s1,2/s0 in Fig. 3(c) indicated a relatively large discrepancy between simulation data on small systems and theoretical predictions. As the bilayer fluctuations are correlated more strongly in small systems, the finite thickness of the simulated system along the bilayer normal could seriously affect the thermal fluctuations. Although Eq. (11) is deduced for an infinitely thick stack of bilayers, in our model an LAM stack of 8 bilayers is sufficiently large to reproduce the behavior described by continuum DH theory. From this, we generated a set of equilibrium flat LAM stack of 8 bilayers with various VH to study the elastic properties in the lamellar phase. 4.2.2. Elastic properties of LAM phase as a function of VH After exemplifying the analysis process for the bilayer fluctuations in a lamellar phase, we can compute the elastic constants for the LAM systems with different VH in the same way. We want to emphasize that even adding VH ¼ 75% homopolymers the lamellar phase formed in our ternary system can still be conceived as a stack of elastic sheets with a well-defined inter-bilayer spacing (the detailed description can be found in Supporting Information). Similar to the before-mentioned sample with VH ¼ 40%, it displays long wave-length thermal fluctuations and can be described successfully by the DH theory. Before presenting our results for the elastic coefficients, as a final consistency test, it is worthwhile to discuss the shape of autocorrelation fluctuation spectra. As mentioned in the theoretical

Z. Bai, H. Guo / Polymer 54 (2013) 2146e2157

2153

Fig. 2. (a) Distribution of inter-bilayer distances in an LAM stack of 8 bilayers at VH ¼ 40%. (b) The mean inter-bilayer distance vs n for the above lamellar stack and its analog configuration with N ¼ 4 at VH ¼ 40%. The dash line is a linear fit. (c) Variance of the distribution of inter-bilayer distances vs n for the above two LAM stacks. The dash line is the prediction of Eq. (13), using Caillé parameter h1 ¼ 0.47.

background section, the DH theory predicts two regimes with a crossover at qc w x1 in the curve of s0. If qt is much larger than qc, fluctuations of different bilayers are incoherent, the fluctuation spectra s0 is proportional to q4 t . Whereas, in the large-wavelength (qt<
proportional to q2.050.05 , which seems to follow the asymptotic t behavior predicted by DH theory for the large-wavelength q2.0 t (qt<x1, s0 of this system does not exhibit the scaling behavior of s0 w q4 t expected from the DH theory for the small-wavelength (qt>>qc) regime. As the in-plane correlation length x of this system is relatively small, the continuum approximation of the DH theory is no longer valid below this correlation length and thus the free membrane regime characterized by a scaling law of s0 w q4 t cannot be observed [69]. While for the LAM system with VH ¼ 55%, the inplane correlation length increases to x ¼ 2.75  0.03rc. We observed that the DH theory describes the data well and the curve of at qt x . Apparently, the layers in the scaling law s0 w q4.030.08 t t LAM system with VH ¼ 55% behave more like free, unconstrained membranes at small wavelength regime but constrained and

Fig. 3. (a) Transbilayers structure factor for an LAM stack of 8 bilayers at VH ¼ 40%. The ratios of s1,2/s0 as a function of qt for (b) the above lamellar stack and (c) its analog configuration with N ¼ 4 at VH ¼ 40%. The dots represent simulation data. The solid lines in (b) and (c) correspond to the predictions of Eq. (11) with the fit parameter x.

2154

Z. Bai, H. Guo / Polymer 54 (2013) 2146e2157

Fig. 5. (a) In-plane correlation length x, (b) Compressibility modulus B, and (c) bending modulus Kc for the LAM stack of 8 bilayers as a function of VH. Fig. 4. The autocorrelation fluctuation spectrum s0 (qt) for the LAM stack of 8 bilayers at (a) VH ¼ 0% and Lx ¼ Ly ¼ 30rc, (b) VH ¼ 55% and Lx ¼ Ly ¼ 50rc, and (c) VH ¼ 75% and Lx ¼ Ly ¼ 50rc.

coherent membranes at the large wavelength regime. As more homopolymers are added to the blend, the in-plane correlation length is further increased, for instance, x ¼ 5.99  0.05rc for the LAM system with VH ¼ 75%. For such a system, the incoherent membrane fluctuations become more pronounced. For example, s0 just exhibits behavior at qt>x1, which approaches the expected a q3.950.06 t 4 scaling law qt for free membranes. We notice that because of the finite box length along (x,y)-plane used in simulations (despite that we have already enlarged the size of simulation box to Lx ¼ Ly ¼ 50rc), we can not get enough points in the regime qt
simulation studies [61,65,68,69,71]. As an example, for our lamellar stack at VH ¼ 20%, the bending modulus and compressibility modulus are 2.01KBT and 0.26KBT$rc4, respectively, very close to the values of Kc ¼ 4.0KBT and B ¼ 0.13KBT$rc4 reported by Loison et al. [69] for the lamellar stack in a binary amphiphile-solvent mixture at a solvent volume fraction of 20%. A further careful analysis of Kc values in Fig. 5(c) and the phase region data in Table 1 indicates that in the ultimate vicinity of the phase transition (e.g. VH ¼ 77.5%), the bending modulus should turn to be smaller than KBT. As a consequence thermal fluctuations could become so large that the highly swollen LAM phase can not sustain the orientation order of the lamellae, and a BmE phase is produced instead. Indeed, our results justify a direct correlation between the vanishing interfacial tension and ultralow bending rigidity of a swollen LAM phase and the proximity to the microemulsion state. Such a correspondence was anticipated by Bates et al. [11]. In particular, our results are also consistent with the picture that the thermal fluctuations bring about the BmE phase [28,32,33]. Additionally, the phase transition from the LAM to the BmE phase can be understood by studying the persistence length of the interface (l) and the lamellar spacing (d). The so-called persistence length is defined as [91]

lzaexpð2pKc =KB TÞ

(24)

where a is a microscopic length scale and in our simulations we chose a of the order of the cut-off distance, i.e. a ¼ 1. As illustrated by Fig. 6, consistent with experimental and theoretical results [11,40,91], the lamellar periodicity increases with VH, but the persistence length decreases with VH. At VH ¼ 75%, l is slightly larger than d. While, if more homopolymers are added to the blend, e.g. VH ¼ 77.5%, l would reduce further. The correlation of the lamellar stack along the normal direction becomes so weak that fluctuations of the local interface position are expected to destroy the long-rang lamellar order, thus the bicontinuous microemulsion structure is formed. In good conformity with a recent fieldtheoretic Monte Carlo study [32], the occurrence of such an LAM to BmE phase transition is also related to matching length scales of the persistence length of composition fluctuations along the microdomain boundaries and of the lamellar spacing. Apparently, the phase transition from the LAM to the BmE phase is expected to occur when the bending modulus turns out to be smaller than KBT or when the persistent length is no larger than the lamellar spacing. Also the two criteria yield the same predictions for the transition point.

Z. Bai, H. Guo / Polymer 54 (2013) 2146e2157

Fig. 6. The persistence length l and the lamellar spacing d for the LAM stack of 8 bilayers as a function of VH.

4.3. Interfacial properties of 2P phase To perform the fluctuation analysis on the interfaces in the 2P phase, we follow the elegant procedure employed by Smit et al. [35,60] in studying the interfacial properties of surfactant monolayers between oil and water. In our DPD simulation, height fluctuation function h(x, y) was defined as the displacement of the local position z(x, y) of each of the two interfaces from its mean position, where z(x, y) was taken to be the weight-average locations of the A (or B) type beads which are the nearest to the interface. The resulting interfacial structure factors S(q) are shown on a double logarithmic plot in Fig. 7(a), for the binary blend (VH ¼ 100%) and the ternary blends with VH ¼ 90% and 85%, respectively. It is worthy emphasizing that even for the 2P system with a homopolymer content as low as FH ¼ 85%, the interfaces display small thermal height fluctuations around their average flat positions (the detailed description can be found in Supporting Information) and still can be regarded as single smoothly undulating surfaces with small curvatures. When plotting 1/(q2S(q)) vs q2 in the inset, we get a straight line at small values of the wave vector q, suggesting that Eq. (4) is indeed a good description of the fluctuation spectrum of the monolayer at long wavelengths. Then by fitting the first four points to y ¼ gþkx, as illustrated in the inserted plot of Fig. 7(a), we get the bending modulus k from the slope and the interfacial tension g on the y-axis by extrapolating the linear part to q ¼ 0. The inferred g and k as a function of VH for 2P phase are shown in Fig. 7(b). As expected, the interfacial tension of the binary blend (VH ¼ 100%) which just contains immiscible homopolymers A and

2155

B is relatively high. Adding diblock copolymers into the binary blends leads to the reduction of interfacial tension, as they are preferentially absorbed at the interfaces and further screen the repulsion between incompatible homopolymers. When the homopolymer content is lowered to VH ¼ 85%, the interfacial tension becomes relatively small, 0.45  0.04KBT/rc2. As more diblock copolymers are added to the blend, i.e. at VH ¼ 82.5%, we find that the interfacial tension becomes vanishingly small. There is no longer a free-energy penalty associated with the formation of extra interfaces. In this case, as shown in Table 1, the thermodynamically stable bicontinuous microemulsions (BmEs) are formed. These results for the interfacial tension are in good agreement with the predictions by experiments [11,15,81] and theories [25,26,28,38]. In particular, our results also support the notion that the phase transition from 2P to BmE coincides with the interfacial tension attaining an ultralow value. Under this circumstance the bending modulus becomes important. As illustrated in Fig. 7(b), in the 2P phase k increases almost monotonically with the addition of diblock copolymers, which confirms that the bending contribution to the interfacial energy becomes more important due to the packing constraints of the copolymers at the interface [25]. Actually, similar behavior has been observed by Smit et al. [35,60]. In contrast to the above results, the opposite trend was reported by Laradji and Mouritsen [38] for a surfactant monolayer on the liquid/ liquid interface. They found that k firstly decreases with increasing density of surfactants for low densities but increases with a further increase of the surfactant density at the interface and attributed this nonmonotonic dependence to a coupling between fluctuations of the interface position and the surfactant orientation. Besides, Smit et al. [35,60] pointed out that for strong or long surfactants local fluctuations in the orientation at the interface, which could decrease the bending modulus from the bare oil/water value, disappears and thus bending the interface becomes increasingly costly with increasing surfactant density, which, of course, can explain the qualitative difference between effect of the copolymer (or surfactant) interface coverage on the bending rigidity. Especially, at VH ¼ 85%, k is increased to 0.85  0.09KBT. A further addition of diblock copolymers to the blend, i.e. at VH ¼ 82.5%, kwould increase to be comparable to KBT. Note that at ultralow or vanishing interfacial tension, a low bending modulus means less resistance of the interface against bending or more proneness to the thermal fluctuations. The fluctuations give rise an entropic repulsion between neighboring monolayers, thus the monolayers at VH ¼ 82.5% become unstable and the blend is expected to undergo a transition from 2P into BmE phases, in which homopolymers form stable interweaving continuous microscopic domains as typically illustrated in Fig. 1(c). Therefore, the transition between 2P and BmE

Fig. 7. (a) lnS(q) vs ln(q) for the binary blend (VH ¼ 100%) and ternary blends with VH ¼ 90% and 85%. The inset shows corresponding 1/(q2S(q)) vs q2. The full lines are fits of the four lowest q values to y ¼ gþkx. The constants k and g are determined by fitting the data. The lowest q is limited by box dimensions, qmin ¼ 2p/Lx ¼ 0.2093rc1. (b) The dependence of interfacial tension g and the bending modulus k of the copolymer monolayer as a function of VH.

2156

Z. Bai, H. Guo / Polymer 54 (2013) 2146e2157

structure not only corresponds to the saturation of the A2/B2 interface with A4B4 diblock copolymers (i.e. attaining a vanishing interfacial tension), which allows a BmE structure to accommodate an extensive amount of copolymer-laden internal interfaces [11,25,26,91], but also is associated with the maintenance of the copolymer monolayers at a low bending modulus k w KBT, which permits thermal fluctuations to stabilize BmE structure [11,25,26,38,40]. 5. Conclusions In the present work, our main objective is to obtain a fundamental understanding of the interfacial properties and their relevance to the phase transitions in ternary symmetric blends. For this, we consider a simple ternary symmetric system, which is composed of one type of short symmetric block copolymers and equal amounts of two short immiscible homopolymers of same lengths, A4B4/A2/B2. We found that in accordance with the experimental [11,14,15,81,87] and theoretical results [30e33] the blend undergoes transitions from the LAM to BmE and then to 2P phase sequentially with increasing the volume fraction of homopolymers VH. For the 2P systems not too close to the phase transition to the BmE, we found that Helfrich model is indeed a good description of the fluctuation spectrum of the copolymer monolayer at long wavelengths. Therefore, the interfacial tension g and the bending modulus k can be straightforward extracted from the undulation spectrum. In good agreement with the predictions by experiments [11,15,81] and theories [25,26,28,38], with the addition of diblock copolymers, g reduces but k increases. In particular, our results support the picture that the phase transition from 2P to BmE not only coincides with attaining a vanishing interfacial tension, but is also associated with keeping the interfacial copolymer monolayers at a low bending modulus of k w KBT. Although the evaluation of g and k via the fluctuation spectrum has now become commonplace for single layer systems, it is exceedingly difficult to measure k for the multilayers (such as those appearing a lamellar stack) in much the same way. As a result, up to now little work has been done to investigate the undulations and thickness fluctuations in the lamellar phase of ternary blends and little attention has been paid to the interfacial properties of these lamellar stacks. Conforming to a recent work on the lamellar phase in binary amphiphileesolvent mixtures [69], we found that the DH theory can describe the fluctuations of the lamellar phase (not too close to the phase transition to the BmE) in ternary blends in a satisfactory way. In particular, we found that in our model an LAM stack of 8 bilayers is sufficiently large to re-produce the behavior (including the distribution of inter-bilayer distances, the transbilayer structure factors) described by continuum DH theory. Moreover, by analyzing the shape of auto-correlation fluctuation spectra s0, we see that as the LAM stack is swollen increasingly with the addition of homopolymers, the bilayers in LAM phases experience a crossover from the single coherent fluctuation with a scaling law s0 w q2 t to the coexisting of coherent fluctuation with s0 w q2 t and incoherent fluctuation with a scaling behavior of s0 w q4 t wherein the incoherent free membrane regime starts at larger length scales and becomes more pronounced. Indeed, all the analysis of the position fluctuations and correlations of bilayers are rather important, as it demonstrates the feasibility of DH theory to calculate the bending modulus Kc and compressibility modulus B of the lamellar stack from the simulation data. Furthermore, with the increase of VH, the in-plane correlation length x increases, but B and Kc decrease, suggesting that adding homopolymers either shorter than or comparable to the diblock copolymer chains can swell the LAM phase and thus fluctuations between adjacent

bilayers become less coherent, bilayers become more flexible, and the swollen LAM system becomes more susceptible to thermal fluctuations. Most notably, our results indicate that the phase transition from LAM to the BmE occurs when the bending modulus turns out to be smaller than KBT or when the persistent length is no larger than the lamellar spacing. Also the two criteria yield the same predictions for the transition point. In a word, our results shed more light on the topic of the correlation of the peculiar changes in interfacial properties to the phase transitions from 2P or LAM to the BmE state, which not only justifies the conjecture of the phase transition from 2P to BmE corresponding to the saturation of interfaces and the transition from LAM to BmE coinciding with the achievement of the flexibility of a copolymer-laden interface above a threshold value (i.e. k  KBT) [11,25,26,38,40,91], but also supports the notion that the thermal fluctuations bring about the BmE phase [11,26,38e40]. Of curse, our present study is restricted to studying on low-MW ternary system. As the thermal fluctuation is strongly affected by the molecular weight, for clarifying the role of chain length systematic studies of several comparison systems with varying chain lengths are underway, and the chain length dependence will be addressed in an upcoming paper. Acknowledgment We thank Dr. Claire Loison for helpful discussions on this research. This work was supported by grants from the NSF of China (21174154, 20874110, and 50930002). We thank the supercomputing center of CAS for supporting computing resources. Appendix A. Supplementary data Supplementary data related to this article can be found at http:// dx.doi.org/10.1016/j.polymer.2013.02.011. References [1] Binder K. Phase-transitions in polymer blends and block-copolymer melts e some recent developments. Berlin: Springer-Verlag; 1994. [2] Gao Y, Ren K, Ning NY, Fu Q, Wang K, Zhang Q. Polymer 2012;53:2792e801. [3] Androsch R, Di Lorenzo ML, Schick C, Wunderlich B. Polymer 2010;51: 4639e62. [4] Anastasiadis SH. Interfacial tension in binary polymer blends and the effects of copolymers as emulsifying agents. Berlin: Springer-Verlag; 2011. [5] Hong KM, Noolandi J. Macromolecules 1981;14:727e36. [6] Hong KM, Noolandi J. Macromolecules 1981;14:736e42. [7] Noolandi J, Hong KM. Macromolecules 1982;15:482e92. [8] Noolandi J, Hong KM. Macromolecules 1984;17:1531e7. [9] Retsos H, Anastasiadis SH, Pispas S, Mays JW, Hadjichristidis N. Macromolecules 2004;37:524e37. [10] Retsos H, Margiolaki I, Messaritaki A, Anastasiadis SH. Macromolecules 2001; 34:5295e305. [11] Bates FS, Maurer WW, Lipic PM, Hillmyer MA, Almdal K, Mortensen K, et al. Phys Rev Lett 1997;79:849e52. [12] Tanaka H, Hasegawa H, Hashimoto T. Macromolecules 1991;24:240e51. [13] Tanaka H, Hashimoto T. Macromolecules 1991;24:5713e20. [14] Stoykovich MP, Edwards EW, Solak HH, Nealey PF. Phys Rev Lett 2006;97: 147802. [15] Hillmyer MA, Maurer WW, Lodge TP, Bates FS, Almdal K. J Phys Chem B 1999; 103:4814e24. [16] Fleury G, Bates FS. Soft Matter 2010;6:2751e9. [17] Jones BH, Lodge TP. J Am Chem Soc 2009;131:1676e7. [18] Jones BH, Lodge TP. Chem Mater 2010;22:1279e81. [19] Zhou N, Bates FS, Lodge TP. Nano Lett 2006;6:2354e7. [20] Gan LM, Chow PY, Liu ZL, Han M, Quek CH. Chem Commun 2005:4459e61. [21] Gan LM, Liu J, Poon LP, Chew CH, Gan LH. Polymer 1997;38:5339e45. [22] Wang LS, Chow PY, Phan TT, Lim IJ, Yang YY. Adv Funct Mater 2006;16:1171e8. [23] Wang ZG, Safran SA. J Phys France 1990;51:185e200. [24] Wang ZG, Safran SA. J Chem Phys 1991;94:679e87. [25] Laradji M, Desai RC. J Chem Phys 1998;108:4662e74. [26] Matsen MW. J Chem Phys 1999;110:4658e67. [27] Chang K, Morse DC. Macromolecules 2006;39:7397e406. [28] Muller M, Schick M. J Chem Phys 1996;105:8885e901.

Z. Bai, H. Guo / Polymer 54 (2013) 2146e2157 [29] [30] [31] [32] [33] [34] [35] [36] [37] [38] [39] [40] [41] [42] [43] [44] [45] [46] [47] [48] [49] [50] [51] [52] [53] [54] [55] [56] [57] [58] [59] [60] [61]

Werner A, Schmid F, Muller M. J Chem Phys 1999;110:5370e9. Kielhorn L, Muthukumar M. J Chem Phys 1997;107:5588e608. Fredrickson GH, Bates FS. J Polym Sci Part B: Polym Phys 1997;35:2775e86. Duchs D, Ganesan V, Fredrickson GH, Schmid F. Macromolecules 2003;36: 9237e48. Duchs D, Schmid F. J Chem Phys 2004;121:2798e805. Senapati S, Berkowitz ML. Phys Rev Lett 2001;87:176101. Rekvig L, Hafskjold B, Smit B. Phys Rev Lett 2004;92:116101. Petrache HI, Gouliaev N, Tristram-Nagle S, Zhang RT, Suter RM, Nagle JF. Phys Rev E 1998;57:7014e24. Bouglet G, Ligoure C. Eur Phys J B 1999;9:137e47. Laradji M, Mouritsen OG. J Chem Phys 2000;112:8621e30. Thompson RB, Matsen MW. Phys Rev Lett 2000;85:670e3. Komura SJ. Phys Condens Matter 2007;19:463101. Groot RD, Warren PB. J Chem Phys 1997;107:4423e35. Hoogerbrugge PJ, Koelman J. Europhys Lett 1992;19:155e60. He LL, Zhang LX, Xia A, Liang HJ. J Chem Phys 2009;130:144907. He LL, Zhang LX, Ye YS, Liang HJ. J Phys Chem B 2010;114:7189e200. He PT, Li XJ, Deng MG, Chen T, Liang HJ. Soft Matter 2010;6:1539e46. He PT, Li XJ, Kou DZ, Deng MG, Liang HJ. J Chem Phys 2010;132:204905. Wu SG, Guo HX. Sci China Ser B-Chem 2008;51:743e50. Groot RD, Madden TJ. J Chem Phys 1998;108:8713e24. Li XJ, Guo JY, Liu Y, Liang HJ. J Chem Phys 2009;130:074908. Liu H, Qian HJ, Zhao Y, Lu ZY. J Chem Phys 2007;127:144903. Soto-Figueroa C, Rodriguez-Hidalgo MDR, Martinez-Magadan JM, Vicente L. Macromolecules 2008;41:3297e304. Zhao Y, You LY, Lu ZY, Sun CC. Polymer 2009;50:5333e40. Zhang ZM, Guo HX. J Chem Phys 2010;133:144911. Zhang ZM, Guo HX. Soft Matter 2012;8:5168e74. Bai ZQ, Xia YZ, Shi SX, Guo HX. Acta Polym Sin 2011:530e6. Gai JG, Li HL, Schrauwen C, Hu GH. Polymer 2009;50:336e46. Rodriguez-Hidalgo MD, Soto-Figueroa C, Martinez-Magadan JM, Vicente L. Polymer 2009;50:4596e601. Huang MX, Li ZQ, Guo HX. Soft Matter 2012;8:2834e45. Zhou Y, Long XP, Zeng QX. Polymer 2011;52:6110e6. Rekvig L, Hafskjold B, Smit B. J Chem Phys 2004;120:4897e905. Rodgers JM, Sorensen J, de Meyer FJM, Schiott B, Smit B. J Phys Chem B 2012; 116:1551e69.

2157

[62] Li XJ, Liu Y, Wang L, Deng MG, Liang HJ. Phys Chem Chem Phys 2009;11: 4051e9. [63] Wu SG, Guo HX. J Phys Chem B 2009;113:589e91. [64] Guo HX, Kremer K. J Chem Phys 2003;118:7714e23. [65] Goetz R, Gompper G, Lipowsky R. Phys Rev Lett 1999;82:221e4. [66] Helfrich W. Z Naturforsch C: Biosci 1973;28:693e703. [67] Shi WX, Guo HX. J Phys Chem B 2010;114:6365e76. [68] Marrink SJ, Mark AE. J Phys Chem B 2001;105:6122e7. [69] Loison C, Mareschal M, Kremer K, Schmid F. J Chem Phys 2003;119:13138e48. [70] Canham PB. J Theor Biol 1970;26:61e81. [71] Lindahl E, Edholm O. Biophys J 2000;79:426e33. [72] Schmid F. Macromol Rapid Commun 2009;30:741e51. [73] Caille A. C R Acad Sci Ser B 1972;274:891e3. [74] Lei N, Safinya CR, Bruinsma RF. J Phys 1995;5:1155e63. [75] Brannigan G, Lin LCL, Brown FLH. Eur Biophys J Biophys Lett 2006;35: 104e24. [76] Lyatskaya Y, Liu YF, Tristram-Nagle S, Katsaras J, Nagle JF. Phys Rev E 2001;63: 011907. [77] Safran SA. Statistical thermodynamics of surfaces, interfaces, and membranes. Reading, Massachusetts: Addison-Wesley; 1994. [78] Loison C, Mareschal M, Schmid F. Comput Phys Commun 2005;169:99e103. [79] Sun DC, Guo HX. Polymer 2011;52:5922e32. [80] Espanol P, Warren P. Europhys Lett 1995;30:191e6. [81] Morkved TL, Stepanek P, Krishnan K, Bates FS, Lodge TP. J Chem Phys 2001; 114:7247e59. [82] Grest GS, Lacasse MD, Kremer K, Gupta AM. J Chem Phys 1996;105: 10583e94. [83] Larson RG. Macromolecules 1994;27:4198e203. [84] den Otter WK, Briels WJ. J Chem Phys 2003;118:4712e20. [85] Guo HX, Kremer K. J Chem Phys 2003;119:9308e20. [86] Jo WH, Kim SH. Macromolecules 1996;29:7204e11. [87] Liu GL, Stoykovich MP, Ji SX, Stuen KO, Craig GSW, Nealey PF. Macromolecules 2009;42:3063e72. [88] Teubner M, Strey R. J Chem Phys 1987;87:3195e200. [89] Zhou N, Lodge TP, Bates FS. J Phys Chem B 2006;110:3979e89. [90] Dadmun MD, Muthukumar M, Schwahn D, Springer T. Macromolecules 1996; 29:207e11. [91] Degennes PG, Taupin C. J Phys Chem 1982;86:2294e304.