Inverse Design of Bulk Morphologies in Multiblock Polymers Using

Aug 28, 2017 - Search; Citation; Subject .... We report on an inverse design strategy to target bulk morphologies utilizing particle swarm optimizatio...
0 downloads 0 Views 1MB Size
Article pubs.acs.org/Macromolecules

Inverse Design of Bulk Morphologies in Multiblock Polymers Using Particle Swarm Optimization Mihir R. Khadilkar,† Sean Paradiso,‡ Kris T. Delaney,† and Glenn H. Fredrickson*,†,‡,§ †

Materials Research Laboratory, ‡Department of Chemical Engineering, and §Materials Department, University of California, Santa Barbara, Santa Barbara, California 93106, United States ABSTRACT: Multiblock polymers are a fascinating class of soft materials owing to their spontaneous self-assembly into a variety of ordered mesophases at the nanoscale. However, controlling the ordering and hence the properties is difficult due to a vast number of design parameters including the choice of monomer species, block sequence, block molecular weights and dispersity, polymer architecture, and binary interaction parameters. Navigating through this design space in order to find an optimal formulation for a target property requires an approach that efficiently searches through the countless parameter combinations in an automated fashion. We report on an inverse design strategy to target bulk morphologies utilizing particle swarm optimization (PSO) as a global optimizer and selfconsistent-field theory (SCFT) as a forward prediction engine. To avoid metastable states in the forward prediction, we utilize pseudospectral variable cell SCFT initiated from a library of defect-free seeds of known block copolymer morphologies. We demonstrate that our approach allows for a robust identification of block copolymers and copolymer alloys that self-assemble into a targeted structure, optimizing parameters such as block fractions, blend fractions, and Flory χ parameters.



optimization29 (PSO), a search and optimization methodology.29 PSO was originally inspired by biological systems that self-aggregate based on some fixed local rules, and it has been previously used for global optimization problems in various fields including network design,30 robotics,31 and scheduling.32 Many technological applications of MBPs including gas separations33 and efficient laser reflectors34 require bulk ordering at the nanoscale, for which such an inverse design scheme has been lacking. This paper extends the same approach of coupling PSO with SCFT to target bulk MBP morphologies. For the inverse design of MBPs to be robust, one has to deal with two separate global optimization problems: first in the forward direction, that of finding a morphology associated with the global minimum of free energy at a given combination of formulation parameters, and the second in the inverse direction, that of optimizing those formulation parameters in order to stabilize a target morphology. The former has been a longstanding problem due to a high density of metastable states around the true free energy minimum.12 Hence, depending on the initial field configuration, SCFT could potentially converge to a multitude of solutions corresponding to competing phases. A successful inverse design method needs to resolve this issue, so that repeated SCFT calculations at a given point in the formulation space yield the same consistent prediction. For the problems addressed by Paradiso et al.,28 this was not an obstacle due to the reduced dimensionality. However, in three-dimensional bulk systems, due to a rough free energy landscape, even the most advanced field relaxation

INTRODUCTION Self-assembly and directed assembly of soft materials have emerged as powerful tools to dictate material structure and properties in a wide range of applications. In particular, multiblock polymers (MBPs), due to their microphase separation, allow the control of structural ordering at length scales not otherwise feasible.1,2 Competition between incompatibility of unlike monomer species and the entropic cost of stretching polymer chains gives rise to a rich range of ordered structures in MBPs. Further, advanced synthesis techniques provide tools for generating a wide array of macromolecules. These developments have opened up numerous opportunities for potential applications including lithography,3−5 membranes,6 photonic crystals,7 and drug delivery.8 Numerical simulations are extremely important and effective in analyzing and understanding MBP phase behavior.1,9−17 A common research paradigm is to start from a given set of design parameters (MBP architecture, block and blend fractions, etc.) and determine the equilibrium morphologies and properties. A more pragmatic approach, however, takes the opposite route: starting from a target morphology and searching for MBP formulations that self-assemble into the same. Such an inverse design, although highly desirable, is challenging since it involves automating a search through an enormous design space12 in a robust, fast, and reliable manner. The inverse design approach has been explored in a wide range of different soft materials.18−23 In MBPs, design via graphoepitaxy or chemoepitaxy prepatterns has been investigated24−27 wherein external templates were optimized to direct self-assembly of thin films of block copolymers into desired targets. On the other hand, a recent effort28 attempted to directly optimize architecture and composition variables in MBP thin films by combining SCFT with particle swarm © 2017 American Chemical Society

Received: June 7, 2017 Revised: August 8, 2017 Published: August 28, 2017 6702

DOI: 10.1021/acs.macromol.7b01204 Macromolecules 2017, 50, 6702−6709

Article

Macromolecules

ϕk, with values smaller than a threshold (in our case 0.005) being set to zero. As mentioned before, the f itness for PSO dynamics is calculated using SCFT, which has long been established as a powerful method to study multiblock polymers.1,9,10,36 Starting from a microscopic model of the polymer chains (typically, but not restricted to continuous Gaussian chains), SCFT uses Hubbard−Stratonovich transforms37 to decouple the manybody interactions between different polymer segments, replacing them with effective one-body interactions involving auxiliary fields. The resulting field theory is then solved selfconsistently using a mean-field approximation. Such a treatment allows for a fast equilibration as well as direct access to the free energy (which is essential in comparing the relative stability of competing phases). Because of the challenges in forward prediction mentioned earlier, we need free energies of several potential candidate phases. We achieve this by starting from a library of unit-cell seed fields and performing separate SCFT calculations for each of the phases. Using a field symmetrizer, symmetry operations of the respective space group are enforced during SCFT. The library of seeds is built on several prior test calculations for each of the individual morphologies in their known regions of stability. After confirming the morphologies visually and through the use of a field symmetrizer, the respective converged fields are used as starting seeds in the simulations for that morphology. The objective function is defined as Ω = FBest alternate − FTarget, where FBest alternate and FTarget are the converged intensive free energies of the best alternate candidate phase and the target phase, respectively. The sign convention dictates that in the favorable scenario Ω is positive and high. SCFT Details. We work with a purely incompressible model of continuous Gaussian chains and perform unit-cell calculations with periodic boundary conditions. By using a variablecell SCFT method,38,39 changes in the size and shape of the cell are allowed (at a fixed concentration) during convergence to a solution. Working in a canonical ensemble is a convenient choice. However, addition of homopolymers to pure MBPs can sometimes lead to a macrophase separation (MPS). In this case, the simulation cell in a canonical ensemble often grows unbounded as it attempts to reduce the free energy cost of the interface between the two macrodomains and, as a result, fails to converge. Working in a grand-canonical ensemble (GCE) allows one to circumvent this problem, since the coexistence line in a canonical ensemble collapses to a single point in GCE as the two coexisting microphases have the same chemical potential (or thermodynamic activity). For blend systems, working in GCE amounts to changing the PSO search variable from blend fraction to thermodynamic activity. Further, due to the incompressibility in our model, one of the activities is redundant and hence in a diblock− homopolymer mixture, we fix the diblock activity to unity and only vary the homopolymer activity. Since the MPS region is absent in GCE (as it collapses to a single line), SCFT runs converge easily and retrieve accurate fitness (which now is the grand potential instead of the free energy). Robustness. To improve SCFT convergence and increase the reliability of the fitness estimation, we employ several robustness checks. First, the converged simulation cell size from the previous PSO step of an agent is used as the starting cell for its current PSO iteration, since it is likely to be close to the ideal cell size at the new location. By seeding this way, we observe that for more than 85% of the total SCFT runs the final

schemes in SCFT require a priori knowledge of the stable phase. We address the challenging forward problem for the bulk 3D case by building a library of defect-free unit-cell seeds for all known ordered morphologies. Separate SCFT simulations are performed for each of the morphology candidates with appropriate symmetry constraints to ensure an accurate estimation of the equilibrium structure. This combined with the PSO search heuristic allows our method to find the formulation parameters stabilizing a given target morphology. We work with single MBP as well as blend systems since alloys of two or more block copolymers often provide an easy and economical route to new materials with improved properties.



METHOD AND IMPLEMENTATION As described earlier, our work combines SCFT with particle swarm optimization29,35 (PSO), a search heuristic in which a collection of optimizers in collaboration via weak local interaction reach a global optimum of the objective function. Using SCFT, the f itness of a phase at a given point in formulation space is calculated and the PSO algorithm is used to find the solutions with improved fitness. In PSO, a set of optimizers (“particles” or agents) move in a D-dimensional search space. An individual agent i is composed of four vectors: its position xi ⃗ = (xi1, xi2, ..., xiD), its velocity vi ⃗ = (vi1, vi2, ..., viD), the best position it has individually found pi ⃗ = (pi1, pi2, ..., piD), and the best position found by any of its neighbors (determined by a network topology), ni ⃗ = (ni1, ni2, ..., niD). The best position found is defined by the highest value of the objective function (defined below), which in our case is related to the intensive free energy of the candidate phases. Agents are initialized at random locations in the parameter space and are updated at each iteration according to the following equations for each agent i and in each dimension d: vidn+ 1 = vidn + χc0ϕ0(pid − xid) + χc1ϕ1(nid − xid) − (1 − χ )vid xidn+ 1 = xidn + vid (1)

where ϕ0 and ϕ1 are independent, uniformly distributed random variables in the interval [0, 1] generated at every update, and c0 and c1 are acceleration coefficients. The parameter χ ∈ [0, 1] is known as the constriction factor. Through the PSO dynamics, the swarm is attracted toward historically fit solutions. For the results shown here, we use between 4 and 20 agents arranged in a ring (each agent connected to two other agents with periodic boundary conditions), χ = 0.8 and c0 = c1 = 2.05 (values close to the “standard” version of the PSO35). While these parameters do influence the PSO trajectory, their effect on the overall performance of the heuristic optimization is rather weak. The relevant variables for our optimization are the polymer architecture parameters, namely chain block fractions f j, blend fractions ϕk (in the case of multiple chains), and interaction strengths. In particular, since block and blend fractions are constrained variables (f j ≥ 0 and ∑f j = 1), we define transformed variables αj = log(f j/f j+1) and βk = log(ϕk/ ϕk+1) to constrict the relevant range. Using this scheme, we have (nb − 1) independent α variables for each chain with nb blocks and (nc − 1) independent β variables for nc distinct chains. Extreme cases where any of the f j or ϕk are zero correspond to an αj or βk being ±∞. However, in practice, our limits on αj and βk are set by our resolution in representing f j or 6703

DOI: 10.1021/acs.macromol.7b01204 Macromolecules 2017, 50, 6702−6709

Article

Macromolecules

examples here, we use between 4 and 12 agents. At each step in the PSO algorithm, the fitness is calculated at each of the agent positions using separate SCFT simulations for the target phase and each of the candidate phases. We included lamellar, hexagonal cylinder, gyroid, BCC, FCC, O70, and disordered phase in our list of candidate morphologies. PSO results for a body-centered-cubic (BCC) sphere phase as a target in a pure AB diblock copolymer using four agents are shown in Figure 1. Initial solutions found by the swarm are in

cell size of the simulation is within 5% of the starting cell size. If the target morphology is not even metastable at a given location (instead converging to a disordered melt), we assign a low fitness to effectively remove its influence on the swarm. A slightly more subtle problem may arise, leading to an erroneous estimation of the fitness. The field symmetrizer used in the respective simulations cannot prevent additional symmetries from being introduced in the fields, which can sometimes lead to convergence to a supergroup (for example, a simulation for the gyroid phase or the Frank−Kasper A15 phase converging to a BCC sphere phase) and has been observed before.40 This is avoided by performing a structure similarity check for the target phase, in which we calculate the absolute value of the inner product of the normalized Fourier component vectors of the seed density fields and the converged fields. The presence of additional symmetries in a supergroup manifests as a different pattern of Fourier modes in the density field. Hence, in the inner product of the seed density fields and the converged fields, a value close to 1 indicates the expected structure, while a value significantly less than 1 may indicate convergence to a supergroup. For diblocks, we observed that the distribution of these inner-product values was strongly bimodal, with most values clustered around 0 or 1. In such a scenario, a threshold value of 0.8 for the inner product was effective in detecting supergroup cases. These supergroup cases were assigned an arbitrary low fitness (−10000) similar to those simulations that converged to a homogeneous disordered state, so that they do not affect the swarm. The distribution is less bimodal for complex sphere phases. For those systems, an a posteriori check for the solutions found by the PSO search in addition to the similarity test is useful to ensure accuracy. An additional issue is the similarity check for A-minority and B-minority seeds for the same ordered morphology. In the current implementation, we do not distinguish between A-minority and B-minority phases of the same symmetry (by taking the absolute value of the innerproduct), however, that could be modified to select only one of the two. In spite of these checks, one is likely to encounter occasional errors. These may arise from unstable SCFT trajectories or simulations that have not fully converged (in case we exhaust the prescribed maximum number of SCFT iterations). These are kept at a minimum by choosing an appropriate time step to ensure SCFT stability, using sufficient SCFT iterations, and a fast cell-updater in variable-cell relaxation.

Figure 1. Results of a PSO search for the BCC phase in pure AB diblock copolymer with 4 agents and 200 iterations varying block fraction fA and interaction strength (χN). (a) A scatterplot of all agent positions throughout the search, overlaid upon the known diblock copolymer phase diagram.42 Different colors represent different agents, and the brightness of the color represents the evolution with early times representing brighter hues and later times representing dark colors. (b) Trajectory of the cumulative best position of the swarm as a whole. Known stable phase regions BCC (S), cylinder (C), gyroid (G), lamellar (L), and disordered (DIS) from the literature are marked accordingly. (c) A plot of best fitness of the swarm as a function of PSO iterations.



RESULTS While classic MBP morphologies like lamellae and hexagonal cylinders are relatively ubiquitous, other bulk morphologies like 3D network phases as well as sphere phases occur in fairly narrow regions in parameter space and often have interesting structural and optical41 properties. Hence, it is instructive to test PSO with these bulk morphologies as targets. Most results shown in this section involve pure diblock copolymer and a diblock−homopolymer blend, which are relatively simple systems that we use as a benchmark. For the diblock examples, the search variables are the block fraction fA of the diblock as well as the interaction strength χN. For the diblock−homopolymer blend, the block fraction fA of the AB diblock and the thermodynamic activity of the homopolymer (A) are PSO variables. The number of PSO agents required depends on the dimensionality of the search space, and in the relatively simple two- and three-dimensional

the disordered region, as in Figure 1b, but soon as the agents collaborate on their fitness information, the swarm finds fitter solutions, eventually converging to the expected BCC region. The fitness evolution of the swarm as plotted in Figure 1c shows how swarm improves its fitness as it evolves under PSO dynamics. Qualitatively, large jumps in fitness are interpreted as the swarm moving away from a local attraction basin to a deeper basin (when potentially one of the agents stumbles upon a new, fitter region), while smaller improvements are associated with optimization within the basin. 6704

DOI: 10.1021/acs.macromol.7b01204 Macromolecules 2017, 50, 6702−6709

Article

Macromolecules

Figure 2. Results of a PSO search for a gyroid (Ia3̅d) phase in pure AB diblock copolymer with four agents and 150 iterations, varying block fraction fA and interaction strength (χN). (a) A scatterplot of all agent positions throughout the search, overlaid upon the known diblock copolymer phase diagram.42 Different colors represent different agents, and the brightness of the color represents the evolution with early times representing brighter hues and later times representing dark colors. (b) Trajectory of the cumulative best position found by the swarm. Known stable phase regions BCC (S), cylinder (C), gyroid (G), lamellar (L), and disordered (DIS) from the literature are marked accordingly.

Figure 3. Results of a PSO search for orthorhombic O70 phase pure AB diblock copolymer with 12 agents and 150 iterations, varying A blockfraction fA and interaction strength χN. (a) A scatterplot of all agent positions during the search, plotted against the known diblock copolymer phase diagram.42 (b) Locations of cumulative best solutions found by the swarm. Known stable phase regions BCC (S), cylinder (C), gyroid (G), lamellar (L), disordered (DIS), and O70 (O) from the literature are marked accordingly.

PSO search required around 30 steps with 12 agents to obtain its first stable solution. A brute force approach in comparison would require at least a 20 × 20 grid for the initial coarse search and a similar grid for a fine search in a smaller region to pinpoint the stable O70 region. Thus, in the present 2D example, the brute force approach would require approximately 800 SCFT simulations for each of the candidate morphologies, whereas the PSO required about 360 SCFT simulations for each competing structure, i.e., a factor of 2 improvement in efficiency. PSO for Blends. Blending different polymers together is one of the most common routes used to tune structural and material properties. Thus, it is greatly beneficial if an inverse design method can work with multicomponent search spaces. As explained earlier, for blends we chose to work in a grandcanonical ensemble (GCE). As an example, we show PSO results for one of the simplest blend systems studied in the literature, that of homopolymer A and AB diblock copolymer. This system has been shown to exhibit a large variety of mesophases in previous studies11,43 and hence is well-suited for a trial investigation. The search space is restricted to only the block fraction fA of the diblock and the homopolymer activity that controls the homopolymer blend fraction (ϕA). Here, the relative length of the homopolymer and χN are also potentially

Similar results are obtained for other targets in AB diblock copolymers. Figure 2 shows results of a PSO search for a gyroid network phase (space group Ia3̅d). Because of the attraction terms in the PSO dynamics, agents tend to oscillate around the local best agent position, until one of the agents to finds a better solution. This balance between exploration and exploitation is how the swarm finds globally optimum solutions. We also tested the PSO algorithm for the orthorhombic O70 phase. This network phase has a minuscule stable region in the AB diblock copolymer phase diagram and was only added after its experimental discovery in linear triblock copolymers.38 PSO search results for this phase are shown in Figure 3, where the PSO agents update their estimate through several basins and ultimately converge to the expected region. We note that for the two-dimensional search examples shown here the computational cost of PSO is comparable to an exhaustive search. However, the exhaustive search quickly becomes impractical as we increase the number of search parameters and hence the dimensionality of the search space. A quantitative measure of efficiency can be developed by comparing the number of PSO steps required to find a solution to the fraction of the area occupied by the target phase in the phase diagram. In the case of O70, the target phase occupied approximately 0.15% of our search region, and the 6705

DOI: 10.1021/acs.macromol.7b01204 Macromolecules 2017, 50, 6702−6709

Article

Macromolecules

Figure 4. (a) Results of a PSO search for the gyroid (Ia3̅d) phase with four agents and 200 iterations in diblock−homopolymer blend (AB + A) at χN = 12.0 and α = 1.0, plotted as a function of the block fraction fA of the diblock and the blend fraction ϕA of the homopolymer, against a known phase diagram from Matsen.43 Although the activity of the homopolymer was used as a PSO variable in the grand-canonical ensemble, we show the homopolymer blend fraction for easy comparison with the known literature. (b) Trajectory followed by the cumulative best positions found by the PSO swarm. Known stable phase regions BCC (S), FCC (F), cylinder (C), gyroid (G), lamellar (L), disordered (DIS), and MPS region are marked.

Figure 5. (a) The scatterplot of agent positions in a PSO search for Frank−Kasper σ phase in SISO at a fixed f S/f I = 0.93, with four agents and 200 iterations, plotted as a function of O-block fraction f O and the S-block asymmetry τ defined as f S1/(f S1 + f S2), against a known phase diagram47 showing respective stable phases. (b) Evolution of swarm fitness as a function of PSO iterations in a three-dimensional PSO search for the same σ phase with 12 agents and 200 iterations allowing f S/f I to also vary. (c) Scatterplot of all agent positions (in magenta) in the 3D search. The projections on three planes are shown in red, green, and blue for clarity. (d) Scatterplot showing only the stable points (in magenta) found in the 3D search, revealing a different, fitter region from that obtained in a 2D slice in (a).

of the PSO agents superimposed with the phase diagram from Matsen.43 As can be seen, after scanning different parts of the parameter space, agents converge on the known stable region for the gyroid phase. Working in the GCE leads to a much more efficient treatment of the inverse problem in the blend system. Although we see a few points in Figure 4a in the MPS region, they are

tunable (making it a four-dimensional search); however, by fixing these, we can compare our results with the twodimensional phase diagrams found in the literature.43 In this diblock−homopolymer blend system, we search for the gyroid phase (with space group Ia3d̅ ) at a fixed χN = 12.0 and a fixed relative length ratio α = 1.0 and include other common phases (namely lamellar, BCC, FCC, and hexagonal cylinders) as alternate candidates. Figure 4 shows the positions 6706

DOI: 10.1021/acs.macromol.7b01204 Macromolecules 2017, 50, 6702−6709

Article

Macromolecules

parameters selected in the implementation including the number of agents, the connection topology, and constriction and inertia coefficients. We have fixed the values of constriction and the inertia coefficients to those used in previous successful PSO studies in the literature.28,35 The efficiency of the search also depends on the number of agents. While using more agents helps to more broadly access the parameter space, the computational expense increases significantly with increasing number of agents and after a point does not result in any reduction in the number of iterations required for successful prediction. For example, the search for the gyroid phase with eight agents requires apprximately 80% more SCFT simulations to reach its first stable solution compared to a search with four agents. The best choice for agent count also depends on the dimensionality of the search space. While four agents were sufficient in a two-dimensional search for the Frank−Kasper σ phase, an effective search in three dimensions for the same target required 12 agents, as shown in Figure 5. The connection topology of the PSO agents is an additional component in the implementation. We tested a ring graph as well as a 2D grid graph (a 2D lattice of m × n agents, each of which is connected to its four nearest neighbors). We observed no significant difference between them in performance, as measured by the number of iterations required to reach the optimal solution. As an example, a 4 × 3 grid of 12 PSO agents required 7% fewer PSO iterations to find the optimal solution than a 12 agent PSO search on a ring topology for the O70 search in the AB diblock copolymer. Admittedly, for a small number of agents (as is the case for most of the examples shown) the two topologies tested are similar; however, the difference between them (qualitatively based on the number of connections) grows as we increase the number of agents. It remains to be seen whether the connection topology has a significant effect on the performance for more complex search spaces using a large number of agents.

due to SCFT simulations that have not converged fully, leading to an inaccurate estimation of the blend fraction ϕA. Frank−Kasper Phases in Tetrablock Terpolymers. As a final test of the method, we applied PSO to an ABA′C tetrablock terpolymer. This terpolymer architecture has recently been shown to exhibit a rich array of low-symmetry Frank−Kasper sphere morphologies in the case of poly(styrene)-b-poly(isoprene)-b-poly(styrene)′-b-poly(ethylene oxide) (SIS′O).44,45 A generic tetrablock terpolymer (a block copolymer with four blocks and three distinct species) has a large number of design parameters, including block ordering (ABCA, ABAC, ACBA, etc.), choice of block chemistry (that sets the binary interaction parameters), and relative block length. Thus, the success of our PSO search in such a case would be an important step toward more complicated search problems. To better match the experimental synthesis strategies used for these tetrablocks,46 we fix the bare χ parameters (by fixing the chemistries and the temperature) and the molecular weight of the parent S1IS2 triblock from which the tetrablocks are synthesized. There are still three free parameters describing the chain: the block fraction of the O-block added, f O, the S-block asymmetry τ = f S1/(f S1 + f S2) in the parent triblock, and the ratio f S/f I (where f S here is the total f S1 + f S2). Without an optimization tool like PSO, even a two-dimensional parameter sweep (by fixing one of these parameters) is relatively laborious. A PSO search to find the Frank−Kasper σ phase was carried out, both in the full three-dimensional search space and in a two-dimensional slice with fixed f S/f I for an easy comparison with a known phase diagram.47 We included other sphere phases, namely A15, C15, BCC, and HCP, as well as lamellar and hexagonal cylinders as candidate phases. Working at a fixed temperature of 180 °C and a fixed molecular weight of 20.8 kDa of the parent triblock S1IS2 (corresponding to N ≈ 299), we replicated the conditions in Chanpuriya et al.47 Figure 5 shows the results for twodimensional and three-dimensional searches for the Frank− Kasper σ phase. Panel a shows a PSO search as a function of f O and τ (at a fixed f S/f I = 0.93), comparing it to the known phase diagram.47 In a full three-parameter search, PSO not only finds the stable region in the original slice but also reveals a different, more stable region at f S/f I = 1.18 (see Figure 5b), which would have been significantly more difficult to find in a comprehensive parameter sweep. As a comparison, even a coarse 20 × 20 × 20 mesh for a parameter sweep would require 8000 SCFT simulations of all competing phases; in contrast, only 2400 SCFT simulations were used in the PSO search. While we have not shown examples in even higher dimensions, the efficiency of the PSO method is likely to improve even further as the dimensionality increases. While no analytical expression exists for the optimum swarm size for a general search problem in higher dimensions, on standard benchmark optimization problems in around 10−30 dimensions, 20−80 agents have routinely been shown to be adequate.48,49 Hence as an example, in a 10-dimensional search even with the larger number of required iterations, PSO is likely to convincingly outperform brute-force search, which would require on the order of 1010− 1013 SCFT simulations per morphology candidate. This underscores an important point made earlier, that as we increase the dimensions in the search, PSO becomes increasingly more effective than a brute-force alternative. Parameter Selection for PSO. Like many global optimization methods, the performance of PSO depends on



CONCLUSION We have presented a computational inverse design method capable of stabilizing a desired bulk multiblock polymer or blend morphology by searching through a large range of design parameters including chain block architecture, blend fractions, and interaction strengths. Combining self-consistent field theory with a PSO search heuristic, this method can efficiently scan and locate optimal parameter combinations, which is particularly useful in high-dimensional search spaces. By combining various robustness checks in SCFT simulations with the collaborative information sharing of PSO agents, the method effectively targets a given morphology. The imperfections in the forward-prediction engine (in determining a global free-energy minimum at a given point in parameter space) present a significant challenge in maintaining accuracy in fitness calculations. However, the checks imposed on the simulations (including the similarity test) greatly improve performance. Nonetheless, we need to include all relevant competing morphologies in our fitness estimation to avoid false-positives. In many cases, the physical intuition about the likely candidates is sufficient to ensure success; for example, when targeting complex Frank−Kasper phases, it is generally sufficient to include all known sphere phases as well as the more common lamellar and cylinder phases.47,50 In problems where an important candidate morphology is missing, one may make a wrong prediction. However, this is where close collaboration between theory and experiments would prove 6707

DOI: 10.1021/acs.macromol.7b01204 Macromolecules 2017, 50, 6702−6709

Article

Macromolecules ORCID

to be a strong asset. Contradictory results from the experiments would help build a richer library of phases and improve reliability in the long term. The algorithm presented is quite flexible when extending to more complicated search spaces. Because of the nature of PSO optimization, one can include all relevant design parameters that would potentially affect the stability and let PSO agents determine the important parameter combinations. For example, in the thin film pattern selection problems considered by Paradiso et al.,28 PSO simulations started from a generic threespecies triblock copolymer sometimes resulted in the elimination of one of the blocks. This removes the burden of selecting factors important for a particular target right at the beginning. Hence, one might start with highest level of block complexity and interaction parameters and let PSO perform a coarse search. Refinement could then be done either with the PSO in a restricted space or with a brute-force mapping of the same subspace. The PSO is also flexible in terms of choice of target. At times, one may only be interested in stabilizing any new phase in a search space. In such a scenario, we are looking for any of a number of different morphologies. With a modification to the fitness function to have a local maximum for each of the individual targets, the same method can be used to search for the combined target. We note that in its conventional form PSO only accommodates a search over continuous variables. However, to search over a discrete degree of freedom, a possible solution is to suitably transform it onto a convenient continuous variable. Still, there are numerous avenues worth exploring in future to improve the method for more complex targets. An alternate chain representation in PSO that accounts for changing architecture (adding/removing branches, permuting the order of blocks) would be beneficial. One possibility is to treat a chain as a graph with variable nodes and connections representing the blocks. The freedom to permute the blocks is still technically available by fully varying χN parameters, but not in a restricted way that one might desire. The structure similarity test is also less robust in the case of more complex phases. To alleviate this problem, tools like machine learning and artificial neural networks could be used to identify “similar” structures. Another avenue is modifying the fitness function to find regions of stability for a target morphology. This would increase the chances of replicating the PSO predictions in experiments where fine control over the block and blend fractions is difficult to achieve. It is also conceivable that for more challenging search spaces PSO dynamics may be modified in a variety of different ways51 to avoid early convergence of the swarm. As more technological applications of MBPs are explored, inverse design methods including the one presented here would be extremely useful to efficiently select formulations for nanomaterials with tailored morphologies. One of the ways of further broadening the use of PSO is by targeting properties directly (for example, photonic bandgap) rather than through structure. Apart from applications, this approach could also be used to better understand structure−property relationships in MBPs at a more fundamental level.



Mihir R. Khadilkar: 0000-0001-7238-2886 Glenn H. Fredrickson: 0000-0002-6716-9017 Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS M.R.K. thanks the Dow Chemical Company for financial support of this research through the Dow Materials Institute at UCSB. K.T.D. was partially supported by the National Science Foundation DMREF-1332842. Extensive use was made of the computational facilities of the Center for Scientific Computing at the California NanoSystems Institute (CNSI) and Materials Research Laboratory (MRL) at U.C. Santa Barbara: an NSF MRSEC (DMR-1121053) and NSF CNS-0960316.



REFERENCES

(1) Fredrickson, G. The Equilibrium Theory of Inhomogeneous Polymers; Oxford University Press: New York, 2006. (2) Bates, F. S.; Fredrickson, G. H. Block Copolymers - Designer Soft Materials. Phys. Today 1999, 52, 32−38. (3) Bates, C. M.; Maher, M. J.; Janes, D. W.; Ellison, C. J.; Willson, C. G. Block copolymer lithography. Macromolecules 2014, 47, 2−12. (4) Luo, M.; Epps, T. H., III Directed block copolymer thin film selfassembly: emerging trends in nanopattern fabrication. Macromolecules 2013, 46, 7567−7579. (5) Ginzburg, V. V.; Weinhold, J. D.; Trefonas, P. Computational modeling of block-copolymer directed self-assembly. J. Polym. Sci., Part B: Polym. Phys. 2015, 53, 90−95. (6) Jackson, E. A.; Hillmyer, M. A. Nanoporous membranes derived from block copolymers: From drug delivery to water filtration. ACS Nano 2010, 4, 3548−3553. (7) Fink, Y.; Urbas, A. M.; Bawendi, M. G.; Joannopoulos, J. D.; Thomas, E. L. Block copolymers as photonic bandgap materials. J. Lightwave Technol. 1999, 17, 1963−1969. (8) Yin, L.; Dalsin, M. C.; Sizovs, A.; Reineke, T. M.; Hillmyer, M. A. Glucose-functionalized, serum-stable polymeric micelles from the combination of anionic and RAFT polymerizations. Macromolecules 2012, 45, 4322−4332. (9) Matsen, M. W.; Schick, M. Stable and unstable phases of a diblock copolymer melt. Phys. Rev. Lett. 1994, 72, 2660−2663. (10) Drolet, F.; Fredrickson, G. H. Combinatorial screening of complex block copolymer assembly with self-consistent field theory. Phys. Rev. Lett. 1999, 83, 4317. (11) Matsen, M. W. New fast SCFT algorithm applied to binary diblock copolymer/homopolymer blends. Macromolecules 2003, 36, 9647−9657. (12) Bates, F. S.; Hillmyer, M. A.; Lodge, T. P.; Bates, C. M.; Delaney, K. T.; Fredrickson, G. H. Multiblock polymers: panacea or Pandora’s box? Science 2012, 336, 434−440. (13) Delaney, K. T.; Fredrickson, G. H. Polymer field-theory simulations on graphics processing units. Comput. Phys. Commun. 2013, 184, 2102−2110. (14) Guo, Z.; Zhang, G.; Qiu, F.; Zhang, H.; Yang, Y.; Shi, A.-C. Discovering ordered phases of block copolymers: new results from a generic Fourier-space approach. Phys. Rev. Lett. 2008, 101, 028301. (15) Matsen, M. W. The standard Gaussian model for block copolymer melts. J. Phys.: Condens. Matter 2002, 14, R21. (16) Shi, A.-C.; Li, B. Self-assembly of diblock copolymers under confinement. Soft Matter 2013, 9, 1398−1413. (17) Xu, W.; Jiang, K.; Zhang, P.; Shi, A.-C. A strategy to explore stable and metastable ordered phases of block copolymers. J. Phys. Chem. B 2013, 117, 5296−5305. (18) Khalatur, P. G.; Khokhlov, A. R.; Krotova, M. K. Evolutionary approach in copolymer sequence design. Macromol. Symp. 2007, 252, 36−46.

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected] (G.H.F.). 6708

DOI: 10.1021/acs.macromol.7b01204 Macromolecules 2017, 50, 6702−6709

Article

Macromolecules (19) Torquato, S. Inverse optimization techniques for targeted selfassembly. Soft Matter 2009, 5, 1157−1173. (20) Cohn, H.; Kumar, A. Algorithmic design of self-assembling structures. Proc. Natl. Acad. Sci. U. S. A. 2009, 106, 9570−9575. (21) Jain, A.; Bollinger, J. A.; Truskett, T. M. Inverse methods for material design. AIChE J. 2014, 60, 2732−2740. (22) Miskin, M. Z.; Jaeger, H. M. Evolving design rules for the inverse granular packing problem. Soft Matter 2014, 10, 3708−3715. (23) Jadrich, R.; Lindquist, B.; Truskett, T. Probabilistic inverse design for self-assembling materials. J. Chem. Phys. 2017, 146, 184103. (24) Latypov, A.; Preil, M.; Schmid, G.; Xu, J.; Yi, H.; Yoshimoto, K.; Zou, Y. Exploration of the directed self-assembly based nanofabrication design space using computational simulations. Proc. SPIE 2013, 868013. (25) Hannon, A. F.; Ding, Y.; Bai, W.; Ross, C. A.; Alexander-Katz, A. Optimizing Topographical Templates for Directed Self-Assembly of Block Copolymers via Inverse Design Simulations. Nano Lett. 2014, 14, 318−325. (26) Khaira, G. S.; Qin, J.; Garner, G. P.; Xiong, S.; Wan, L.; Ruiz, R.; Jaeger, H. M.; Nealey, P. F.; de Pablo, J. J. Evolutionary Optimization of Directed Self-Assembly of Triblock Copolymers on Chemically Patterned Substrates. ACS Macro Lett. 2014, 3, 747−752. (27) Ouaknin, G.; Laachi, N.; Delaney, K.; Fredrickson, G.; Gibou, F. Shape optimization for DSA. Proc. SPIE 2016, 97770Y−97770Y. (28) Paradiso, S. P.; Delaney, K. T.; Fredrickson, G. H. Swarm Intelligence Platform for Multiblock Polymer Inverse Formulation Design. ACS Macro Lett. 2016, 5, 972−976. (29) Eberhart, R. C.; Kennedy, J. Proceedings of the sixth international symposium on micro machine and human science 1995, 1, 39−43. (30) Juang, C.-F. A hybrid of genetic algorithm and particle swarm optimization for recurrent network design. IEEE Transactions on Systems, Man, and Cybernetics, Part B (Cybernetics) 2004, 34, 997− 1006. (31) Pugh, A.; Martinoli, J. Inspiring and Modeling Multi-Robot Search with Particle Swarm Optimization. 2007 IEEE Swarm Intelligence Symposium 2007, 332−339. (32) Pandey, S.; Wu, L.; Guru, S. M.; Buyya, R. A Particle Swarm Optimization-Based Heuristic for Scheduling Workflow Applications in Cloud Computing Environments. 2010 24th IEEE International Conference on Advanced Information Networking and Applications; 2010; pp 400−407. (33) Choi, S.; Drese, J. H.; Jones, C. W. Adsorbent materials for carbon dioxide capture from large anthropogenic point sources. ChemSusChem 2009, 2, 796−854. (34) Yoon, J.; Lee, W.; Thomas, E. L. Optically pumped surfaceemitting lasing using self-assembled block-copolymer-distributed Bragg reflectors. Nano Lett. 2006, 6, 2211−2214. (35) Bratton, D.; Kennedy, J. Defining a standard for particle swarm optimization. 2007 IEEE swarm intelligence symposium; 2007; pp 120−127. (36) Leibler, L. Theory of microphase separation in block copolymers. Macromolecules 1980, 13, 1602−1617. (37) Stratonovich, R. L. Sov. Phys. Solid State 1958, 2, 1824. (38) Tyler, C. A.; Morse, D. C. Orthorhombic Fddd Network in Triblock and Diblock Copolymer Melts. Phys. Rev. Lett. 2005, 94, 208302. (39) 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. (40) Tyler, C. A.; Qin, J.; Bates, F. S.; Morse, D. C. SCFT study of nonfrustrated ABC triblock copolymer melts. Macromolecules 2007, 40, 4654−4668. (41) Chan, C. T.; Ho, K.; Soukoulis, C. Photonic band gaps in experimentally realizable periodic dielectric structures. EPL (Europhysics Letters) 1991, 16, 563. (42) Matsen, M. W. Effect of Architecture on the Phase Behavior of AB-Type Block Copolymer Melts. Macromolecules 2012, 45, 2161− 2165.

(43) Matsen, M. Phase behavior of block copolymer/homopolymer blends. Macromolecules 1995, 28, 5765−5773. (44) Lee, S.; Bluemle, M. J.; Bates, F. S. Discovery of a Frank-Kasper σ phase in sphere-forming block copolymer melts. Science 2010, 330, 349−353. (45) Zhang, J.; Bates, F. S. Dodecagonal Quasicrystalline Morphology in a Poly (styrene-b-isoprene-b-styrene-b-ethylene oxide) Tetrablock Terpolymer. J. Am. Chem. Soc. 2012, 134, 7636−7639. (46) Touris, A.; Chanpuriya, S.; Hillmyer, M. A.; Bates, F. S. Synthetic strategies for the generation of ABCA’type asymmetric tetrablock terpolymers. Polym. Chem. 2014, 5, 5551−5559. (47) Chanpuriya, S.; Kim, K.; Zhang, J.; Lee, S.; Arora, A.; Dorfman, K. D.; Delaney, K. T.; Fredrickson, G. H.; Bates, F. S. Cornucopia of scale ordered phases in sphere-forming tetrablock terpolymers. ACS Nano 2016, 10, 4961−4972. (48) Shi, Y.; Eberhart, R. C. Empirical study of particle swarm optimization. Evolutionary Computation, 1999. CEC 99. Proceedings of the 1999 Congress, 1999; pp 1945−1950. (49) Nápoles, G.; Grau, I.; Bello, R. Constricted Particle Swarm Optimization based algorithm for global optimization. Polibits 2012, 46, 5−11. (50) Xie, N.; Liu, M.; Deng, H.; Li, W.; Qiu, F.; Shi, A.-C. Macromolecular metallurgy of binary mesocrystals via designed multiblock terpolymers. J. Am. Chem. Soc. 2014, 136, 2974−2977. (51) Nakisa, B.; Nazri, M. Z.; Rastgoo, M. N.; Abdullah, S. A survey: Particle swarm optimization based algorithms to solve premature convergence problem. J. Comput. Sci. 2014, 10, 1758.

6709

DOI: 10.1021/acs.macromol.7b01204 Macromolecules 2017, 50, 6702−6709