Optimized Phase Field Model for Diblock Copolymer Melts

Mar 15, 2019 - Expressions for these objects are given in the Supporting Information. In our previous procedure, we set the multiplicative constant by...
1 downloads 0 Views 1MB Size
Article Cite This: Macromolecules XXXX, XXX, XXX−XXX

pubs.acs.org/Macromolecules

Optimized Phase Field Model for Diblock Copolymer Melts Jimmy V. Liu,†,‡ Carlos J. García-Cervera,§,⊥ Kris T. Delaney,‡ and Glenn H. Fredrickson*,†,‡,∥ †

Department of Chemical Engineering, ‡Materials Research Laboratory, §Department of Mathematics, and ∥Materials Department, University of California, Santa Barbara, Santa Barbara, California 93106, United States ⊥ Visiting Professor at BCAMBasque Center for Applied Materials, Mazarredo 14, E48009 Bilbao, Basque Country, Spain

Macromolecules Downloaded from pubs.acs.org by KEAN UNIV on 03/28/19. For personal use only.

S Supporting Information *

ABSTRACT: Self-consistent field theory (SCFT) is a versatile framework for investigating self-assembly and thermodynamic properties in broad classes of polymer systems. Nonetheless, it is expensive to implement in higher space dimensions d as propagator fields of d + 1 dimensions must be resolved in a numerical SCFT simulation. In principle, SCFT can be exactly reduced to a phase field (or density functional) theory involving only d-dimensional species density fields, but the density functional of the reduced theory is not known analytically. Although asymptotic approximations to the functional are available in the literature, e.g., the Ohta−Kawasaki (OK) model for diblock copolymers, these functionals tend to produce inaccurate and unreliable predictions. Here, we present a phase field mapping technique adapted from force-matching concepts in the particle coarse-graining literature to best-fit coefficients in a trial density functional using data obtained from low-dimensional SCFT simulations. Such a mapping scheme results in an optimized phase field (OPF) model with improved accuracy over ad hoc choices of the fit parameters. We demonstrate the method by generating an OPF model for diblock copolymer melts, either symmetric or asymmetric, and we study its accuracy and transferability. We show that compared to the OK model, the OPF model makes significantly more accurate predictions of microphase domain periods. The OPF model also semi-quantitatively captures densities and energies, making it useful for computationally intensive applications, either on its own or in a multiscale approach alongside SCFT.



6(Ns M log M ), whereas a field update in a phase field simulation is a d-dimensional calculation with a complexity of only 6(M log M ). Here, Ns is the number of grid points used to resolve the contour variable and M is the number of grid points used to resolve the d spatial dimensions. Furthermore, in phase field models, stable and metastable states correspond to pure minima. This enables the use of efficient nonlinear conjugate gradient optimizers that require a convex free energy landscape (in SCFT, the presence of pressure-like fields leads to mixed saddle points and precludes conjugate gradient optimization). Finally, phase field models are purely realvalued. This eliminates the sign problem present in the original field theory, where complex-valued fields and energies lead to difficulties in sampling the probability distribution (beyond mean-field theory or SCFT). Because an accurate SCFT simulation typically requires Ns of about 100, we expect that a phase field simulation can be run to convergence faster than a comparable SCFT simulation by roughly two orders of magnitude.33 Historically, phase field models have been inaccurate compared to SCFT. We previously showed that well-known phase field models for the symmetric diblock copolymer melt fail to capture the correct shape or amplitude for the density profile and make poor predictions of length scales at moderate

INTRODUCTION Polymer field theory is a valuable tool for studying the behavior of dense polymer systems on length scales ranging from nanometers to microns. In particular, self-consistent field theory (SCFT) simulations have been used to predict phase behavior,1−15 to guide directed self-assembly schemes for patterning electronic devices,16−25 to design materials with desirable mechanical properties,26,27 and to study dynamical behavior of polymers in processes like solvent vapor annealing.28,29 SCFT has been successful at predicting density configurations, length scales, and energetics when thermal fluctuations can be neglected (that is, when the polymers are long and strongly overlapped). At an increased computational cost, we can also perform field-theoretic simulations using the complex Langevin method to sample these thermal fluctuations.30−32 Nonetheless, there are applications that are too computationally intensive even for SCFT, examples being extremely large cell simulations, large parameter sweeps, or non-equilibrium simulations where thermodynamic forces must be computed at each time step. Phase field (or density functional) models offer a more computationally efficient route than SCFT to tackle such problems. These models do not explicitly encode the contour variable s that indicates the position along a polymer chain in SCFT. This benefits numerical calculations of phase field models in several ways. For a simulation in d spatial dimensions conducted using an efficient pseudo-spectral scheme, a field update in SCFT is a d+1-dimensional calculation with a computational complexity of © XXXX American Chemical Society

Received: January 28, 2019 Revised: March 15, 2019

A

DOI: 10.1021/acs.macromol.9b00194 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules segregation strengths.34 Consequently, these models cannot reproduce defect states or defect formation energies known from SCFT. These failures are largely due to the approximations used to derive such models: the Landau−Brazovskii model35 and the Ohta−Kawasaki (OK) model36 both rely on the random phase approximation (RPA) first applied to the diblock copolymer by Leibler.37 The RPA is an asymptotic approximation that decreases in accuracy as χN increases above the weak segregation limit, χN → χNs, where χN is the segregation strength characterizing the repulsion between unlike segments and χNs is the mean-field spinodal value signaling the onset of microphase separation. In the intermediate and strong segregation regimes where χN > χNs, the assumed small parameter in the RPA (the amplitude of density variations) is large and the RPA is not strictly valid, but the OK model is often used as a crude substitute for SCFT. The OK model invokes two additional approximations. First, it replaces the quadratic coefficient in the RPA expansion with its limiting behavior for small and large wavevectors. This preserves the spinodal value χNs but leads to an incorrect critical wavevector k*. Second, the OK model makes a crude approximation of assuming the third- and fourth-order expansion coefficients are fully local. Later, Uneyama and Doi generalized the OK model to broader classes of copolymers and blends, but this generalized OK model, like the original OK model, makes significant errors relative to SCFT in predicting domain periods, density profiles, and energies.38 To address the inaccuracies of these asymptotic and phenomenological phase field models, we recently introduced a new approach to developing phase field models, which we refer to as phase field mapping.34,39 The method is based on the coarse-graining technique of force-matching for particle theories40 that was subsequently adapted to field theories.41 Like the coarse-graining method, phase field mapping parametrizes a cheaper model by minimizing the difference in thermodynamic forces between two modelsin this case, SCFT and a phase field model. In principle, the mapping can be exact: no accuracy is necessarily lost in switching to a phase field description of the theory. We previously outlined a procedure to parametrize an optimized phase field (OPF) model using phase field mapping and a small number of inexpensive SCFT calculations and reported the parametrization for a melt of symmetric diblock copolymers.39 Here, we describe in detail an improved and simplified procedure, and for the first time, provide results for the general case of asymmetric block fractions. We further investigate the performance of the model with respect to the accuracy of density profiles, domain periods, and defect energies relative to SCFT, and transferability to microphases other than lamellae.



R g = b N/6 , and all energies in units of kBT, the product of Boltzmann’s constant and the temperature. The canonical partition function for this system is often written in terms of two auxiliary chemical potential fields wA(r) and wB(r), eliminating explicit dependence on the species density fields (volume fractions) ϕA(r) and ϕB(r). To elucidate the connection between this model and the phase field model, however, we instead write an equivalent density-explicit form A=

∫ +ϕ ∫ +wA ∫ +wB e−H[w ,w ,ϕ] A

B

(1)

so that a phase field Hamiltonian (see eq 8 below) can emerge by the elimination of the chemical potential degrees of freedom in eq 1. We have abbreviated ϕ(r) = ϕA(r) and used incompressibility to eliminate ϕB(r) = 1 − ϕ(r). The action functional or effective Hamiltonian of this system is given by H[wA , wB , ϕ] = − χN C

∫ dr ϕ(r)2 − i ∫ dr[wA(r) − wB(r)]ϕ(r)

− V ln Q [iwA , iwB]

(2)

n V

Here, C = is the chain density given in units of chain partition function Q is given by

Q [iwA , iwB] =

1 V

R−3 g .

∫ dr q(r, 1; [iwA , iwB])

The single-

(3)

which in turn requires the solutions to the Fokker−Planck equation for a chain propagator q

∂q(r, s ; [iwA , iwB]) = ∇2 q(r, s ; [iwA , iwB]) ∂s − iw(r, s)q(r, s ; [iwA , iwB])

(4)

subject to the initial condition q(r, 0) = 1 and w(r, s) = wA(r) for contour positions s within the A block (0 ≤ s ≤ f) and w(r, s) = wB(r) otherwise ( f < s ≤ 1). These equations can be studied at a mean-field level, or we can use the complex Langevin technique to account for thermal fluctuations. In this study, we use mean-field theory to compute a representative configuration ϕ* for phase field mapping. The stable configuration ϕ* is a stationary saddle point in H that satisfies the self-consistent field equations: δH[wA , wB , ϕ] δϕ(r) δH[wA , wB , ϕ] δwA(r) δH[wA , wB , ϕ] δwB(r)

=0 (5)

ϕ = ϕ*, wA = wA* , wB = wB*

=0 (6)

ϕ = ϕ*, wA = wA* , wB = wB*

=0 (7)

ϕ = ϕ*, wA = wA* , wB = wB*

In computing ϕ*, we also apply the variable cell shape method to relax the dimensions of the simulation cell to a stress-free configuration reflecting the commensurate domain spacing L0. Pseudo-spectral methods are employed in all SCFT simulations, as described in ref 31. Phase Field Model. By formally evaluating the integrals over w configurations in eq 1, we can in principle express the partition function in terms of an exact reduced action He[ϕ]: 43

METHODS

Self-Consistent Field Theory. To obtain training sets for phase field mapping, we use a standard SCFT model for an incompressible AB diblock copolymer melt with continuous Gaussian chains and equal statistical segment lengths b for both A- and B-type monomers.42 In a system of volume V, each of the n chains is composed of N monomers, of which a fraction f belong to the A block. A Flory χ parameter characterizes the strength of the repulsive interactions between A and B monomers. Throughout this work, we express all lengths in units of the unperturbed radius of gyration

e−He[ϕ] = A=

∫ +wA ∫ +wB e−H[w ,w ,ϕ] A

∫ +ϕ e−H [ϕ] e

B

(8) (9)

The functional He rigorously defines a phase field model whose statistics are consistent with the field theory of eq 1, although no B

DOI: 10.1021/acs.macromol.9b00194 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules

As χN increases, other terms in eq 13 begin to play a role in setting L0, and the single-mode approximation breaks down. Equation 13 defines the functional form of the optimized phase field model, but we have not yet specified values for the coefficients {ci}, which should depend on f and χN. The parametrization is where the OPF model and the OK model differ. In the OK model, three of the coefficients are prescribed from the RPA by matching the secondorder vertex function in the limit of small and large wavevectors, yielding c6 and c5, respectively, then choosing c2 to reproduce the RPA spinodal χNs. In the present notation, these expressions are:

closed analytical expression for it exists. Conceptually, however, the form of He grants it several numerical advantages over the original complex-valued Hamiltonian H. Because He is purely real-valued, the Boltzmann weight e−He is positive definite, circumventing the sign problem associated with the original complex weight e−H. This would allow efficient Monte Carlo simulations of the theory in eq 9 and the use of advanced techniques such as flat histogram sampling not applicable to the complex-valued theory. At the mean-field level, because equilibrium solutions are pure minima in He, we can perform phase field calculations using efficient numerical methods such as the nonlinear conjugate gradient method. Finally, because the single-chain partition functions have been integrated out, the phase field theory does not require resolving the chain contour variable s, nor generating solutions to the modified diffusion equation in eq 4, reducing the theory from d + 1 internal dimensions to d spatial dimensions. In the following sections, it will be useful to express the thermodynamic force for the phase field model in terms of the original complex theory. By explicitly differentiating eq 8, we find that δHe[ϕ] δϕ(r)

= ϕ = ϕ0

δH[wA , wB , ϕ] δϕ(r)

ϕ0

c 2,OK = − χN + χNs − c5,OKk*2 − c6,OKk*− 2

(10)

∫ +wA ∫ +wB G[wA , wB , ϕ0] e−H[wA , wB , ϕ0] ∫ +wA ∫ +wB e−H[wA , wB , ϕ0] (11)

Like the original field theory, the model of eq 9 can be studied either at the mean-field level or accounting for fluctuations in ϕ. For this study, we focus exclusively on mean-field solutions. For such purposes, He can be viewed as a Helmholtz free energy functional whose minima represent mean-field solutions. The functional form of He is not analytically accessible. We instead approximate it using a linear expansion in a set of basis functionals {hi}: He[ϕ] ≈ Hr[ϕ] =

∑ cihi[ϕ] i

(12)

Here, the coefficients {ci} are left as parameters to be optimized by the mapping procedure. There are many possible choices for the basis set {hi}; choosing one that will yield an accurate, transferable phase field model is not straightforward. Previous studies using basis sets from the Landau−Brazovskii35 and Ohta−Kawasaki36 models have shown that the latter was better able to capture SCFT density profiles and length scales at intermediate segregation strengths.34,44 For the current extended study, we thus use a basis set drawn from the Ohta− Kawasaki model Hr[ϕ] = C

∫ dr{c2δϕ(r)2 + c3δϕ(r)3 + c4δϕ(r)4 + c5 |∇ϕ(r)|2 + c6

∫ dr′ G(r − r′)δϕ(r)δϕ(r′)}

(13)

where δϕ = ϕ − f is an order parameter with no k = 0 Fourier mode. The polynomial terms in the expansion produce a double well in the free energy that controls the amplitude of oscillations in the density field, in the manner of Landau theory. The square gradient term penalizes the formation of interfaces. In the final term, G is the Green’s function for Laplace’s equation: it satisfies ∇2G(r − r′) = −δ(r − r′). This term reflects the connectivity of the diblock copolymer and prevents macrophase separation. In the weak segregation limit, δϕ can be approximated by a single Fourier mode of wavenumber k* determined by the balance of the last two terms. 2π This sets a preferred length scale for domain sizes L0 = k * (valid when χN ≈ χNs) according to

k* = (c6/c5)1/4

c5,OK =

1 4f (1 − f )

(16)

c6,OK =

3 4f 2 (1 − f )2

(17)

The cubic and quartic polynomial coefficients c3 and c4 are not specified; the authors36 only note that in the strong segregation limit, the polynomial in δϕ must have two degenerate minima corresponding to the densities in the A-rich and B-rich phases. Other authors have assumed phenomenological values that produce the expected double well in Hr.45,46 To represent the OK model as fairly as possible in comparisons to our model, we parametrize c3 and c4 using phase field mapping, while retaining Ohta and Kawasaki’s original expressions for the remaining quadratic coefficients. Instead of this combination of asymptotic and phenomenological arguments, to obtain the five coefficients for our optimized phase field model, we apply the following phase field mapping method using a density configuration ϕ* obtained from SCFT. Critical to the utility of the method is that the SCFT reference configuration can be developed at a low computational cost (e.g., using a small 1D simulation), whereas the optimized phase field model can be applied to large-scale simulations in two or three dimensions. Phase Field Mapping. Phase field mapping is a new technique for producing accurate and transferable phase field models. It systematically parametrizes an OPF model using only a small number of inexpensive SCFT calculations combined with a set of basis functionals and associated coefficients. The method draws inspiration from force-matching, a technique for coarse-graining particle and field theories. Given an ideal (“complete”) set of basis functionals, we expect that the reduced action with optimized coefficients would agree perfectly with SCFT. Hence, the extent to which perfect correspondence with SCFT is not achieved is a measure of the quality of the basis set. We shall find that even a simple basis such as the OK form in eq 13 is capable of accurately predicting diblock copolymer structure and thermodynamics across a range of χN and f values. In this sense, phase field mapping has much in common with electronic density functional theory, where an exact functional is not known, but various approximate forms have been shown to provide good numerical agreement with more accurate theories. Given the Hamiltonian H, a representative SCFT configuration ϕ*, and a reduced action Hr with undetermined parameters {ci}, our phase field mapping procedure involves minimizing the following cost metric Ψ quantifying the disagreement in force and stress between the two models ÄÅ ÉÑ2 ÅÅ ÑÑ ÅÅ δHr[ϕ] ÑÑ δ [ ϕ ] H w , w , 1 A B ÑÑ Ψ({ci}) = − drÅÅÅ ÑÑ ÅÅ δϕ(r) δϕ V ( r ) Ñ ÅÅÇ ϕ* Ñ ÑÖ ϕ=ϕ* ÄÅ ÉÑ2 ÅÅÅ ∂Hr,V[ϕ] ÑÑÑ ∂HV[wA , wB , ϕ] ÑÑ + λÅÅÅÅ − ÑÑÑ ÅÅ ∂ϵ ∂ ϵ ϕ* Ñ ÅÇ ÑÖ ϕ=ϕ* (18) where we have defined the volume-intensive quantities Hr,V = Hr/V and HV = H/V. For the linear basis expansion form of Hr from eq 12, Ψ is a quadratic form in {ci} and minimization amounts to solving a

where the average over w field configurations at a constrained value of ϕ = ϕ0 is defined as ⟨G[wA , wB , ϕ]⟩ϕ0 =

(15)



(14) C

DOI: 10.1021/acs.macromol.9b00194 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules

Figure 1. From left to right: OPF model parameters −c2, −c3, c4, c5, c6, and the ratio c6/c5 as functions of χN for block fractions f from 0.5 (darkest) to 0.2 (lightest) obtained from phase field mapping (points), the regression (solid lines), and values from the Ohta−Kawasaki model (dashed lines). fix c2 = −1 before minimization, then determine the multiplicative factor by matching the elastic modulus from SCFT

system of linear equations in those coefficients. The linear system corresponding to eq 13 is provided in the Supporting Information. The first line of eq 18 represents force-matching. It penalizes differences between the thermodynamic forces in the two models. This penalty is zero in the ideal case where the forces exactly match, such as when Hr is exactly equal to He, as shown in eq 10, and positive otherwise. The second line of eq 18 matches stresses between the two models. The derivatives

∂Hr,V ∂ϵ

and

∂HV ∂ϵ

∂ 2Hr,V[ϕ] ∂ϵ

L − L0 , L0

ϕ=ϕ*

∂ 2HV[wA*, wB*, ϕ] ∂ϵ 2

ϕ=ϕ*

(19)

Expressions for these objects are given in the Supporting Information. In our previous procedure, we set the multiplicative constant by computing optimized values of Hr for several values of ϵ evenly spaced between −0.5 and 0.5.39 Equation 19 eliminates the need to compute additional density states by computing the necessary derivatives at ϕ*.

represent the elastic stress

required to produce an isotropic strain ϵ =

=

2

where L is a scale

factor for the system size associated with state ϕ* and L0 is a reference value associated with the commensurate (zero-stress) state from SCFT. The dimensionless parameter λ controls the balance between force- and stress-matching. As it is increased, the optimization coefficient values change and lead to improved predictions of the microdomain period until the values eventually saturate. We found that setting λ = 100 was sufficient to produce parameters within 1% of their saturated values. Finally, it should be noted that the constrained averages appearing in eq 18 are evaluated in the mean-field/SCFT approximation by simply replacing wA and wB in the force and stress terms by their saddle point values of w*A and w*B . If ϕ* is chosen according to the conditions described in the previous sections (if it is a stationary point in H that reflects the commensurate domain spacing L0), then the SCFT force and stress terms in the above expression are both zero. Minimization of eq 18 then amounts to minimizing the squares of the phase field force and stress in configuration ϕ*, with the relative weighting of the stress determined by λ. The linear system becomes homogeneous, so the mapping admits a family of solutions that differ only by a multiplicative factor on all of the {ci}. To break this degeneracy, we



RESULTS Optimized Parameters. Using a commensurate onedimensional (lamellar) density profile from SCFT as ϕ*, we performed phase field mapping to obtain model parameters c2, ..., c6 that represent the optimal values according to the metric of eq 18. Together with the OK basis set in eq 13, these parameters define the optimized phase field model for the diblock copolymer melt. We repeated the mapping process at various points in the diblock copolymer phase space (0.2 ≤ f ≤ 0.5 and 10.6 ≤ χN ≤ 35) for which the lamellar phase is at least metastable. Throughout this range, all of the parameters appear to be smooth functions of both variables. Figure 1 compares the mapped parameters and fitted curves to the Ohta−Kawasaki model parameters. In the weak segregation regime (near the spinodal), the mapped parameters recover similar values to the RPA and the OK model. As χN increases away from the spinodal, the optimized parameters rapidly D

DOI: 10.1021/acs.macromol.9b00194 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules diverge from the RPA values. In particular, we note that the ratio c6/c5, which plays a dominant role in setting the domain spacing, adjusts to reflect the information obtained from SCFT. For convenience and future reference, we fit the OPF model parameters to simple polynomial functions of f and χN using least-squares regression to the intermediate segregation data and taking advantage of the symmetry about f = 0.5: 2

Table 2. Diblock Copolymer OPF Model Parameters b(2) jk

b(3) jk

2

j 2k c 2 = − ∑ ∑ b(2) jk x g

b(4) jk

(20)

j=0 k=0 2

2

j 2k + 1 c3 = − ∑ ∑ b(3) jk x g

b(5) jk

(21)

j=0 k=0 2

c4 =

2

b(6) jk

j 2k ∑ ∑ b(4) jk x g

(22)

j=0 k=0

j j j j j j j j j j j j j j j

= = = = = = = = = = = = = = =

0 1 2 0 1 2 0 1 2 0 1 2 0 1 2

k=0

k=1

k=2

5.920 2.025 0.005522 9.741 9.224 0.06433 9.686 36 0.02068 0.7853 0.2356 0.1185 0.5 0.2357 0.0006666

−14.31 −4.285 −0.1915 −39.46 14.69 −1.281 53.00 − −0.2385 −5.654 −1.170 −0.7423 − −1.956 −0.02858

398.5 39.47 −1.003 −999.0 −510.3 15.70 −1775.0 − −0.4559 −16.22 3.659 5.481 − 7.147 0.05316

2

c5 =

2k (5) ∑k = 0 (b0(5) k + b1k x)g 2

2k 1 + ∑k = 0 b2(5) k xg 2

(23)

ance of the OPF model with regards to accuracy and transferability by comparing its predictions to those of SCFT in various contexts. A key property for a phase field model to reproduce is the commensurate domain spacing L0. This length scale describes the preferred periodicity for lamellar ordering at a given value of f and χN. From previous studies using SCFT, L0/Rg is known to scale as Nα−0.5, where α is approximately 1 just above the weak segregation limit and gradually decreases toward 2 as 3 χN increases.47 Note that this scaling is given for L0 in units of R g = b N/6 , as throughout the rest of this work; if lengths are instead written in dimensional units, the relation becomes L0 ∼ Nα. The scaling exponent α has historically been incorrect for phase field models neglecting vertex functions beyond quadratic order.48 Phase field mapping can correct this issue systematically even with the simple OK basis set. We computed equilibrium solutions that minimize HV and Hr,V with respect to the simulation cell size to obtain predictions of L0 for SCFT, for the OPF model, and for the OK model throughout the range of f and χN described above. We then computed the signed fractional error η made by each phase field model relative to SCFT, given by L0,PF − L0 η= L0 (25)

2

j 2k c6 = −c 2/ ∑ ∑ b(6) jk x g

(24)

j=0 k=0

Here, we have defined the block fraction asymmetry g = 0.5 − f and the distance from the spinodal x = χN − χNs, where χNs is computed numerically from the RPA37 for each value of g (some values are shown in Table 1). Terms with exponents j = Table 1. RPA Spinodal Conditions for the Diblock Melt |g|

χNs

0.00 0.05 0.10 0.15 0.20 0.25 0.30 0.35

10.495 10.698 11.344 12.562 14.635 18.172 24.613 38.038

0 or k = 0 are treated as 1 (constants). The parameters b(i) jk appearing in eqs 20−24 are listed in Table 2. These expressions are valid in the region 0.2 ≤ f ≤ 0.8 and 10.6 ≤ χN ≤ 35 (provided χN > χNs). Values indicated by a dash were set to zero manually, because they did not significantly impact the quality of the fit. As Figure 1 shows, the regression is in good agreement with the original data throughout, particularly in the intermediate segregation regime, where the error made by the regression is within 2%. Although the individual c5 and c6 coefficients show some discrepancy in the weak segregation regime, the ratio c6/ c5 is still accurate to within 5%. The regression thus retains the advantages of the OPF model over the OK model, as we discuss next. Predicted Domain Spacing. The phase field mapping procedure described above provides a method for generating an optimized phase field model, but the quality of the mapped model is largely dependent on the choice of basis functionals and the training configuration ϕ*. We evaluated the perform-

Figure 2 shows that the OK model fails to capture the scaling of L0 as χN increases. We attribute this to several approximations in the model. The asymptotic expansion of the RPA is not appropriate at intermediate segregation, where density perturbations are no longer small, and neglecting the nonlocal character of the third- and fourth-order RPA vertex functions will also distort predictions of L0. Furthermore, approximating the second-order vertex function by its limiting behavior for small and large wavevectors causes k* to be incorrect (according to eq 14), so the OK model overestimates L0 even in the weak segregation limit. On the other hand, by modifying the ratio c6/c5 as χN increases, the OPF model provides qualitative and sometimes quantitative agreement with SCFT, especially for block fractions near symmetry where the lamellar phase is globally stable. It also correctly captures E

DOI: 10.1021/acs.macromol.9b00194 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules

Figure 2. (a) Predicted domain spacing L0 in units of Rg for SCFT at f = 0.5 (solid line) and at f = 0.3 (dashed line), for the OPF model at f = 0.5 (closed circles) and at f = 0.3 (open circles), and for the OK model at f = 0.5 (closed triangles) and at f = 0.3 (open triangles). (b, c) Error η in the lamellar domain spacing predicted by two phase field models, the OPF model (center) and the Ohta−Kawasaki model (right), as functions of χN for various values of f.

the scaling of L0 with increasing N: the error quickly settles to a constant value as χN increases above the spinodal. Density Fields. Next, we compare the shape of the density fields predicted by SCFT and by the OK and OPF models. Figure 3 shows a representative example of equilibrium density

partition these domains in asymmetric cases. We expect that resolving this issue would require expanding the basis set to include additional terms (other than the c3 term) that vanish for f = 0.5. Defect Annihilation and Kinetic Barriers. For some applications, a pure ordered state serves only as a reference, whereas the main focus of study is on defect states and their associated energies. In a previous study, we showed that the OPF model can consistently predict the same types of lamellar defect states as SCFT.39 Müller and Orozco Rey also found qualitative success in the study of lamellar defects and kinetic effects using the OK model.44 Both studies suggest that caution should be exercised in drawing conclusions with the OK model. For example, we also showed that at f = 0.5 and χN = 25, the OPF model predicts defect formation energies that are within 20% of the SCFT values, whereas the OK model predicts energies with the wrong sign.39 Here, we extend our analysis of the OPF model to examine kinetic pathways for defect annihilation. A number of studies have investigated kinetic pathways in SCFT by applying the string method.21,23,49−51 To do this, we use the mean-field approximation to construct a free energy landscape F[ϕ]= H[w A*, w B*, ϕ] that allows arbitrary configurations for the density field ϕ while relaxing the remaining fields to values stationary in H at the prescribed ϕ. For the OPF model, the analogous free energy functional is simply F[ϕ] = Hr[ϕ]. For such functionals, we use the string method52 to compute a minimum energy path (MEP) connecting two states, such as a metastable defect ϕd and the stable perfect lamellar configuration ϕp. This MEP consists of a continuous curve of configurations ϕ(α) beginning at ϕ(0) = ϕd and ending at ϕ(1) = ϕp, such that the entire curve is always parallel to the force ∇F. That is

Figure 3. Representative stress-free equilibrium density fields for SCFT (solid), the OPF model (dashed blue), and the Ohta− Kawasaki model (dotted orange) for the lamellar phase at f = 0.4 and χN = 35.

profiles for the lamellar phase at f = 0.4 and χN = 35; density profiles for other combinations of f and χN are provided in the Supporting Information. Because we know that the OK model makes poor predictions of L0 here, we allowed each model to use a periodic cell with its own commensurate domain spacing. Thus, each curve corresponds to an equilibrium state free of internally generated stresses, and each curve terminates at the model’s respective value for L0 or L0,PF. We see in Figure 3 that, overall, the OPF model follows the SCFT profile more closely than the OK model does, a trend that is consistent for all values of f and χN tested. On symmetry (when f = 0.5), the OPF result is almost indistinguishable from SCFT for all values in 10.4 ≤ χN ≤ 35. As f deviates from 0.5, though, the OPF profile begins to deviate from the SCFT result. Note that both the OK and OPF model profiles tend to overestimate the width of the minority A domain and underestimate the width of the majority B domain. This indicates that the basis set from eq 13 struggles to correctly

(∇F )⊥ (ϕ(α)) = 0

(26)

where the arc length α is a distance defined by the L norm of ϕ, normalized to 1, and (∇F)⊥ represents the component of the force normal to ϕ(α). Along the MEP, the free energy is maximized at an intermediate barrier state ϕ† with energy F[ϕ†]. Assuming isotropic mobilities so that the MEP is the highest flux kinetic pathway, we can then estimate the barrier crossing frequency as proportional to e−Eb/kBT, where the barrier height is Eb = F[ϕ†] − F[ϕd]. 2

F

DOI: 10.1021/acs.macromol.9b00194 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules

Figure 4. Top: minimum energy pathway for the annihilation of a disclination defect in a lamella-forming diblock copolymer melt with χN = 25 and f = 0.5 as found in SCFT (diamonds) and the OPF model (circles). Bottom, left to right: snapshots of density profiles along the MEP from the OPF model (visually indistinguishable from SCFT) corresponding to the local extrema of F.

slightly with the ends of the cell, as evidenced by a slight tilting of the lamellae most distant from the asymmetric defect. Although the qualitative shape of the MEP is captured by the OPF model, the OPF model underestimates the formation energy for the defect F(0) by about 20%, just as in our previous work. We also find that the OPF model significantly overestimates the barrier heights for both the disclination and dislocation. We further investigated these trends by repeating the string calculation for different values of χN for which the disclination remained metastable (down to χN = 20). Figure 5

To test the accuracy and transferability of the OPF model for such an application, we replicated an SCFT calculation from Takahashi et al. for the annihilation of a disclination defect in a lamellar diblock copolymer melt.21 We chose to study the most relevant melting pathway (the one predicted to have the lowest energy barrier), an asymmetric pathway that passes through a second metastable state, a dislocation defect. To avoid the use of confining walls, as were used in the Takahashi study, we imposed periodic boundary conditions as usual but surrounded the defect region with additional nondefective lamellae. These serve to stabilize the system against finite size effects by isolating the defect from the edges of the cell. Our version of this system contains seven commensurate lamellar periods at f = 0.5 and has dimensions of 6 L0 × 7 L0. At χN = 25, this corresponds to 25.7 Rg × 30.0 Rg. As done by Takahashi et al., to obtain extensive energies, we assume a film thickness of 4 Rg and a chain density of 3.1 R−3 g , corresponding to a polystyrene-b-poly(methyl methacrylate) diblock copolymer with a molecular weight of 72 kg/mol, a monomer density of 1 g/cm3, and a radius of gyration Rg = 7.2 nm. Figure 4 compares the minimum energy pathways predicted by SCFT and the OPF model for asymmetric melting of the disclination defect at χN = 25. Qualitatively, the two models are in strong agreement for the entire MEP. Both models predict MEPs containing two barriers (local maxima), and snapshots of the density fields along the string are indistinguishable between the two models. At the first transition state, the polymer chains rearrange to form a narrow bridge across the break in the lower lamella (for SCFT: Eb = 2.7 kBT; for OPF: Eb = 10 kBT). Next, we find a local minimum representing a dislocation defect. A second bridge then forms to repair the remaining broken lamella (for SCFT: Eb = 5.0 kBT; for OPF: Eb = 20 kBT). This second transition state gradually heals to the perfect lamellar state at α = 1. We see a slight finite size error in the intermediate states: in spite of the extra lamellae surrounding the defect, the defect still interacts

Figure 5. Energy barriers Eb for the melting of disclination (open shapes) and dislocation (closed shapes) defects versus χN as predicted by SCFT (diamonds) and the OPF model (circles).

shows the resulting defect formation energies and barrier heights as a function of χN. We observe that the OPF model correctly predicts increases in Eb with χN, but that it consistently overestimates its magnitude. It may also allow both defects to remain metastable for smaller values of χN where SCFT will indicate that they are unstable. Although the OPF model is seen to fail quantitatively in this test, it is important to note that the model was trained by a 1D SCFT simulation in a perfect lamellar phase, yet the string calculation G

DOI: 10.1021/acs.macromol.9b00194 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules probes its ability to reproduce complex free energy pathways in rugged high-dimensional landscapes. That it works even qualitatively is frankly quite remarkable. Transferability Across Different Morphologies. We have demonstrated that the OPF model produced by mapping from the commensurate lamellar choice for ϕ* can be used for qualitative investigations of defects and transition states. Another important question is whether the lamellar mapping is still valid if the resulting OPF model is instead used to study a different morphology. To address this issue, we performed calculations for the hexagonally packed cylinder phase (HEX) using both SCFT and the OPF model parametrized from the lamellar phase (with coefficients summarized in Table 2). We again minimized HV and Hr,V with respect to the simulation cell size and obtained predictions for the commensurate domain spacing for the HEX phase L0,HEX, defined as the side length of a rhombic unit cell. The fractional error ηHEX made by the OPF model (mapped from the lamellar phase) relative to SCFT (for the HEX phase) is shown in Figure 6a. As in the lamellar case, the OPF model qualitatively captures the magnitude of L0,HEX both in the weak segregation and intermediate segregation regimes. It performs best for block fractions near f = 0.3, where L0,HEX is accurate to within 4%. This is especially relevant since HEX is the globally stable phase at these conditions.

An alternate approach to studying the HEX phase is to redo the mapping from a cylindrical SCFT reference state ϕ*HEX. For points in the phase space where the HEX phase was at least metastable, we used phase field mapping to obtain another set of values for the model parameters that we will refer to as OPF-HEX. Note that the OPF-HEX mapping extends to f = 0.15, since the HEX phase remains metastable for stronger asymmetries than the lamellar phase, but it could not be parametrized for the symmetric case f = 0.5. The coefficient values differ somewhat from the lamellar mapping but still follow similar trends in χN and f; the data are shown in the Supporting Information. We computed values of L0,HEX predicted by the HEX mapping. The resulting error (again with respect to SCFT) is shown in Figure 6b. As expected, because the OPF-HEX model was trained in the HEX geometry, it tends to make slightly better predictions of L0,HEX than the lamellar mapping. Still, the original mapping obtained simply from a onedimensional SCFT calculation performs comparably well for block fractions near f = 0.30.



DISCUSSION Our results show that phase field mapping can systematically parametrize a useful model using an inexpensive onedimensional SCFT calculation. Even from this limited information, the mapping procedure will automatically adjust model parameters to capture the appropriate length scale for ordering. Furthermore, the mapping is highly transferable: the OPF model that we obtained by training from a nondefective lamellar state ϕ* can be used to study defect states, transition states, and the cylindrical morphology with qualitative and sometimes quantitative accuracy. Using tabulated data like that provided in eqs 20−24 and Table 2, the OPF model parameters c2−c6 can be used as direct replacements in existing numerical codes using the Ohta−Kawasaki model to improve the simulation accuracy at no additional computational cost. As for computational efficiency relative to SCFT, we found that OPF model calculations optimized to use a nonlinear conjugate gradient solver could be completed roughly two orders of magnitude faster than the corresponding SCFT calculations on the same hardware33 (our SCFT implementation employs the semi-implicit Seidel scheme detailed in ref 53. All simulations were performed using pseudo-spectral methods that efficiently handle operators local in real space or reciprocal space). For applications in which calculation time is a limiting factor, then, the OPF model is strongly advantaged relative to SCFT. The OPF model has another computational advantage: it is less memory intensive. SCFT calculations require the storage of 6(Ns) fields corresponding to the propagator q(r, s) at different contour points s, whereas OPF calculations are independent of Ns. In our implementation, we found that a bulk SCFT simulation with Ns = 100 contour points must allocate 109 fields, whereas a corresponding OPF simulation requires only 13 fields. On any given hardware with a fixed amount of memory, then, a phase field simulation can access system volumes 9 times larger than SCFT simulations. This makes the OPF model an appealing choice for repetitive simulations of large 3D structures, including studies of systems with large unit cells such as the Frank−Kasper sigma phase. We caution the reader, however, that the approach has yet to be tested on conformationally asymmetric block copolymers. It

Figure 6. Error ηHEX in the cylindrical domain spacing L0,HEX as a function of χN for various values of f, relative to SCFT, predicted by the OPF model (a) when mapped from the lamellar phase and (b) when mapped from the HEX phase. H

DOI: 10.1021/acs.macromol.9b00194 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules

models for polymer systems other than the diblock copolymer melt, such as diblock/homopolymer blends or complex architecture block and graft polymers, for which the mapping method is yet untested. In these cases, we must track more than one order parameter. Identifying an appropriate basis set, perhaps using an asymptotic approach built on the foundations of work by Uneyama and Doi,38 will be crucial.

remains to be seen whether OPF models can discriminate between complex phases possessing small differences in free energy. In terms of accuracy, we found that the OPF model can qualitatively reproduce density fields computed from SCFT and typically captures the relevant length scales to within 5%. Energy calculations are also qualitatively accurate. Unlike the OK model, the OPF model can consistently predict the correct sign for defect formation energies, but it also tends to underestimate their magnitudes. The same OPF model also overestimates the sizes of the energetic barriers for escaping metastable defects. If quantitative accuracy in such applications is a concern, we suggest a multiscale approach that first computes an approximate MEP pathway for the density using the OPF model, then uses it as an initial condition for a short SCFT relaxation step to reduce the overall cost of a full SCFT calculation. The mapping that we obtained is not transferable across temperatures: all of the OPF model parameters vary strongly with χ. In practice, this is straightforward to address. Because the parameters in Figure 1 are all smooth functions of χN, it is inexpensive to perform phase field mapping at a few representative values of χN and then use interpolation to accurately model the polymer system throughout the range of interest. In principle, however, the χ-dependence of the exact reduced action He in eq 9 should be fully contained in the quadratic term. This can be seen by tracing the first term in eq 2 into the integrand of eq 8. The χ term is independent of both wA and wB and can be moved outside the integrals, so it should appear in the same form in He: c 2 = −χN + c 2S



ASSOCIATED CONTENT

S Supporting Information *

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.macromol.9b00194. Force and stress calculations for the cost function Ψ in eq 18; statistics for the parameter regression; HEX phase coefficients; and representative stress-free equilibrium density fields for SCFT, the OPF model, and the OK model (PDF)



AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. ORCID

Jimmy V. Liu: 0000-0002-3880-8772 Kris T. Delaney: 0000-0003-0356-1391 Glenn H. Fredrickson: 0000-0002-6716-9017 Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS We are grateful to Samsung Electronics and the National Science Foundation CMMT Program under grant DMR1822215 for partial support of this research. We also acknowledge support from the Center for Scientific Computing from the CNSI, MRL: an NSF MRSEC (DMR-1720256).

(27)

and the remaining coefficients c2S, c3, ..., c6 should depend on f but not on χN. One can then ask why we observe χ dependence in model parameters other than c2. We attribute this to a deficiency in the basis set: eq 13 is an incomplete approximation for the exact action functional He. The mapping procedure can compensate for physics absent from the basis set by redistributing the effects of missing terms through the other coefficients, but the redistribution is not consistent across different values of χN. In this sense, the extent to which c2S, c3, ..., c6 are independent of χ is a measure of the quality of the basis set. Phase field mapping will unequivocally produce the best possible model for a given basis set according to the cost metric in eq 18, but the accuracy and transferability of the model still rest on the selection of this basis set. Because we do not have a complete set of basis functionals that spans the space of He, in this study, the basis set was chosen intuitively from past experience, drawing on work done using asymptotic methods and Landau theory. An improved basis set would likely introduce nonlocal terms beyond quadratic order in δϕ to correct density profiles at asymmetric block fractions. Alternately, we envision a combinatorial approach, where phase field mapping is applied to basis sets drawn from a large library of basis functionals. The library can be tested for quality by metrics such as the value of the cost function Ψ, the domain period error of eq 25, or by the transferability of the basis coefficients across χ. Directions for future work include studies of the transferability of a bulk parametrization to confined systems, where density gradients may be large and chain stretching effects significant. Future work will also include developing optimized



REFERENCES

(1) Matsen, M. W.; Schick, M. Stable and unstable phases of a diblock copolymer melt. Phys. Rev. Lett. 1994, 72, 2660−2663. (2) Matsen, M. W. Phase behavior of block copolymer/ homopolymer blends. Macromolecules 1995, 28, 5765−5773. (3) Tang, P.; Qiu, F.; Zhang, H.; Yang, Y. Morphology and phase diagram of complex block copolymers: ABC linear triblock copolymers. Phys. Rev. E 2004, 69, No. 031803. (4) Tyler, C. A.; Morse, D. C. Orthorhombic Fddd network in triblock and diblock copolymer melts. Phys. Rev. Lett. 2005, No. 208302. (5) Suo, T.; Yan, D.; Yang, S.; Shi, A.-C. A theoretical study of phase behaviors for diblock copolymers in selective solvents. Macromolecules 2009, 42, 6791−6798. (6) Qin, J.; Bates, F. S.; Morse, D. C. Phase behavior of nonfrustrated ABC triblock copolymers: weak and intermediate segregation. Macromolecules 2010, 43, 5128−5136. (7) Matsen, M. W. Effect of architecture on the phase behavior of AB-type block copolymer melts. Macromolecules 2012, 2161. (8) Shi, A.-C.; Li, B. Self-assembly of diblock copolymers under confinement. Soft Matter 2013, 9, 1398−1413. (9) Xie, N.; Li, W.; Qiu, F.; Shi, A.-C. σ phase formed in conformationally asymmetric AB-type block copolymers. ACS Macro Lett. 2014, 3, 906−910. (10) Mester, Z.; Lynd, N. A.; Delaney, K. T.; Fredrickson, G. H. Phase coexistence calculations of reversibly bonded block copolymers: A unit cell Gibbs ensemble approach. Macromolecules 2014, 47, 1865−1874.

I

DOI: 10.1021/acs.macromol.9b00194 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules (11) Li, W.; Delaney, K. T.; Fredrickson, G. H. Fddd network phase in ABA triblock copolymer melts. J. Polym. Sci., Part B: Polym. Phys. 2016, 1112. (12) Tsai, C. L.; Delaney, K. T.; Fredrickson, G. H. Genetic Algorithm for Discovery of Globally Stable Phases in Block Copolymers. Macromolecules 2016, 6558. (13) Khadilkar, M.; Paradiso, S.; Delaney, K.; Fredrickson, G. Inverse Design of Bulk Morphologies in Multiblock Polymers Using Particle Swarm Optimization. Macromolecules 2017, 6702. (14) Chang, A. B.; Bates, C. M.; Lee, B.; Garland, C. M.; Jones, S. C.; Spencer, R. K. W.; Matsen, M. W.; Grubbs, R. H. Manipulating the ABCs of self-assembly via low-χ block polymer design. Proc. Natl. Acad. Sci. U.S.A. 2017, 114, 6462−6467. (15) Kim, K.; Arora, A.; Lewis, R. M.; Liu, M.; Li, W.; Shi, A.-C.; Dorfman, K. D.; Bates, F. S. Origins of low-symmetry phases in asymmetric diblock copolymer melts. Proc. Natl. Acad. Sci. U.S.A. 2018, 115, 847−854. (16) Bosse, A. W.; García-Cervera, C. J.; Fredrickson, G. H. Microdomain ordering in laterally confined block copolymer thin films. Macromolecules 2007, 40, 9570−9581. (17) Tang, C.; Lennon, E. M.; Fredrickson, G. H.; Kramer, E. J.; Hawker, C. J. Evolution of block copolymer lithography to highly ordered square arrays. Science 2008, 429. (18) Detcheverry, F. A.; Nealey, P. F.; de Pablo, J. J. Directed assembly of a cylinder-forming diblock copolymer: Topographic and chemical patterns. Macromolecules 2010, 43, 6495−6504. (19) Detcheverry, F. A.; Liu, G.; Nealey, P. F.; de Pablo, J. J. Interpolation in the directed assembly of block copolymers on nanopatterned substrates: simulation and experiments. Macromolecules 2010, 43, 3446−3454. (20) Nagpal, U.; Müller, M.; Nealey, P. F.; De Pablo, J. J. Free energy of defects in ordered assemblies of block copolymer domains. ACS Macro Lett. 2012, 1, 418−422. (21) Takahashi, H.; Laachi, N.; Delaney, K. T.; Hur, S. M.; Weinheimer, C. J.; Shykind, D.; Fredrickson, G. H. Defectivity in laterally confined lamella-forming diblock copolymers: Thermodynamic and kinetic aspects. Macromolecules 2012, 45, 6253−6265. (22) Liu, C.-C.; Ramírez-Hernández, A.; Han, E.; Craig, G. S. W.; Tada, Y.; Yoshida, H.; Kang, H.; Ji, S.; Gopalan, P.; de Pablo, J. J.; Nealey, P. F. Chemical Patterns for Directed Self-Assembly of Lamellae-Forming Block Copolymers with Density Multiplication of Features. Macromolecules 2013, 46, 1415−1424. (23) Laachi, N.; Iwama, T.; Delaney, K. T.; Shykind, D.; Weinheimer, C. J.; Fredrickson, G. H. Directed self-assembly of linear arrays of block copolymer cylinders. J. Polym. Sci., Part B: Polym. Phys. 2015, 53, 317−326. (24) Ouaknin, G.; Laachi, N.; Bochkov, D.; Delaney, K.; Fredrickson, G. H.; Gibou, F. Functional level-set derivative for a polymer self consistent field theory Hamiltonian. J. Comput. Phys. 2017, 207. (25) Carpenter, C. L.; Nicaise, S.; Theofanis, P. L.; Shykind, D.; Berggren, K. K.; Delaney, K. T.; Fredrickson, G. H. Orientational Preference in Multilayer Block Copolymer Nanomeshes with Respect to Layer-to-Layer Commensurability. Macromolecules 2017, 8258. (26) Drolet, F.; Fredrickson, G. H. Optimizing chain bridging in complex block copolymers. Macromolecules 2001, 5317. (27) Lynd, N. A.; Oyerokun, F. T.; O’Donoghue, D. L.; Handlin, D. L.; Fredrickson, G. H. Design of soft and strong thermoplastic elastomers based on nonlinear block copolymer architectures using self-consistent-field theory. Macromolecules 2010, 3479. (28) Fraaije, J. G. E. M. Dynamic density functional theory for microphase separation kinetics of block copolymer melts. J. Chem. Phys. 1993, 9202. (29) Paradiso, S. P.; Delaney, K. T.; García-Cervera, C. J.; Ceniceros, H. D.; Fredrickson, G. H. Block copolymer self assembly during rapid solvent evaporation: Insights into cylinder growth and stability. ACS Macro Lett. 2014, 16.

(30) Fredrickson, G. H.; Ganesan, V.; Drolet, F. Field-theoretic computer simulation methods for polymers and complex fluids. Macromolecules 2002, 35, 16−39. (31) Fredrickson, G. H. The Equilibrium Theory of Inhomogeneous Polymers; Oxford University Press, 2006; pp 1−452. (32) Delaney, K. T.; Fredrickson, G. H. Recent Developments in Fully Fluctuating Field-Theoretic Simulations of Polymer Melts and Solutions. J. Phys. Chem. B 2016, 120, 7615−7634. (33) Liu, J.; Delaney, K. T.; Fredrickson, G. H. Optimized Phase Field Models in Confinement: Fast and Accurate Simulations of Directed SelfAssembly, Proceedings of the SPIEAdvances in Patterning Materials and Processes XXXIV; SPIE, 2017; Vol. 10146, 101460Z. (34) Liu, J.; Laachi, N.; Delaney, K. T.; Fredrickson, G. H. In Advantages and Limitations of Density Functional Theory in Block Copolymer Directed Self-Assembly, Proceedings of the SPIE Alternative Lithographic Technologies VII; SPIE, 2015; Vol. 9423, 94231I. (35) Fredrickson, G. H.; Helfand, E. Fluctuation effects in the theory of microphase separation in block copolymers. J. Chem. Phys. 1987, 87, 697. (36) Ohta, T.; Kawasaki, K. Equilibrium morphology of block copolymer melts. Macromolecules 1986, 19, 2621−2632. (37) Leibler, L. Theory of microphase separation in block copolymers. Macromolecules 1980, 13, 1602−1617. (38) Uneyama, T.; Doi, M. Density Functional Theory for Block Copolymer Melts and Blends. Macromolecules 2005, 38, 196−205. (39) Liu, J.; Delaney, K. T.; Fredrickson, G. H. In Phase Field Mapping for Accurate, Ultrafast Simulations of Directed Self-Assembly, Proceedings of the SPIEAlternative Lithographic Technologies VIII; SPIE, 2016; Vol. 9779, p 977920. (40) Izvekov, S.; Voth, G. A. A multiscale coarse-graining method for biomolecular systems. J. Phys. Chem. B 2005, 109, 2469−2473. (41) Villet, M. C.; Fredrickson, G. H. Numerical coarse-graining of fluid field theories. J. Chem. Phys. 2010, 132, No. 034109. (42) Matsen, M. W. The standard Gaussian model for block copolymer melts. J. Phys.: Condens. Matter 2002, 14, R21. (43) Barrat, J. L.; Fredrickson, G. H.; Sides, S. W. Introducing variable cell shape methods in field theory simulations of polymers. J. Phys. Chem. B 2005, 109, 6694−6700. (44) Müller, M.; Orozco Rey, J. C. Continuum models for directed self-assembly. Mol. Syst. Des. Eng. 2018, 3, 295−313. (45) Yoshimoto, K.; Taniguchi, T. Large-Scale Simulations of Directed Self-Assembly with Simplified Model. J. Photopolym. Sci. Technol. 2013, 26, 809−816. (46) Muramatsu, M.; Nakano, T.; Tomita, T.; Yamamoto, K.; Matsuzaki, K.; Kitano, T. In Simulation Analysis of Directed SelfAssembly for Hole Multiplication in Guide Pattern, Proceedings of the SPIEAlternative Lithographic Technologies VI; SPIE, 2014; Vol. 9049, 904921. (47) Matsen, M. W.; Bates, F. S. Unifying weak- and strongsegregation block copolymer theories. Macromolecules 1996, 29, 1091−1098. (48) McMullen, W. E. A fourth-order, density-functional, randomphase approximation study of monomer segregation in phaseseparated, lamellar, diblock-copolymer melts. Macromolecules 1993, 26, 1027−1036. (49) Cheng, X.; Lin, L.; E, W.; Zhang, P.; Shi, A.-C. Nucleation of ordered phases in block copolymers. Phys. Rev. Lett. 2010, 104, No. 148301. (50) Kim, B.; Laachi, N.; Delaney, K. T.; Carilli, M.; Kramer, E. J.; Fredrickson, G. H. Thermodynamic and kinetic aspects of defectivity in directed self-assembly of cylinder-forming diblock copolymers in laterally confining thin channels. J. Appl. Polym. Sci. 2014, No. 40790. (51) Iwama, T.; Laachi, N.; Kim, B.; Carilli, M.; Delaney, K. T.; Fredrickson, G. H. Computational study of directed self-assembly in neutral prepatterns for a graphoepitaxial pitch-multiplication application. Macromolecules 2015, 48, 1256−1261. (52) E, W.; Ren, W.; Vanden-Eijnden, E. String method for the study of rare events. Phys. Rev. B 2002, 66, No. 052301. J

DOI: 10.1021/acs.macromol.9b00194 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules (53) Ceniceros, H. D.; Fredrickson, G. H. Numerical Solution of Polymer Self-Consistent Field Theory. Multiscale Model. Simul. 2004, 2, 452−474.

K

DOI: 10.1021/acs.macromol.9b00194 Macromolecules XXXX, XXX, XXX−XXX