Self-Assembly of Cylinder-Forming Diblock Copolymer Thin Films

Aug 5, 2013 - We study the self-assembly of cylinder-forming diblock copolymers confined in ultrathin films by means of dissipative particle dynamics ...
0 downloads 0 Views 2MB Size
Article pubs.acs.org/Macromolecules

Self-Assembly of Cylinder-Forming Diblock Copolymer Thin Films Arash Nikoubashman,* Richard A. Register, and Athanassios Z. Panagiotopoulos Department of Chemical and Biological Engineering, Princeton University, Princeton, New Jersey 08544, United States ABSTRACT: We study the self-assembly of cylinder-forming diblock copolymers confined in ultrathin films by means of dissipative particle dynamics simulations paired with an elaborate polymer model that mimics the experimental reference system (PS−PHMA) as faithfully as computationally feasible. In particular, the effects of the polymer composition, film thickness, and surface interactions are investigated systematically, where we achieve quantitative agreement with experiments. We observe that the cylindrical domains change their orientation from perpendicular to parallel with respect to the substrate twice as the film thickness is increased. Furthermore, the onset thickness and the sharpness of this transition strongly depend on the microscopic details of the copolymers and their interaction with the substrate. fraction f was increased.10,11 We achieve almost perfect quantitative agreement with these experiments and provide a detailed explanation for the observed phenomena. Such surface-induced phase transitions have previously been studied through self-consistent field theory (SCFT)9,23 and dynamic density functional theory (DDFT).14−16 Here, the seminal work by Horvat et al.16 is especially noteworthy, where the phase behavior of cylinder-forming triblock copolymers was extensively analyzed via DDFT calculations and good qualitative agreement with experiments was achieved. However, this and other prior simulation studies lack a direct quantitative comparison with specific experiments. In contrast to these mean-field approaches, we employ a particle-based simulation model that can describe irregular structures and capture fluctuations. Because of the microscopic nature of our system and the large system sizes (32 beads per chain and up to 32 400 polymers), we are also able to study the effects of polydispersity, which has not been done before. Furthermore, we introduce a novel quantification of order, which facilitates direct comparisons with experimental data. Such large-scale simulations are only feasible by harnessing the power of modern massively parallel graphics processing units (GPUs). As a comparison, the presented simulations would have taken about 95 processor years on a conventional computer cluster. In addition to the simulations, we also develop here a simple theoretical model that successfully captures the most distinctive feature of the system, i.e., the sequential turning of the microdomain orientation with increasing film thickness. Here, our approach distinctly differs from previous models,9,23 since it is based solely on the microscopic properties of the system and plain geometric considerations. From our model calculations we were able to conclude that the driving force behind the || to

1. INTRODUCTION Thin block copolymer films are highly relevant for many scientific and industrial applications due to their ability to form uniform domains of controllable shape at nanometer length scales. Significant efforts have been devoted to understanding, predicting, and controlling the underlying self-assembly mechanisms.1−6 For diblock copolymer melts in the bulk, the emerging morphologies are dictated almost exclusively by the volume fraction of the minority block f and the strength of monomer interactions, quantified by the Flory−Huggins parameter χ.7 However, this picture changes completely for ultrathin films, where the film thickness H and the detailed copolymer−substrate interactions play a major role. Furthermore, one also needs to consider the relative orientation of the microdomains with respect to the substrate. In fact, it has been shown in both experiments8−11 and simulations12−16 that, for lamella- and cylinder-forming copolymers, the orientation of the microdomains switches from perpendicular (⊥) to parallel (||) (and vice versa) as the film thickness or the surface properties are changed. From a technological point of view, films with parallel cylinders may be of interest in the fabrication of in-plane nanowire arrays,17−19 while upright cylinders and spheres could have potential applications in the patterning of hexagonal arrays for data storage.20,21 In this paper, we report on the phase behavior of confined cylinder-forming diblock copolymers using large-scale computer simulations. We employ a sophisticated coarse graining protocol22 in combination with dissipative particle dynamics (DPD) to study ultrathin films of polystyrene-b-(n-hexyl methacrylate) (PS−PHMA). These diblock copolymers are promising candidates for the fabrication of polarizing grids in the ultraviolet spectrum and show a very rich phase behavior.17−19 Recent experiments systematically analyzed the phase behavior of these systems as a function of film thickness and polymer composition, revealing the existence of || to ⊥ transitions, which became less prominent as the volume © XXXX American Chemical Society

Received: April 26, 2013 Revised: June 6, 2013

A

dx.doi.org/10.1021/ma400867s | Macromolecules XXXX, XXX, XXX−XXX

Macromolecules

Article

⊥ transition must be enthalpic in nature and that entropy plays only a subordinate role. The rest of this article is organized as follows. First, we present our copolymer model and a detailed description of the coarse-graining procedure. Additionally, we briefly introduce the employed simulation technique. In the subsequent section, we present our simulation results and compare them to experimental findings; first we study the effects of the copolymer composition and size and then discuss the influence of polydispersity. In the second half of the Results section, we introduce our theoretical model and compare its predictions with the results from the experiments and simulations. Finally, we summarize the findings and draw our conclusions. We also provide an outlook on how our approach can be used for sheared systems.

Here, Re is the end-to-end distance of the copolymer in the melt, which can be estimated via the freely jointed chain model as Re2 = ∑N−1 bi2.32 Furthermore, √N̅ = ρRe3/N denotes the i so-called interdigitation number, which is a measure for the strength of fluctuations.33 The incompatibility of the type A and B monomers is quantified by the Flory−Huggins parameter χ, while κ sets the compressibility of the melt.22 Lastly, ω(r) is the density cloud around each particle’s center. In addition to these intramolecular forces, the bead−surface interactions need to be taken into account as well. In our setup, we have applied periodic boundary conditions in the x and y directions of the system and placed walls at ±(H + σ)/2 in the xy-plane. The direction of the exerted force is then perpendicular to this plane, and its magnitude depends only on the position of the bead relative to each of the walls. In our simulations, we employ the following wall potential:

2. MODEL AND SIMULATION METHOD For studying the behavior of block copolymers in thin films and to make meaningful comparisons to real experiments, we need a theoretical model that is as detailed as computationally feasible. There is a wide variety of different approaches to choose from (each with its own advantages and drawbacks), such as the cell dynamics (CD)3,24,25 and molecular dynamics (MD) techniques.26,27 In CD simulations, the time evolution of an order parameter is evaluated on a lattice, where the driving forces are in general due to chemical potential gradients and diffusive dynamics. This drastic coarse-graining enables simulations over long time scales and for large system sizes but prohibits the detailed analysis of microscopic properties, such as the polymer stiffness and polydispersity. MD simulations of bead−spring models, on the other hand, retain the microscopic details of the system to a much higher degree but quickly become computationally infeasible. To bridge the gap between the two aforementioned approaches, we have opted to employ the so-called theoretically informed coarsegrained (TICG) model22 in conjunction with the dissipative particle dynamics simulation technique;28−31 the bead−spring model enables us to study the influence of microscopic properties, such as polydispersity, while the inherent softness of the DPD particles allows for rather large time steps. In the TICG model, each of Np polymer chains is discretized into N beads of diameter σ, where we denote the position and velocity of the ith particle as ri and vi, respectively. The interaction between individual beads can be split up into a contribution due to the chemical bonding between beads of the same chain (Ub) and into a second part covering the interactions between all particles (Umm). In this model, the bonds between two adjacent beads (with indices i and j) are described by Hookean springs: Ub(rij) kBT

2 3 ⎛ rij ⎞ = ⎜ ⎟ 2 ⎝ bi ⎠

⎛ z ⎞−6 Uwall(z) = εS⎜ ⎟ kBT ⎝ σS ⎠

where z denotes the (shortest) distance to each wall. In eq 3, the strength and the range of the repulsion can be controlled through the parameters εS and σS. Whereas the above-described rules governing the polymer model are general, the simulation of specific copolymer systems requires special care and will be discussed in what follows. In this work, we follow a bottom-up coarse-graining approach to model thin films consisting of linear PS−PHMA diblock copolymers, for which experimental data are available.10,11 More specifically, we will focus on the cylinder-forming diblocks PS−PHMA 31/90 and PS−PHMA 24/89, where in this notation the first number signifies the molar mass (in kg/ mol) of the first block and so forth. Table 1 shows the relevant Table 1. Molecular Weight of the Repeat Unit M0, Molecular Weight of the Kuhn Segment MK, Mass Density ρ, and Kuhn Length b of the Constituent Styrene (S) and Hexyl Methacrylate (HMA) Moleculesa

kBT

=

N̅ [κN + χN (1 − δij)] R e3

(r − rj) dr

monomer

M0 [kg/mol]

MK [kg/mol]

ρ [g/cm3]

b [nm]

S HMA

0.104 15 0.170 25

0.720 1.62

0.969 0.95

1.8 2.44

a

The properties of the S and HMA monomers were obtained from ref 32 and refs 34−37, respectively.

properties of the constituent monomers (styrene and hexyl methacrylate). We can then use this information to determine several physical properties of the studied copolymers, such as the number of monomer units N and the volume fraction of the minority block f (see Table 2). For the computational model, we will refer to the different particle types simply as A and B in lieu of S and HMA for the

(1)

Table 2. Monomer Number of Type A, NA, and of Type B Particles, NB, Volume Fraction of the Minority Block f, Endto-End Length Re, and Interdigitation Number √N̅ 1/2 of the Copolymers

where kB denotes the Boltzmann constant, T the temperature, rij the distance between the two particles, and bi the bond length of the connecting spring. The monomer−monomer interaction Umm between two particles i and j of either type A or B is given by22 Umm(rij)

(3)

∫V ω(r − ri)ω (2) B

copolymer

NA

NB

f

Re

√N̅ 1/2

PS−PHMA 24/89 PS−PHMA 31/90 A−B 7/25 A−B 8/24

230 298 7 8

523 529 25 24

0.21 0.25 0.22 0.25

20.8 nm 21.7 nm 7.2σ 7.1σ

46.1 48.4 46.1 48.4

dx.doi.org/10.1021/ma400867s | Macromolecules XXXX, XXX, XXX−XXX

Macromolecules

Article

sake of clarity and to underline the general nature of our findings.14−16 With the information on the monomeric level at hand, we can derive the properties of the copolymer as a whole and subsequently perform the coarse-graining. For the latter step, we have chosen to use N = 32 beads and to keep both the parameters f and √N̅ 1/2 invariant during the mapping procedure. Furthermore, we set in our model the bond length between A beads to bA = σ and as a result the bond length between B particles to bB = 1.4σ. Table 2 shows the properties for the experimental copolymers as well as for the coarsegrained ones. From these data, we get the conversion factor σ ≈ 2.9 nm for A−B 7/25 and σ ≈ 3.1 nm for A−B 8/24. Furthermore, we can conclude that we need to conduct our simulations at a particle density of ρ = √N̅ N/R3e ≈ 4.5 in order to faithfully represent the original experimental systems.22,38 Finally, we still need to derive the remaining interaction parameters presented in eqs 2 and 3, namely χ, κ, and ω(r). For the Flory−Huggins parameter χ, the following relationship has been determined for PS−PHMA block copolymers in the temperature range 423 K < Tabs < 488 K (with absolute temperature Tabs):39 χ ≈ 0.035 + 3.93 K/Tabs

Figure 1. Pair potential between particle pairs of the same type UAA and dissimilar type UAB as a function of interparticle distance r.

(4)

The reference experiments for the PS−PHMA copolymers have been conducted at Tabs = 423 K,10,11 and therefore we have used χ ≈ 0.044. For the compressibility term in Umm, we have chosen κN = 15 according to the arguments brought forward in refs 22 and 38. The last parameter that needs to be specified is the functional form of ω(r), i.e., the density cloud around each particle. Here, we chose spherical steplike functions, where ω(r) = 1 for r < σ/2 and ω(r) = 0 elsewhere, so that the beads can be considered as soft spheres. Hence, the integral in eq 2 reduces to the normalized overlap volume of two spheres at a distance of r. We then utilized the parameters derived from PS−PHMA 24/89 for all our simulations, instead of using an extra set of parameters for modeling PS−PHMA 31/90. This step is justified, since the relevant parameters (and thus the resulting pair potentials) differ only marginally between both copolymers. The advantage of this approach is that the number of model parameters is reduced significantly, and thereby the comparison between the simulations is facilitated. The resulting pair potentials between particle pairs of the same type (AA or BB) and dissimilar type (AB) are shown in Figure 1. Now, in order to make a (quantitative) study of the behavior in thin films, we need to replicate the appropriate copolymer− surface interactions. In our model, the interaction between each bead type (A or B) and each surface (lower or upper boundary) is characterized by two parameters (cf. eq 3), resulting in a total of eight free parameters that need to be determined. In the experiments, the films are spin-coated atop of a silicon substrate, and it was observed that the PHMA (matrix) block wets both the substrate and the air surface.19 Thus, we can reasonably impose the same interaction potential Uwall for the upper and lower boundaries, thereby cutting the number of free parameters in half. We set εS,B = 3.0 kBT for the majority block and εS,A = 8.0 kBT for the minority block to replicate the wetting behavior observed in the experimental setup. Moreover, we decided to utilize the same length scale (reflected by the parameter σS) for both particle types, and chose σS = σ for the sake of simplicity. Figure 2 shows the resulting potentials.

Figure 2. Wall potential for A and B particles as a function of wall distance z.

The time evolution of the particle trajectories was obtained by the standard velocity Verlet algorithm,40,41 where we have chosen a value of Δt = 0.02 for the time step. In our simulations, NVT conditions were established via a DPD thermostat. The DPD model consists of the summation of conservative forces (FC), random noise (FR) and dissipative forces (FD). The conservative forces can be derived from eqs 1−3 through the relationship FC = −∇U. The random and dissipative forces on each particle are given by28−31 FR (rij) = −ξij 3

2γkBT ωD(rij)riĵ Δt

FD(rij) = −γωD2(rij)(vij·rij)riĵ

(5) (6)

where ξij denotes a uniformly distributed random number with zero mean and unit variance, γ the friction coefficient, r̂ij the unit vector between the beads i and j, and vij = vi − vj. Furthermore, the functional form of the weight function ωD(rij) has been chosen as C

dx.doi.org/10.1021/ma400867s | Macromolecules XXXX, XXX, XXX−XXX

Macromolecules

Article

Figure 3. Film morphologies obtained at H = 30 nm (a) and H = 40 nm (b). The left half of each image shows the experimental results for PS− PHMA 24/89, while the right half shows the findings from the simulations of A−B 7/25. Each snapshot shows a 400 nm × 400 nm section of the film, and the scale bars indicate 100 nm.

Figure 4. Fraction of parallel cylinders p(C∥) as a function of film thickness H retrieved from simulations (red, solid) and experiments (blue, dashed) for cylinder-forming diblock copolymer with f ≈ 0.21 (a) and f ≈ 0.25 (b).

⎧ ⎪1 − rij / σ , rij < σ ωD(rij) = ⎨ ⎪ rij ≥ σ ⎩ 0,

Ntotal

ρ (x , y ) =

∑ δ(r − ri)e(H /2 − z)/l

d

i

(7)

(8)

In eq 8, ld is an arbitrary decay length which we set to ld = 5σ by visually matching the resulting density profiles with the TMAFM images. We then measured the area taken up by perpendicular cylinders A⊥ (circular areas in Figure 3b) and defined p(C||) = 1 − A⊥/A, where A denotes the total area covered by styrene beads. Here we decided that a continuous area corresponded to the C⊥ phase, when its size Ai was below 100 σ2 and its circularity e = 4πAi/Si 2 exceeded a value of 0.7 (with perimeter Si and e = 1.0 for a perfect circle). This procedure was then repeated five times for each state point to achieve proper sampling of the data.

We then conducted our computer experiments by means of the HOOMD simulation package,42−44 which is highly optimized for GPUs. The system size was set to L × L × H, with L = 80σ and 6σ ≤ H ≤ 36σ, which corresponds to a periodicity of about 250 nm and film thicknesses between 20 and 110 nm. For such system sizes, the number of total interaction sites ranged between 172 800 ≤ Ntotal ≤ 1 036 800. We started all our simulations from a completely disordered, random configuration and equilibrated the system over 8 × 105 time steps. The microstructures were then analyzed following a scheme that mimics the experimental measurement protocol as closely as possible, in order to facilitate direct comparisons; in the experiments, the thin copolymer films were analyzed via tapping mode atomic force microscopy (TM-AFM), yielding two-dimensional images that typically covered 2 μm × 2 μm and contained a sizable number of microdomains (see Figure 3). The resulting images were then utilized to compute the fractional coverage of perpendicular cylinders.10,11 Accordingly, we first calculated the two-dimensional density distribution in the xy-plane, ρ(x,y), for our simulation data. Here, we replicated the finite depth perception capabilities of the TMAFM by weighting the contribution of each bead to ρ(x,y) on the basis of its distance to the surface:

3. RESULTS We begin our discussion with the detailed analysis of the cylinder orientation as a function of film thickness. In the experiments,10,11,17−19 thin films of PS−PHMA block copolymers were spin-coated on silicon substrates. The ensuing macroscopic film thickness H was measured via ellipsometry, while the microscopic morphologies were captured through TM-AFM. In our simulations, we can precisely control H and capture the shape of the microdomains, since we have access to the position of each single bead at any given time. Figure 3 shows a comparison between the morphologies obtained from experiments and simulations, and we can see good qualitative agreement between them; the simulations do not only predict D

dx.doi.org/10.1021/ma400867s | Macromolecules XXXX, XXX, XXX−XXX

Macromolecules

Article

Figure 5. p(C∥) as a function of film thickness H for various volume fractions f at fixed polymerization N = 32 (a) and various chain lengths N at fixed f = 0.25 (b).

near H ≈ 1.5Dhex and show almost full recovery to a parallel configuration at H ≈ 2Dhex. Following this argument, a subsequent C∥ → C⊥ transition can be expected at H = 2.5Dhex and so forth. Indeed, we observed a second minimum for p(C∥) located at H ≈ 70 nm, which is however much less pronounced than the preceding one. We can observe a very similar overall behavior for copolymers with f ≈ 0.25 in Figure 4b, where we can identify some small, but significant, differences: in this case, the first C∥ → C⊥ transition occurs at slightly thinner films, and the emerging parallel cylinders are stable over a somewhat larger thickness range. Furthermore, the subsequent minimum of p(C∥) is shifted to H ≈ 45 nm and is significantly more narrow. Finally, these copolymer films did not exhibit the second C∥ → C⊥ transition, observed for the f ≈ 0.21 case. These findings can be explained through the increased volume fraction f, which leads to an elevated number of A particles at the surfaces in the C⊥ phase and a concomitant energy penalty due to the selective nature of the copolymer−surface interactions (cf. eq 3 and Figure 2). In order to better understand the influence of f, we conducted additional simulations of cylinder-forming A−B 9/ 23 copolymers with f ≈ 0.28 and compared its behavior with the previous diblocks in Figure 5a. Here, we can immediately notice that both the width and the depth of the minimum decreased significantly as we increased the volume fraction to f ≈ 0.28. Furthermore, we observed predominantly parallel cylinders for very thin films and were not able to identify any minimum beyond the first one. These results are in line with the findings so far and corroborate our interpretation of the C∥ → C⊥ transition. We then studied the dependence of the transition on the chain length N, as shown in Figure 5b, where we can identify two distinct trends which occur upon decreasing the degree of polymerization: a shift of the minimum to smaller H and an overall increase of stable parallel cylinders for all values of H. The first effect can be easily understood, since shorter copolymers develop smaller microdomains, reflected in smaller values of d and D. Therefore, following the previous arguments, the C∥ → C⊥ transition and the ensuing recovery to the C∥ phase have to occur at smaller film thicknesses. The second effect is more difficult to comprehend, since there is no direct

the appropriate phases but also reproduce a similar size and shape distribution of the microdomains. From our simulation data, we determined both the diameter of the cylindrical domains and their periodicity as d = 17 ± 1 nm and D = 36 ± 1 nm for A-B 7/25 and d = 18 ± 1 nm and D = 39 ± 1 nm for A-B 8/24. These values are in good agreement with the experimental findings of d = 16 nm and D = 37 nm for PS−PHMA 24/89 and d = 17 nm and D = 42 nm for PS− PHMA 31/90, which have been determined via small-angle Xray scattering.10,11 The small, but sizable, discrepancies in D are most likely due to the employed value for the Kuhn length of HMA bB = 2.44 nm (cf. Table 1), which has been determined at Tabs = 306 K, while the experiments and simulations have been conducted at Tabs = 423 K.34−37 In order to make a more quantitative analysis, we measured the fraction of parallel and perpendicular cylinders, p(C||) and p(C⊥) = 1 − p(C||), respectively (see preceding section for details). Figure 4 shows a comparison between experiments and simulations for this quantity (cf. Table 2); we see remarkable quantitative agreement between the data up to H ≈ 45 nm. Past this point, the simulations underestimate the fraction of parallel cylinders consistently but still reproduce the correct qualitative behavior. Figure 4a shows that diblocks with f ≈ 0.21 formed preferably perpendicular cylinders in very thin films with H ≤ 20 nm. The orientation of the microdomains then changed from perpendicular to parallel as the film thickness was increased, until the system consisted almost exclusively of parallel cylinders at H = D. However, upon further increasing H, we observed a very sharp transition to upright cylinders at H ≈ 40 nm (cf. Figure 3). This phenomenon can be attributed, for the most part, to the incompatibility between the lattice spacing of the microdomains and the film thickness: the natural thickness of a layer of hexagonally packed cylinders is given by Dhex = D√3/2, which leads to a value of Dhex ≈ 32 nm for the experimental system and Dhex ≈ 31 nm for the model system. In the transient regime D < H < 2Dhex, the C|| phase is not necessarily the one with the lowest free energy, and other phases (or combinations of them) might be favorable. This geometric explanation for the C∥ → C⊥ → C∥ transition is consistent with the observation from experiments as well as simulations, since the respective p(C∥) have their minimum E

dx.doi.org/10.1021/ma400867s | Macromolecules XXXX, XXX, XXX−XXX

Macromolecules

Article

connection between the degree of polymerization and the relative amount of perpendicular microdomains. To better understand this behavior, we studied the time evolution of the films at those thicknesses where the C⊥ phase was expected to be preferential; at the beginning of the simulations, the systems were in a completely disordered state. After a very short time, perpendicular cylinders formed spontaneously. These domains, however, did not remain well separated from each other, and a significant amount of neighboring sites fused during the equilibration. We attribute this behavior to the decreased PHMA buffer size between the PS microdomains, which facilitates such merging processes. This interpretation is corroborated by the fact that the effect is amplified as N is decreased (see Figure 5b). As a consequence, mean-field calculations should not be able to predict such a behavior, since this effect is due to the finite size of the chains. Next we studied the influence of polydispersity, since in real experiments the degree of polymerization is not exactly the same for all copolymers in the system. Here, we focus on the A−B 8/24 copolymer and its experimental counterpart PS− PHMA 30/91, for which a polydispersity index of 1.08 was reported.11 We expect similar behavior for the other copolymer compositions, if polydispersity is taken into account. Figure 6

Figure 7. Simulation results for cylinder-forming A−B 8/25 diblock copolymers at film thicknesses H = 11σ (top row) and H = 14σ (bottom row) with different strengths of the bead−surface interaction. The position of the boxes indicates the parameter where the simulation has been conducted and its gray level and distribution the phases that were observed at that point. The snapshots show typical configurations of the disordered (D), perpendicular cylinder (C⊥), parallel cylinder (C∥), and perforated lamella (PL) structures. Type B particles are not shown here for the sake of clarity.

at H = 11 σ (top row) and H = 14 σ (bottom row) as a function of ΔεS = εS,A − εS,B. We chose these two film thicknesses, since for the experimental system at hand (cf. Figure 4a) they correspond to configurations of fully aligned parallel and perpendicular cylinders, respectively. At ΔεS = 0, the system formed a disordered phase (D) where no clear patterns were discernible. For slightly selective surfaces (ΔεS > 0), we were able to identify a well-ordered hexagonal array of perpendicular cylinders (C⊥). As ΔεS was increased even further, the cylindrical domains changed their orientation to parallel with respect to the substrate. Finally, the copolymers eventually formed perforated lamellae (PL) for the largest ΔεS values studied here. We observed a very similar sequence of phases for H = 14σ, albeit shifted to higher values of ΔεS. This progression is in very good qualitative agreement with the results for cylinder-forming triblock copolymers published in ref 16. To better understand the peculiar phase behavior of these thin films, we also develop here a simple model that captures the salient features of the system at hand. Our approach is solely based on geometric considerations and microscopic properties and is therefore far less involved than either the main DPD simulations of the present work or earlier SCFT calculations.9,23 The advantage of such a simple model is that it is very intuitive and that it highlights the prominent features of the problem. In what follows, we will focus on the C∥ to C⊥ transitions that occur as the film thickness H is increased. We first examine the state C∥ at H = D, shown schematically in Figure 8a. In this case, the thickness of the monolayer film matches exactly the periodicity of the microdomains and parallel cylinders with a circular cross section of diameter d

Figure 6. Comparison of p(C∥) between a monodisperse and a polydisperse system for f = 0.25 as a function of film thickness H. The inset shows the experimental probability distribution of the molecular weights M (solid line) and the corresponding values of N (symbols).

shows p(C∥) as a function of the film thickness H for both the mono- and polydisperse case, and besides some minor quantitative differences, which are well within the error bars, the two plots are indistinguishable. This behavior becomes more comprehensible if we consider the size distribution of the diblocks, p(N), which is shown in the inset of Figure 6; the chain lengths follow almost perfectly a Gaussian distribution with a mean of N = 32 and a standard deviation of σ ≈ 7. Because of the rather narrow and symmetric shape of p(N), the behavior of the polydisperse copolymer melt is dominated by the N = 32 population, whereas the contributions from the remaining polymer chains, which equally consist of more or fewer monomers than the average, cancel out. Finally, we performed a qualitative study on the effect of surface wetting by systematically varying the interaction strength between the walls and the particles. Figure 7 shows the emerging phases for confined A−B 7/25 copolymer melts

Figure 8. Schematic representation of three (distorted) cylindrical microdomains: (a) parallel microdomain with a circular cross-section at H = D, (b) parallel microdomain with an ellipsoidal cross section at H′ > H, and (c) perpendicular domain with circular cross section. The rectangular areas at the top and bottom surfaces in panels a and c highlight the difference between the two configurations (see text). F

dx.doi.org/10.1021/ma400867s | Macromolecules XXXX, XXX, XXX−XXX

Macromolecules

Article

parallel alignment was restored at H = 16.8σ for A−B 7/25, H = 15.5σ for A−B 8/24, and H = 14.5σ for A−B 9/23. This is in remarkable agreement with the findings from the experiments and simulations, where we also witnessed a strong fdependence (cf. Figures 4 and 5a). For film thicknesses beyond H = 25σ this simplified model starts to break down, since irrespective of the polymer composition, it predicts further C∥ to C⊥ transition sequences, although such a transition was found only for the f ≈ 0.21 case. This discrepancy is partly due to the fact that this model is not capable of capturing the regions where C∥ and C⊥ structures coexist. Furthermore, more complex structures than perfectly parallel and perpendicular cylinders can be stable in thick films, which we did not explicitly account for in our modeling. Nonetheless, this simplistic approach reproduces the alternation pattern observed in both experiments and simulations with surprising accuracy for thin films. Therefore, we can conclude that the alternation of the C∥ and C⊥ structures is mostly due to enthalpic contributions and that entropy only assumes a minor role.

form. We can estimate the diameter of the microdomain as d = 2.0RAe , where RAe denotes the end-to-end distance of the copolymer’s A block. The numerical value of D can then be derived from space-filling conditions, and we get d = 16 nm and D = 33 nm for the A−B 7/25 composition, which is in good agreement with the results from both experiments and simulations. Upon increasing the film thickness of a monolayer film to H′ > D, the additional A particles can be incorporated into the system in essentially three ways, by either (i) forming an additional layer of parallel cylinders, or (ii) growing the existing parallel microdomains along the z-axis, thereby preserving the in-plane spacing D, or (iii) changing the microdomain orientation to perpendicular. The latter two pathways are illustrated schematically in panels b and c of Figure 8, and we assume that any energy difference between the two configurations must be dominated by the surface terms (the remarkable agreement with experimental and simulation results offers a posteriori justification for this assumption). We can discern two competing contributions, where the first one stems from the incompatibility of A and B particles and the second one is due to the selectivity of the bead−surface interaction. For a fixed cross-sectional area, microdomains with an ellipsoidal cross section exhibit a significantly larger interface with the surrounding copolymer matrix than domains with a circular shape. If we take the C∥ configuration as our reference, we can estimate the resulting energy difference per particle as

Δu1 = f (1 − f )Δγ ΔS

4. CONCLUSIONS We performed dissipative particle dynamics simulations with a sophisticated coarse-graining protocol to study the phase behavior of PS−PHMA diblock copolymers confined in thin films. We found that the minority block (PS) of the copolymers formed cylindrical microdomains, in order to minimize the contact surface with the majority block (PHMA). As the film thickness was increased, orientation of these domains sequentially changed from perpendicular to parallel with respect to the substrate. We studied the influence of the volume fraction f and found that the thickness window in which the perpendicular configuration was preferable shrank as the PS content of the copolymers was increased. Furthermore, we investigated the effects of the copolymer size and found that the parallel configuration became increasingly favored for shorter chains. Lastly, the effects of polydispersity were studied, revealing no significant differences to the monodisperse case. Our simulation results agree in an almost perfect quantitative fashion with the findings from recent experiments,10,11 rendering the employed technique as a valuable asset in polymer research; following this route, the development of optimized and tailored materials can be accelerated significantly, since the simulation scheme presented here offers a very fast and cost-efficient alternative to traditional laboratory prototyping. In the second part of our study, we introduced a geometric model that accurately captured the sequence of C∥ and C⊥ structures observed in both experiments and simulations. The model calculations revealed that the surface terms are the prevalent driving forces behind the C∥ ↔ C⊥ transition and that entropy plays only a minor role. Because of its simple nature, this model can easily be applied to estimate the phase behavior of other copolymer melts confined in thin films. Future research is planned on the behavior of these thin films under steady shear, a method which is often employed to achieve long-range ordering of the microdomains. The particlebased design of the employed simulation method makes the implementation of shear flow with its ensuing hydrodynamic interactions rather straightforward compared to mean-field approaches such as SCFT.

(9)

where Δγ is the excess surface energy per unit length (surface tension), and ΔS is the discrepancy between the circumference of an ellipse and a circle enclosing the same area. In conjunction with eq 2, the surface tension can be estimated by σ Δγ = ⟨UAB⟩ − ⟨UAA ⟩ ≈

N̅ χN 2R e 3

(10)

If we now compare the C⊥ configuration with the undistorted C∥ case, it is readily visible that the two states differ only at the top and bottom surfaces of the film (highlighted by the small rectangular areas in Figure 8), where we find a sizable amount of A particles in the C⊥ structure. We can estimate the resulting energy difference as Δu 2 =

2σ f ΔεS H

(11)

where ΔεS denotes the energy penalty for placing an A particle at the surface. For the wall potential employed, and our set of parameters (see preceding section), we can estimate this quantity as ΔεS = 5kBT. Furthermore, we impose the constraint that every microdomain has to be surrounded by at least one layer of B particles, in order to avoid the interpenetration of two neighboring sites. We can then use eqs 9 and 11 to predict the orientation of the emerging cylinders, which are aligned parallel in the case of Δu1 < Δu2 and perpendicular otherwise. For all studied volume fractions (0.2 ≤ f ≤ 0.3), we identified the C⊥ structure for ultrathin films with H ≲ 7.5σ. At H ≈ 7.5σ, the orientation of the microdomains changed to parallel with respect to the substrate and remained so until the film thickness reached H ≈ 12.2σ, where the system went through a C∥ → C⊥ transition. The thickness range over which this C⊥ structure was preferred strongly depended on the exact value of f, and full G

dx.doi.org/10.1021/ma400867s | Macromolecules XXXX, XXX, XXX−XXX

Macromolecules



Article

(27) Chremos, A.; Margaritis, K.; Panagiotopoulos, A. Z. Soft Matter 2010, 6, 3588−3595. (28) Hoogerbrugge, P. J.; Koelman, J. M. V. A. Europhys. Lett. 1992, 19, 155. (29) Koelman, J. M. V. A.; Hoogerbrugge, P. J. Europhys. Lett. 1993, 21, 363. (30) Español, P.; Warren, P. B. Europhys. Lett. 1995, 30, 191. (31) Groot, R. D.; Warren, P. B. J. Chem. Phys. 1997, 107, 4423− 4435. (32) Rubinstein, M.; Colby, R. H. In Polymer Physics; Rubinstein, M., Colby, R. H., Eds.; Oxford University Press: New York, 2003. (33) Daoulas, K.; Müller, M.; de Pablo, J. J.; Nealey, P. F.; Smith, G. Soft Matter 2006, 2, 573−583. (34) Fetters, L. J.; Lohse, D. J.; Colby, R. H. In Physical Properties of Polymers Handbook, 2nd ed.; Mark, J. E., Ed.; Springer: New York, 2007; pp 447−454. (35) Mays, J. W.; Hadjichristidis, N. J. Macromol. Sci., Rev. Macromol. Chem. Phys. 1988, 28, 371−401. (36) Ferry, J. D.; Child, W. C., Jr. J. Colloid Sci. 1957, 12, 389−399. (37) Chinai, S. N. J. Polym. Sci., Part A: Polym. Chem. 1957, 25, 413. (38) Peters, B. L.; Ramírez-Hernández, A.; Pike, D. Q.; Müller, M.; de Pablo, J. J. Macromolecules 2012, 45, 8109−8116. (39) Ruzette, A. V. G.; Banerjee, P.; Mayes, A. M.; Pollard, M.; Russel, T. P.; Jerome, R.; Slawecki, T.; Hjelm, R.; Thiyagarajan, P. Macromolecules 1998, 31, 8509−8516. (40) Allen, M. P.; Tildesley, D. J. In Computer Simulation of Liquids; Oxford University Press: New York, 1989. (41) Frenkel, D.; Smit, B. In Understanding Molecular Simulation: From Algorithms to Applications; Academic Press: New York, 2001. (42) http://codeblue.umich.edu/hoomd-blue. (43) Anderson, J. A.; Lorenz, C. D.; Travesset, A. J. Comput. Phys. 2008, 227, 5342−5359. (44) Phillips, C. L.; Anderson, J. A.; Glotzer, S. C. J. Comput. Phys. 2011, 230, 7191−7201.

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected] (A.N.). Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS The authors thank Paul M. Chaikin and Raleigh L. Davis for many helpful discussions. This work has been supported by the Princeton Center for Complex Materials (PCCM), a U.S. National Science Foundation Materials Research Science and Engineering Center (Grant DMR-0819860). We also acknowledge the computational resources from the Texas Advanced Computing Center (TACC) at The University of Texas at Austin, which were provided through the Extreme Science and Engineering Discovery Environment (XSEDE), a National Science Foundation initiative (Grant OCI-1053575).



REFERENCES

(1) Colby, R. H. Curr. Opin. Colloid Interface Sci. 1996, 1, 454−465. (2) Bates, F. S.; Fredrickson, G. H. Phys. Today 1999, 52, 32−38. (3) Hamley, I. W. Macromol. Theory Simul. 2000, 9, 363−380. (4) Segalman, R. A. Mater. Sci. Eng., R 2005, 48, 191−226. (5) Darling, S. B. Prog. Polym. Sci. 2007, 32, 1152−1204. (6) Hamley, I. W. Prog. Polym. Sci. 2009, 34, 1161−1210. (7) Matsen, M. W.; Bates, F. S. Macromolecules 1996, 29, 1091−1098. (8) Knoll, A.; Horvat, A.; Lyakhova, K. S.; Krausch, G.; Sevink, G. J. A.; Zvelindovsky, A. V.; Magerle, R. Phys. Rev. Lett. 2002, 89, 035501. (9) Khanna, V.; Cochran, E. W.; Hexemer, A.; Stein, G. E.; Fredrickson, G. H.; Kramer, E. J.; Li, X.; Wang, J.; Hahn, S. F. Macromolecules 2006, 39, 9346−9356. (10) Davis, R. L.; Michal, B. T.; Adamson, D. H.; Chaikin, P. M.; Register, R. A. Proc. ACS Div. Polym. Mater. Sci. Eng.; 244th ACS National Meeting, Philadelphia, PA, August 19−23, 2012. (11) Davis, R. L.; Register, R. A. Personal communications, 2013. (12) Brown, G.; Chakrabarti, A. J. Chem. Phys. 1994, 101, 3310− 3317. (13) Brown, G.; Chakrabarti, A. J. Chem. Phys. 1995, 102, 1440− 1448. (14) Huinink, H. P.; Brokken-Zijp, J. C. M.; van Dijk, M. A.; Sevink, G. J. A. J. Chem. Phys. 2000, 112, 2452−2462. (15) Huinink, H. P.; van Dijk, M. A.; Brokken-Zijp, J. C. M.; Sevink, G. J. A. Macromolecules 2001, 34, 5325−5330. (16) Horvat, A.; Lyakhova, K. S.; Sevink, G. J. A.; Zvelindovsky, A. V.; Magerle, R. J. Chem. Phys. 2004, 120, 1117−1126. (17) Pelletier, V.; Asakawa, K.; Wu, M.; Adamson, D. H.; Register, R. A.; Chaikin, P. M. Appl. Phys. Lett. 2006, 88, 211114. (18) Hong, Y.-R.; Asakawa, K.; Adamson, D. H.; Chaikin, P. M.; Register, R. A. Opt. Lett. 2007, 32, 3125−3127. (19) Papalia, J. M.; Adamson, D. H.; Chaikin, P. M.; Register, R. A. J. Appl. Phys. 2010, 107, 084305. (20) Naito, K.; Hieda, H.; Sakurai, M.; Kamata, Y.; Asakawa, K. IEEE Trans. Magn. 2002, 38, 1949−1951. (21) Ruiz, R.; Kang, H. M.; Detcheverry, F. A.; Dobisz, E.; Kercher, D. S.; Albrecht, T. R.; de Pablo, J. J.; Nealey, P. F. Science 2008, 321, 936−939. (22) Pike, D. Q.; Detcheverry, F. A.; Müller, M.; de Pablo, J. J. J. Chem. Phys. 2009, 131, 084903. (23) Matsen, M. W. Macromolecules 2010, 43, 1671−1674. (24) Vega, D. A.; Harrison, C. K.; Angelescu, D. E.; Trawick, M. L.; Huse, D. A.; Chaikin, P. M.; Register, R. A. Phys. Rev. E 2005, 71, 061803. (25) Pinna, M.; Zvelindovsky, V. Eur. Phys. J. B 2012, 85, 210. (26) Horsch, M. A.; Zhang, Z.; Iacovella, C. R.; Glotzer, S. C. J. Chem. Phys. 2004, 121, 11455−11462. H

dx.doi.org/10.1021/ma400867s | Macromolecules XXXX, XXX, XXX−XXX