Quantifying Fluctuation/Correlation Effects on the Order-Disorder Transition of Symmetric Diblock Copolymers

0.2 MB·English

by  Jing Zong

Quantifying Fluctuation/Correlation Effects on the Order-Disorder Transition of
Symmetric Diblock Copolymers
Jing Zong and Qiang Wang∗
Department of Chemical and Biological Engineering,
Colorado State University, Fort Collins, CO 80523-1370
(Dated: January 24, 2012)
Using fast off-lattice Monte Carlo simulations with experimentally accessible fluctuations, we re-
port thefirst systematicstudyunambiguously quantifyingtheshift of theorder-disorder transition
(ODT) χ∗ of symmetric diblock copolymers from the mean-field prediction χ∗MF. Our simulations
2 are performed in a canonical ensemble with variable box lengths to eliminate the restriction of pe-
1 riodic boundary conditions on the lamellar period. Exactly the same model system (Hamiltonian)
0 is used in both our simulations and mean-field theory; the ODT shift is therefore due to the fluc-
2 tuations/correlations neglected by thelatter. While χ∗/χ∗MF 1 ¯−k is found with ¯ denoting
the invariant degree of polymerization, k decreases around th−e ¯∝-vNalue corresponding Nto the close
n N
packingof polymer segments as hard spheres, indicating theshort-range correlation effects.
a
J
PACSnumbers: 64.60.Cn,64.70.km
2
2
] Owing to the well developed polymer self-consistent vs. theinvariantdegreeofpolymerizationN¯ ≡nRe3,0/V,
t field(SCF)calculations[1],goodorevenquantitativeun- where n denotes the number of copolymer chains, Re,0
f
o derstanding has been achieved for the self-assembly of the end-to-enddistanceofanidealchain,andV the sys-
s flexiblelineardiblockcopolymer(DBC)meltsinbulk[2]. tem volume; these are the most accurate data obtained
t. Due to its mean-field approximation, however,SCF the- using each method as explained below.
a
orygivesqualitativelyincorrectpredictionsinthe region Beardsley and Matsen performed conventional Monte
m
near the order-disorder transition (ODT) where the sys- Carlo (MC) simulations on a face-centered cubic (FCC)
- tem fluctuations it neglects become important. In par- lattice in a canonical ensemble with replica-exchange to
d
n ticular, it fails to capture the fluctuation-induced first- study symmetric DBC of N = 30 segments modeled
o orderphasetransitionforODTofsymmetricDBC[3]and by the self- and mutual-avoiding walks with nearest-
c the direct transition between the gyroid and disordered neighbor repulsion between A and B segments ǫ (in
AB
[ phases[4]. Although the former is a classic problem in units of kBT, where kB is the Boltzmann constant and
1 polymer science and has been extensively studied by ex- T the thermodynamic temperature); about 20% lattice
v periments, theories and simulations, no quantitative un- sites were unoccupied, treated as anathermal solvent.[5]
5 derstandingofthefluctuationeffectsonODTisachieved. From the peak of constant-volume heat capacity, they
3 Fig. 1 summarizes recent simulation results on the ODT determined ODT to be χ∗N zǫ∗ N = 40.5, where
5 shift of symmetric DBC from the mean-field prediction z =12 is the lattice coordinati≡on nuAmBber.[5] The mean-
4
field ODT was determined to be χ∗ N = 14.654 using
. MF
1 lattice SCF calculations based on the same Hamiltonian
0 as used in their simulations.[6]
2 2
Using fast off-lattice Monte Carlo (FOMC) simula-
1
: 1 tions in an isothermal-isobaric ensemble, de Pablo and
v
co-workers studied compressible symmetric DBC melts
i
X 1 modeled by discrete Gaussian chains (DGC) of N = 64
r MF with a position-independent but anisotropic pair poten-
a tial of cubic symmetry.[7] They determined χ∗N by cal-
Ref. [5]
Ref. [7] culating the chemical potentials of the disordered and
0.1 Ref. [9] lamellarphases.[7]Notethattheirχ∗N wasrenormalized
Ref. [12]
to take into account the short-rangecorrelations(due to
Ref. [14]
0.04 the finiteinteractionrangeusedintheirsimulations)ap-
1 2 3 4 5 6
10 10 10 10 10 10 proximately;theirχ∗ N wasthereforedeterminedusing
MF
the random-phaseapproximationfor DGC with Diracδ-
function interactions.[7, 8]
FIG. 1: (Color online) Some literature results on the ODT Mu¨ller and Daoulas performed single-chain-in-mean-
shift. See text for details. field(SCMF)simulationstostudycompressiblesymmet-
ric DBC melts modeled by DGC of N = 32, and deter-
mined χ∗N =13.65 0.1 by comparing the free energies
±
of the disordered and lamellar phases.[9] Their SCMF
∗Electronicaddress: [email protected] simulation is similar to canonical-ensemble FOMC sim-
2
ulation with a spatial discretization scheme (i.e., with
20000
a position-dependent and anisotropicpair potential)[10],
PS-P2VP
except that a second-order term in the energy differ- PS-PI
10000
ence due to MC trial moves is neglected by the quasi- PS-PB
PS-PMMA
instantaneous field approximation[11]. They also esti-
PEP-PEE
mated χ∗ N = 11.7 0.5 by fitting the lamellar com-
MF ±
position profiles (averagedover directions parallelto the
lamellar interfaces) of a system at ¯ 2.5 109 at var-
N ≈ ×
ious χN to a sinusoidal function and then extrapolating
the square of the so-obtained amplitude of the composi- 1000
tion profiles to 0.[9]
Using field-theoretic simulations in a canonical en- 500
10 20 40 60 80 100
semble, Fredrickson and co-workers studied incompress-
ible DBC melts modeled by continuous Gaussian chains L0 (nm)
with Dirac δ-function interactions.[12] They determined
χ∗N = 12.2 (for slightly asymmetric DBC with a com- FIG.2: (Coloronline)Theinvariantdegreeofpolymerization
position f = 0.49) by comparing the free energies of the ¯ vs. thebulklamellarperiodL0 fornearlysymmetricDBC
N
disordered and lamellar phases.[12] For this “standard” melts commonly used in experiments.[17, 18]
model (which cannot be directly used/tested in molec-
ular simulations), χ∗ N = 10.495 for f = 0.5 is the
MF
well-known Leibler’s result[13] and is used in Fig. 1. longer holds and N becomes a chain discretization pa-
Finally, the prediction of the fluctuation theory[14] by rameter that does not correspond to the actual chain
FredricksonandHelfand(FH),χ∗/χ∗ 1=3.91 ¯−1/3, length used in experiments. This point is crucial for un-
MF− N
is also shown as a straight line in Fig. 1. This theory is derstanding coarse-grained models with soft potentials,
basedontheOhta-KawasakieffectiveHamiltonian[15]for the use of which is the basic idea of the recently pro-
the above“standard”modelandthe Hartreeanalysisby posed fast MC simulations[10, 19].
Brazovskii[16]. Althoughthelatterisrigorouslyaccurate Inthis work,we performextensiveFOMC simulations
only for ¯ &1010[14], the FH prediction has often been to quantify the ODT shift of symmetric DBC melts
comparedNwith simulation results at much smaller ¯ in with the following Hamiltonian: = C + E,
the literature. With the different models and metNhods where C = (3k T/2a2) n N−H1(R H RH)2
H B Pk=1Ps=1 k,s+1− k,s
used and the scarce simulation data in Fig. 1, it is clear is the Hamiltonian due to chain connectivity, and
thattheODTshiftofsymmetricDBCduetofluctuations E = (1/2κρ ) drdr′[ρˆ (r)+ρˆ (r)]u (r
0 A B 0
H R | −
is a problem far from being well understood. r′ )[ρˆ (r′)+ρˆ (r′)] + (χ/ρ ) drdr′ρˆ (r)u (r
A B 0 A 0
The system fluctuations are controlled by ¯, as sug- r′|)ρˆ (r′) is the Hamiltonian duRe to non-bonde|d in−-
B
N |
gested by the FH theory. For monodisperse DBC melts teractions; here a is the effective bond length, R
k,s
where each copolymer chain consists of N monomers denotes the spatial position of the sth segment on
A
of type A followed by N monomers of type B, with the kth chain, κ and χ are the generalized Helfand
B
R2 = N a2 + N a2 and nN = ρ¯ V, we have ¯ = compressibility[20] for the copolymer melts and the
e,0 A A B B 0 N
Nρ¯2a¯6; here N =N +N , ρ¯ φ ρ +(1 φ )ρ , generalized Flory-Huggins interaction parameter for
a¯ 0 fa2 +(1 fA)a2, Bf 0 ≡N /AN0,,Aa (a−) dAeno0t,Bes the repulsion between A and B segments, respectively,
the≡stpatistAical seg−mentBlengt≡h ofAthe A A(B) bBlock, ρ0,A ρ0 ≡ nN/V, ρˆA(r) ≡ Pnk=1PNs=A1δ(r − Rk,s) and
(ρ ) the A (B) monomer number density, and φ the ρˆ (r) n N δ(r R ) are the microscopic
0,B A B ≡Pk=1Ps=NA+1 − k,s
overall volume fraction of the A block. Fig. 2 shows number density of A and B segments at spatial position
¯ vs. the bulk lamellar period L for nearly sym- r, respectively, and u (r) = (15k T/2πσ3)(1 r/σ)2
0 0 B
N −
metric DBC melts commonly used in experiments[17]; for r < σ and 0 otherwise is a normalized (i.e.,
the range of L0 = 10 100nm, which is of interest dru0(r) = kBT) isotropic pair potential depending
∼ R | |
for most applications of DBC, roughly corresponds to only on the distance r between two segments with σ
¯ =500 20,000[18]. denoting the finite interaction range. Note that u (r) is
0
N ∼
In all conventional molecular simulations with hard essentially the same as the potential for the conservative
excluded-volume interactions (e.g., the Lennard-Jones force used in dissipative particle dynamics[21]; its
potential or the self- and mutual-avoiding walks), ¯ is Fourier transform is positive definite, thus avoiding the
N
on the same order of magnitude as N, which is at most formation of clustered crystals at large ρ0[22]. Taking
one hundred or so for concentrated polymer solutions Re,0 = √N 1a as the length scale, we have five
or melts (e.g., the leftmost data point in Fig. 1 is at parameters in−our model: ¯, χN, N/κ, N, and σ/a.
¯ 139). Using soft potentials thatallow particle over- OursimulationsareperfNormedinacanonicalensemble
N ≈
lapping is therefore the only way to study DBC melts withtrialmovesofrandomhopping,reptation,pivot[23],
within the above ¯-range (i.e., with experimentally ac- and box-length change. Note that the highly efficient
N
cessible fluctuations) at present, where nN = ρ¯ V no pivot algorithm cannot be used in multi-chain simula-
0
3
tions with hard excluded-volume interactions due to its
0.5 0.05
extremely small acceptance rates. With soft potentials, 8
however,wecanachieveanacceptancerateofabout70% 6
here. We also use replica-exchange[24] at different χN 0.4 ) 0.04
(with acceptance rates of 50 80%) to further improve P(b4 )
∼ 2 N
our sampling efficiency. 0.3 0.03
In simulations of periodic structures such as lamellae, 00.0 0.1 0.2( )0.3 0.4 d(
the periodic boundary conditions limit the allowed ori- 0.2 0.02
entations of the structure and thus its period. In par- d
ticular, for lamellae with a normal direction n in a sim- 0.1 d d( N) 0.01
ulation box with length L in the j(= x,y,z) direction,
j
L j·n=n L(n) must be satisfied, where n is the num- N
j j j 0.0 0.00
ber of periods contained in the box along the j direc- 50 60 70 80 90 100
tion (which could be 0) and L the lamellar period; this
N
gives L(n) = 1 (n /L )2.[25] For a fixed-length
.qPj j j
boxascommonlyusedincanonical-ensemblesimulations, FIG.3: (Coloronline) Ensembleaverageoftheorderparam-
eter, Ψ , and d Ψ /d(χN). The error bar of Ψ is taken
boththelamellarorientationanditsperiodaretherefore h i h i h i
as three times the standard deviation with sample correla-
discretized (i.e., L can hardly be L ); this problem is
0 tions taken into account. The vertical line marks the ODT
the most severe for cubic boxes.[25] To eliminate it, we ∗
χ N determined using the equal-weight criterion[27], which
change box lengths at a fixed V (i.e., using rectangular
isveryclosetothelocation ofthemaximumofd Ψ /d(χN)
boxes) with segmental positions re-scaled according to (i.e., another way to determine ODT[28]). The hinsiet shows
the Metropolis acceptance criterion. Compared to sim- the reweighted histogram of Ψ at χ∗, Pb(Ψ;χ∗), where the
ulations in an isothermal-isobaric ensemble, this has the vertical line marks Ψ (χ∗). N = 10, σ/a = 0.3, N/κ = 50,
advantage that ¯ is constant during each simulation. and ¯ =10,000. h i
N N
We use a new order parameter for ODT, which char-
acterizesthedegreeofpositionalorderinlamellae. Fora
∗
given direction n, we calculate the volume fraction of TABLE I: Mean-field ODTχMFN.
A segments φA(t) ρˆA(t) [ρˆA(t) + ρˆB(t)] as a func- σ/a=0.3 σ/a=0.1√19 σ/a=2/√3
≡ (cid:14)
tion of position t along n in a collected configura- N =10 10.047 11.427
tion (averaged over directions perpendicular to n), and N =20 10.405 10.462 11.102
ψ(n) dtexp[4πit/L(n)]f(t) dtf(t) with f(t)
≡ (cid:12)R (cid:12)(cid:14)R ≡
1 2φ (cid:12)(t) 1 and i √ 1.(cid:12) The order parameter
A
− | − | ≡ −
Ψ ∈ [0,1] for a collected configuration is defined as the sults are independent of N/κ and ¯, and can be under-
largest ψ over all possible lamellar orientations in the N
stood based on our previous work[10, 29].
box. Fig. 3 shows the ensemble average of Ψ, Ψ , and
h i Collectingallthedata,Fig.4showsourresultsofODT
dhΨi/d(χN)vs. χN. Wefurtherusemultiplehistogram shift vs. ¯ at various N, σ/a and N/κ.[30] We find
reweighting[26] to accurately locate ODT according to χ∗/χ∗ N1 ¯−k in all the cases. While this is con-
the equal-weight criterion[27] described below. MF− ∝ N
sistent with the functional form of FH prediction, our
The inset of Fig. 3 shows the histogram of Ψ at a
ODTshift islargerthantheir predictionforallthe cases
1
given χ, Pb(Ψ;χ), with R0 dΨPb(Ψ;χ) = 1. The double we have studied. We also find a decrease in the nega-
peak near ODT is a signature of the first-order phase tive exponent k around ¯ = 2(N 1)3 N2(σ/a)6 in
cp
transition. The ODT is determined as χ∗ at which N − (cid:14)
the cases of σ/a = 0.3 and 0.1√19 0.44, which corre-
R0hΨidΨPb(Ψ;χ∗) = Rh1ΨidΨPb(Ψ;χ∗). To estimate the spondsto the FCC closepackingofp≈olymersegmentsas
statistical error of χ∗, we calculate ODT using the first- hard spheres. This is therefore due to the local packing
and second-half of our samples collected after equilibra- ofsegments(i.e.,theshort-rangecorrelations). Denoting
tion, respectively, and take three times their largest de- the k-value for ¯ < ¯ by k and that for ¯ > ¯
cp 1 cp
viation from χ∗ determined using all the samples as the by k , we find thNat kN k decreases with incrNeasingNN
2 1 2
−
error bar. Note that, even with the same order param- (i.e., in the case of σ/a=0.3), which leads to decreasing
eter, there are other ways to determine the transition σ (thus weaker correlations) at constant R . One may
e,0
point[28],whichgiveslightlydifferentODTforfinitesys- therefore expect a single power-lawdecayof χ∗/χ∗ 1
tems as shown in Fig. 3. Systematic study of the finite- with increasing ¯ in the limit of N (i.e., σMF−0),
N → ∞ →
sizeeffectsrequiresmuchmoreextensivesimulationsand which has no such correlation and is the case studied by
is out of the scope here. FredricksonandHelfand[14]. Ontheotherhand,increas-
WiththesameHamiltonianasusedinourFOMCsim- ing σ/a at constantN leads to strongercorrelationsand
ulations, we determine the mean-field ODT χ∗ using thus increasing k k (i.e., in the case of N =20). We
therandom-phaseapproximation[10]. TableIlisMtsFχ∗ N also note that va1ry−ing2N at constant σ/a exhibits both
MF
for various N and σ/a. We note that the mean-field re- thecorrelationandchaindiscretizationeffects. Toexam-
4
100. Thislargeσ/a-valuealsomakessimulationsat ¯ >
8 N
2000 expensive to do. Nevertheless, we see from Fig. 4
that increasing N increases the ODT shift in this case,
1.3
in contrast to that of σ/a = 0.3. On the other hand, as
0.75
1 0.8 we increase N/κ from 0 to 50, k decreases here and the
0.6 2
1 0.75 two lines with the same N in Fig. 4 cross at ¯ 1300
1.2 cr
MF 1 0.5 (for N = 10) or 2500 (for N = 20). TherNefore≈, χ∗N
0.65 0.6 decreases for ¯ < ¯ and increases for ¯ > ¯ with
cr cr
N N N N
increasing N/κ. While the latter is consistent with the
0.45 0.5
a 0.3 0.44 1.15 case of σ/a = 0.3, we find that both k1 and k2 increase
N 10 Open: N 50 with increasing N/κ there.
N 20 Filled: N 0
0.1 To summarize, we have performed extensive FOMC
2 3 4 5
10 10 10 10 simulations with experimentally accessible fluctuations
to systematically and unambiguously quantify the ODT
shift of symmetric DBC from the mean-field prediction.
Exactly the same model system (Hamiltonian) is used
FIG. 4: (Color online) Our results on the ODT shift. The
negativeslopekisgivennexttoeachline. Seetextfordetails. in both our simulations andmean-field theory; the ODT
shift is therefore due to the fluctuations/correlationsne-
glected by the latter. We hope this work will stimulate
inethelatteralone,onemaycomparethecaseofN =10 the development of advanced theories better describing
and σ/a = 0.3 with that of N = 20 and σ/a = 0.1√19 the fluctuation/correlationeffects identified here.
(both of which have σ = 0.1Re,0). We expect that the WethankProf.MarcusMu¨ller,Prof.MarkStoykovich,
chaindiscretizationeffectsdiminishquicklywithincreas- and Dr. Kostas Daoulas for helpful discussions. This
ingN,assupportedbyFig.6(a)ofRef.[10]showinghow work was supported by the U.S. Department of En-
χ∗ N varies with N. ergy, Office of Basic Energy Sciences, Division of Mate-
MF
In the case of σ/a = 2/√3 1.15, ¯ < 15 for both rials Sciences and Engineering under Award DE-FG02-
cp
N-values, and we did not per≈form simNulations at ¯ < 07ER46448.
N
[1] G. H. Fredrickson, The Equilibrium Theory of Inhomo- and A.Skoulios, ibid.15, 258 (1982); P.F. Green et al.,
geneous Polymers, Oxford University Press, 2006. ibid.22,2189(1989); C.M.Papadakiset al.,J.Phys.II
[2] M.W.Matsen andM.Schick,Phys.Rev.Lett.72,2660 France7,1829(1997);K.Almdaletal.,Phys.Rev.Lett.
(1994); C. A. Tyler and D. C. Morse, ibid. 94, 208302 65, 1112 (1990).
(2005). [18] ¯ is quite sensitive to a; see Supplemental Material for
[3] F. S.Bates et al.,Phys. Rev.Lett. 61, 2229 (1988). Nthe valuesof a and ρ0 for each block.
[4] F. S.Bates et al.,Faraday Discuss. 98, 7 (1994). [19] Q. Wang, Soft Matter 5, 4564 (2009); ibid. 6, 6206
[5] T.M.BeardsleyandM.W.Matsen,Eur.Phys.J.E32, (2010).
255 (2010). [20] E.HelfandandY.Tagami,J.Polym.Sci.,PartB:Polym.
[6] O.N.Vassiliev andM.W.Matsen, J.Chem.Phys.118, Lett. 9, 741 (1971).
7700 (2003). [21] R.D.GrootandP.B.Warren,J.Chem.Phys.107,4423
[7] F. A. Detcheverry et al., Phys. Rev. Lett. 102, 197801 (1997).
(2009). [22] C. N.Likos et al.,Phys.Rev. E63, 031206 (2001).
[8] Their result of χ∗MFN = 10.3 is incorrect and should be [23] S. D. Stellman and P. J. Gans, Macromolecules 5, 516
10.481 instead, which is used in Fig. 1. (1972).
[9] M. Muller and K. C. Daoulas, J. Chem. Phys. 128, [24] R. H. Swendsen and J. S. Wang, Phys. Rev. Lett. 57,
024903 (2008). 2607 (1986).
[10] Q.WangandY.Yin,J.Chem.Phys.130,104903(2009). [25] Q. Wang et al.,J. Chem. Phys.112, 450 (2000).
[11] K. C. Daoulas and M. Muller, J. Chem. Phys. 125, [26] A. M. Ferrenberg and R. H.Swendsen, Phys. Rev. Lett.
184904 (2006). 63, 1195 (1989).
[12] E.M.Lennonetal.,Phys.Rev.Lett.101,138302(2008). [27] M. Muller and N. B. Wilding, Phys. Rev. E 51, 2079
[13] L. Leibler, Macromolecules 13, 1602 (1980). (1995).
[14] G. H. Fredrickson and E. Helfand, J. Chem. Phys. 87, [28] See, for example, A. M. Ferrenberg and D. P. Landau,
697 (1987). Phys. Rev. B 44, 5081 (1991); M. S. Challa et al., ibid.
[15] T. Ohta and K. Kawasaki, Macromolecules 19, 2621 34, 1841 (1986).
(1986). [29] Q. Wang, J. Chem. Phys.129, 054904 (2008).
[16] S.A. Brazovskii, Sov.Phys. JETP 41, 85 (1975). [30] See Supplemental Material for a list of our simulation
[17] Y.Matsushitaetal.,Macromolecules23,4313(1990);T. parameters and results.
Hashimotoetal.,ibid.13,1237(1980);G.Hadziioannou