A framework for bioinspired search

Suppose a material can be in some configuration, x, and that there is a set of parameters, ${{\rm{\lambda }}_{i,}}$ with which we can tune the material’s properties, averaged over all configurations. A bioinspired evolutionary materials design process then starts with a number of trial designs that sample the possible configurations, x, and proceeds by changing the ${{\rm{\lambda }}_i}$ in an iterative manner such that those configurations become selected whose properties are closest to the design target. We can express this by:

(1)$${\rm{\dot \rho }} = {\rm{\rho }}\left[ {f\left( x \right) - \left\langle {f(x)} \right\rangle } \right].$$

Here, ${\rm{\rho }}\left( {x|{{\rm{\lambda }}_i}} \right)$ is the chance of finding the material in configuration x and $f\left( x \right)$ is a fitness function that measures how close a configuration x is to the desired target. The dot over a variable indicates a time derivative and angle brackets indicate the average over configurations. With each time step or iteration cycle (i.e., with each generation of new trials), configurations that are better (fitter) than the configuration average $\left\langle {f(x)} \right\rangle$ increase in likelihood, while those that are worse decrease. In evolutionary biology and game theory, Equation 1 is the replicator equation.

In the form written, in Equation 1, the dependence on control parameters ${{\rm{\lambda }}_{i,}}$ is implicit. However, the time dependence in ${\rm{\rho }}\left( {x|{{\rm{\lambda }}_i}} \right)$ enters through the time dependence in ${{\rm{\lambda }}_i}$. Using this fact we have:

(2)$${\rm{\partial}} _{{\rm{\lambda }}_i } log\left( \rm{\rho} \right){\rm{\dot \lambda }}_i = \left[ {f\left( x \right) - \left\langle {f\left( x \right)} \right\rangle } \right].$$

Here, ${\partial _{{{\rm{\lambda }}_i}}}$ indicates a partial derivative with respect to ${{\rm{\lambda }}_i}$. If this is multiplied by ${\rm{\rho }}{\partial _{{{\rm{\lambda }}_j}}}\log \left( {\rm{\rho }} \right)$ and integrated over configuration space, we arrive at:

(3)$$\left\langle {\rm{\partial}} _{{\rm{\lambda }}_j } log \left( {\rm{\rho }} \right){\rm{\partial}} _{{\rm{\lambda }}_i } log \left( {\rm{\rho }} \right)} \right\rangle {\rm{\dot \lambda }}_i = \left\langle ({\rm{\partial}} _{{\rm{\lambda }}_j } log \left( {\rm{\rho }} \right)\left[ {f\left( x \right) - \left\langle {f\left( x \right)} \right\rangle } \right]} \right\rangle .$$

The first term on the left-hand side is the Fisher information metric. Inverting it leads to:

(4)$${\rm{\dot \lambda }}_i = \left\langle {\partial _{{\rm{\lambda }}_j } log\left( {\rm{\rho }} \right)\partial _{{\rm{\lambda }}_i } log\left( {\rm{\rho }} \right)} \right\rangle ^{ - 1} \left\langle {\partial _{{\rm{\lambda }}_j } log\left( {\rm{\rho }} \right)\left[ {f\left( x \right) - \left\langle {f\left( x \right)} \right\rangle } \right]} \right\rangle {\rm{.}}$$

This equation provides an explicit prescription for the evolution of the parameters, ${{\rm{\lambda }}_i}$. ^{Reference Miskin13} It propagates the set of ${{\rm{\lambda }}_i}$ into the direction of fitness increase, while simultaneously, the Fisher information metric serves to constrain the associated changes in entropy (i.e., it prevents large changes in the diversity between successive time steps so as to prevent getting stuck prematurely). In the optimization literature, adding this constraint is known as natural gradient descent.

Note that in this way of treating the evolution of configurations, each parameter ${{\rm{\lambda }}_i}$ is a deterministic variable following a path prescribed by Equation 4. The random variations enter at the level of the different configurations, x, which are sampled by the population of trial designs at every step. The algorithm collects information from random samples in configuration space to direct an optimization in parameter space.

For a specified distribution, ${\rm{\rho }},$ the right-hand side depends only on configuration averages. Equation 4, therefore, acts as a ready-made algorithm—sampling ${\rm{\rho }}$ provides the information for evaluating the right-hand side of Equation 4, which then can be integrated to update the parameters ${{\rm{\lambda }}_i}$ for the next iteration step.^{Reference Miskin, Khaira, de Pablo and Jaeger11}

The particular advantage of this formalism is that the search process is able to exploit information about the distribution of material configurations ${\rm{\rho }}\left( {x|{{\rm{\lambda }}_i}} \right)$. In many situations, this distribution is known a priori. Furthermore, with distributions ${\rm{\rho }}$ that are exponential functions, the right-hand side of Equation 4 simplifies and reduces to covariances. For example, if ${\rm{\rho }}\left( {x|{{\rm{\lambda }}_i}} \right)$ is given by the canonical ensemble ${e^{ - {{\rm{\lambda }}_i}{h_i}\left( x \right)}}$, where the ${h_i}\left( x \right)$ are components of an energy functional, executing the derivatives in Equation 4 gives the parameter evolution:^{Reference Miskin, Khaira, de Pablo and Jaeger11}

(5)$${{\rm{\dot{\lambda }}}_i} = - Cov{\left[ {{h_i}\left( x \right),\,{h_j}\left( x \right)} \right]^{ - 1}}Cov\left[ {{h_j}\left( x \right),\,f\left( x \right)} \right].$$

Here, Cov stands for the covariance matrix. Alternatively, one could simply assume a multivariate Gaussian distribution for ${\rm{\rho }}$, without reference to any problem-specific configuration statistics. In this implementation, the right-hand side can still be expressed in terms of covariances and a version of Equation 4 leads to update rules for the parameters defining the Gaussian, namely, the means of the different variables and the elements of the covariance matrix. Structurally, this becomes equivalent to the differential equations describing CMA-ES.^{Reference Miskin13}

At first glance, it might seem that the assumption of a multivariate Gaussian distribution for ${\rm{\rho }}$ would mitigate many of the advantages of evolution strategies over more conventional approaches such as adaptive simulated annealing, in particular, because the latter makes explicit use of the statistical mechanics governing the distribution of configurations. However, for reasons of computational efficiency, CMA-ES remains highly competitive (see Reference Reference Miskin13), and has been shown to outperform other algorithms across a wide range of different optimization problems.

Consider, for example, the problem of finding the most compact configuration of four spheres of equal size (i.e., the configuration that minimizes the sum of the distances to the center of mass). For objects made with less than 13 equal spheres, the correct solutions are deltahedra (i.e., polyhedra with triangular faces).^{Reference Sloane, Hardin, Duff and Conway14} Finding these solutions is remarkably hard for search algorithms and provides a good benchmark test. In Figure 1, we show this for the problem of four spheres, comparing the performance of CMA-ES against a Monte Carlo simulation and adaptive simulated annealing. In the following sections, we illustrate the power of evolution strategies for materials design with two examples, one based on standard CMA-ES and the other using the full information contained in Equation 4 for generating improved, problem-adapted design engines.

Figure 1.

Comparison of algorithms for the task of finding the most compact configuration of four equal spheres. Plotted is relative distance of the best solution f _min from the optimal solution f _opt, as achieved after a given number of function calls, for three different algorithms: covariance matrix adaptation evolution strategy (CMA-ES), random search using a Monte Carlo algorithm, and simulated annealing. For CMA-ES the performance does not require tuning—since the CMA-ES adjusts its covariance matrix, and thus rate of convergence, automatically, we do not need to worry about whether parameters were initially set correctly. Conversely, parameters had to be carefully tuned to achieve acceptable performance from simulated annealing.^{Reference Miskin and Jaeger10}

Granular materials optimization and discovery

Granular materials are inherently far from equilibrium in the (statistical mechanics) sense that, left on their own, they are unable to reach the lowest energy configuration of the constituent grains and instead get trapped, glass-like, in one of a multitude of metastable configurations. If each of these metastable configurations corresponds to a certain value of a measured property of the granular material as a whole, how can one identify those particle-level parameters that will produce a particular, desired value? This constitutes a materials design problem that is particularly challenging, because we cannot limit the search to the subset of lowest-energy “ground-state” configurations.

As a concrete example, consider the case of pouring, under gravity, a large number of identical particles into a box so as to generate an amorphous, glassy aggregate. Tapping the box will then densify this aggregate until some limit is reached asymptotically. A long-standing question has been how this asymptotic packing density depends on particle shape and, in particular, which shape will produce the densest random packing. In the past, the standard approach has been trial and error (i.e., guessing a shape and evaluating the packing properties). The problem has resisted conventional search methods because of its complexity and the ruggedness of the search landscape. Particle shape allows, of course, for endless variations and any search therefore requires some a priori restrictions in terms of the types of solutions to be considered.

Suppose we limit the search to the family of shapes that can be generated by gluing together N spheres, each of variable radius. The only rule here is that all N spheres form one rigid compound particle (“granular molecule”). As the number of degrees of freedom increases with N, the options for placing adjacent spheres within each molecule grow exponentially. Already for relatively small N, say 5–10, the previously discussed build rule enables a myriad of different compact, extended, or even branched granular molecules. As Figure 1 already hints, CMA-ES turns out to be highly efficient at parsing through this large number of options and finding a shape that brings the aggregate properties close to the desired target. Figure 2 shows this for a situation where the optimizer was given N = 10 spheres per particle. After a few hundred iterations and, along the way, setting seven of the radii to essentially zero, it arrived at a trimer molecule where two smaller spheres are positioned roughly 75° apart on a larger central sphere (upper left of top panel). During this shape evolution toward the final trimer molecule, the algorithm discovered and exploited the size distribution of voids in a random packing of spheres, as well as the spatial correlations between voids, such that the two “ears” can fit into the interstices and thereby increase the overall packing fraction to ∼68%, compared to ∼64% for spheres.^{Reference Miskin and Jaeger15} The top panel of Figure 1 shows shapes and packing densities of trial solutions sampled during the evolution toward the densest-packing (trimer) and the loosest-packing molecule (rod in lower right).

Figure 2.

Optimizing the packing density of granular material comprised of particles whose shape can be varied. Each particle consists of 10 constituent spheres with variable radii. (a) Spheres bonded without overlap. In searching for maximum and minimum packing density, the algorithm explored a wide range of different shapes labeled here along the horizontal axis by a shape index that is the sum of all i = 10 sphere radii r _i in the particle (blue crosses). The two shapes emerging after optimization for the density extrema are indicated by arrows.^{Reference Miskin and Jaeger15} (b) Finding the maximum density for particles comprised of bonded spheres with variable overlap. Data points show the performance of the median of the population for each generation, with representative shapes indicated; shaded area indicates one standard deviation. Inset shows an aggregate packing of the optimized trimers.^{Reference Roth and Jaeger16}

Different types of granular molecules become possible once the parameter space is opened up to allow, in addition, for overlap between the spheres constituting an individual particle. This enables the optimizer to explore shapes with surfaces that are less corrugated. This same task now leads to an oblate granular molecule with an essentially equilateral triangular shape. Remarkably, under pouring and tapping this simple particle achieves a random packing density of ∼73%, almost reaching the 74% for a perfect crystal made from equal spheres.^{Reference Roth and Jaeger16}

We have so far discussed searches for new material properties, where the design target is a particular desired property, such as packing density or material stiffness, under fixed boundary or loading conditions. However, the same evolution-based strategies can also be applied when the design target is a particular, desired functional dependence under load. Furthermore, the same search process can be used not only to optimize known materials properties, such as finding shapes that pack more densely, but can be tasked with discovering novel behavior not previously known to occur. For example, all ordinary granular materials held under constant confining pressure weaken when compressed. Using CMA-ES, a specifically shaped small granular molecule was discovered^{Reference Miskin and Jaeger10} that generates the opposite behavior: the aggregate becomes more resistive with increasing compression (i.e., it exhibits strain stiffening).

From statistical mechanics to automated design

As we noted earlier, black box algorithms do not exploit the information contained in the statistical distribution of microstates or configurations. Methods such as CMA-ES and adaptive simulated annealing deal with the actual configurations only through the objective function, effectively condensing a vast amount of information into one parameter. When this extra information is taken into account via the specific form of ${\rm{\rho }}\left( {x|{{\rm{\lambda }}_i}} \right)$ in Equation 4, significant computational speedup can arise.^{Reference Miskin, Khaira, de Pablo and Jaeger11} For example, we discuss here the problem of directing the self-assembly of diblock-copolymer thin films into an ordered array of stripe domains, where the design target is a particular, desired stripe orientation.

Diblock copolymers consist of two, essentially immiscible polymer blocks that are covalently bonded. During annealing, they spontaneously phase separate into patterns of alternating domains comprised of one or the other block.^{Reference Bates, Hillmyer, Lodge, Bates, Delaney and Fredrickson17} These patterns have tunable (via the molecular weight of each block) domain spacing on the scale of 10–100 nm. This has applications as large-area, self-assembled nanotemplates for semiconductor devices or high-density storage media, where the patterns are transferred into an underlying substrate by exploiting the difference between the two copolymer blocks in their resistance to etching. One way of guiding the domains into following a particular direction along the substrate is to implant (chemical) markers that attract one of the two blocks so that it registers spatially with the marker, inducing a surrounding, local patch of the domain pattern to register as well.

The materials design problem now becomes the following: How should the interaction energies between the two blocks and also the interactions with the markers be tuned such that quick and reliable alignment is achieved even if the markers are spaced far apart? In a typical optimization of this type, more than a dozen material and geometry parameters are being tuned.^{Reference Qin, Khaira, Su, Garner, Miskin, Jaeger and de Pablo8,Reference Khaira, Qin, Garner, Xiong, Wan, Ruiz, Jaeger, Nealey and de Pablo9,Reference Jaeger and de Pablo18}

In Figure 3, we show the simulated evolution during annealing from initially disordered domains to a final stripe configuration that is highly ordered and aligned with the markers in the substrate. In the top panel, we plot the performance of black box CMA-ES for a case where the markers were placed a distance corresponding to three times the inherent domain spacing of the copolymer blocks. In the bottom panel, we contrast this with results from an algorithm directly based on Equation 4. Here, ${\rm{\rho }}\left( {x|{{\rm{\lambda }}_i}} \right)$ is implemented as canonical ensemble and the ${{\rm{\lambda }}_i}$ take on the role of setting the strength, relative to the thermal energy kT (where k is the Boltzmann constant and T is temperature), of the various interaction energies in the copolymer Hamiltonian.^{Reference Miskin, Khaira, de Pablo and Jaeger11} With this ${\rm{\rho }}\left( {x|{{\rm{\lambda }}_i}} \right),$ the algorithm becomes significantly more powerful. While standard, black box CMA-ES took 75–100 generations to converge, with the improved algorithm convergence is reached after only 15–20 generations, even when applied to the more challenging problem where the markers are placed a distance of six domain spacings apart.

Figure 3.

Optimizing directed self-assembly of diblock-copolymer thin films. The design goal is to anneal the copolymer into an array of regularly spaced copolymer domains comprised of alternating blocks (shown in red and blue). The insets show snapshots of the domain pattern during the evolution. The data show the evolution of the objective function during optimization as a function of the number of iterations, or generations. (a) Standard, black box covariance matrix adaptation evolution strategy (CMA-ES) with substrate markers, indicated by slight vertical protrusions in the gray substrate and placed a distance of three domain repeat spacings apart.^{Reference Jaeger and de Pablo18} (b) Advanced evolution strategy using problem-specific ${\rm{\rho }}\left( {x|{{\rm{\lambda }}_i}} \right)$ in Equation 4, with substrate markers of adjustable width w placed a distance 6 m of six domain repeat spacings apart. Inset: Composition profile ${\rm{\psi }}$ after optimization, demonstrating the nanoscale phase segregation of the blocks (${\rm{\psi }}$ here indicates the fraction of blue block at different positions along the substrate) as well as the desired alignment with the markers.^{Reference Miskin, Khaira, de Pablo and Jaeger11}

Modularity and adaptability

Biologists have long observed that biological systems are sometimes modular—that is, made of clearly delineated modules that each solve a subtask and whose outputs are put together to perform a more complex task. For example, different cell types in an embryo use similar signal amplification pathways connected in different ways to different cell-specific receptors to generate a large variety of input–output relationships.^{Reference Kashtan and Alon20} Similarly, a great diversity of proteins can be built from the same basic protein domains, each of which might perform a task such as binding and allosteric regulation.^{Reference Bhattacharyya, Remenyi, Yeh and Lim23}

Evolution in time-varying environments is thought to be a significant factor in the emergence of modularity.^{Reference Kashtan and Alon20,Reference Lipson, Pollack and Suh24,Reference Parter, Kashtan and Alon25} Time-varying environments favor modular designs because such designs are quickly adapatable to new environments since new function only requires changing the relationship between existing modules, with no changes inside each module. In an influential work, Kashtan et al.^{Reference Kashtan and Alon20} performed in silico evolutionary simulations to evolve circuits of logical gates in order to compute a particular Boolean function, f. The resulting circuit, shown in Figure 4, is not modular in any sense. They then simulated evolution with a target that was periodically switched between the function, f, and another function, g. As seen in Figure 4, this procedure also results in a circuit equally competent at computing f, but dramatically different in structure—the circuit is highly modular unlike the nonmodular circuit evolved in a constant environment.

Figure 4.

Varying goals during design parameter optimization can select for highly adaptable materials whose function is switched with minimal changes in design parameters. (a) Traditionally, design optimization algorithms are run with a fixed goal until completion. However, biology inspires an alternative approach where the goal is switched between two or more distinct goals as a function of algorithm time, even if optimization is not completed in early iterations. (b) Such time-varying goals are thought to discover parameter sets for each goal that are close to each other in design parameter space. In contrast, fixed-goal evolution typically identifies equivalent parameters for the two goals that are far from each other. (c) A computer evolution study found benefits of such time-varying environments in circuit design. Fixed environments led to a nonmodular complex circuit to compute a given function, f. Varying the goal between two functions, f and g, led to modular circuits that could easily switch between the two functions through a few small changes. (c) Reprinted with permission from Reference Reference Kashtan and Alon20. © 2005 National Academy of Sciences.

Similar results have been seen in other biological systems,^{Reference Kussell, Leibler and Grosberg26,Reference Kussell and Leibler27} with a recent promising example in immunology.^{Reference Wang, Mata-Fink, Kriegsman, Hanson, Irvine, Eisen, Burton, Wittrup, Kardar and Chakraborty28} The human immune system creates antibodies that can specifically detect any particular foreign protein (“antigen”) through a rapid evolutionary process. However, to deal with rapidly mutating viruses such as HIV, the body must evolve “generalist” antibodies that can detect structural motifs common to several multiple antigens instead of being highly specific to just one of them. Recent work^{Reference Wang, Mata-Fink, Kriegsman, Hanson, Irvine, Eisen, Burton, Wittrup, Kardar and Chakraborty28} has shown that using a temporal sequence of vaccination that cycles the antigens shown can lead to generalist antibodies (also called broadly neutralizing antibodies), even if such generalists are rare and on average, less fit, than the specialists. Thus, evolution in time-varying environments is thought to naturally lead to adaptable and often modular designs in diverse biological systems.

Time-varying design pressure for materials

These biological works suggest that time-varying goals in evolutionary algorithms can help identify design parameters that correspond to multifunctionality or adaptability, but are hard to find by static protocols. However, discovering multifunctional design parameters requires reconceptualizing a parameter search in design and evolutionary algorithms as a dynamic process. In this language, traditional design approaches correspond to equilibrium dynamics in parameter space. In contrast, the biological examples previously discussed point at an alternative strategy—we can exploit nonequilibrium effects during parameter search to dynamically stabilize adaptable and multifunctional parameter choices.

The biological examples point to some clear lessons for what kinds of materials problems will benefit from time-varying design:

1. 1. High-dimensional parameter space. The materials design problem must involve a search over many parameters (e.g., tens or hundreds of parameters). Common reasons for large numbers of materials design parameters are the spatial extent of systems (e.g., elastic properties in different parts of a solid) or the complexity of constituent elements of a material (e.g., many parameters required to describe the shape of a particle).

2. 2. Degeneracy of parameter choices. Many different parameter choices should result in the same material property. Time-varying environments exploit this degeneracy to pick out parameter sets that are adaptable or multifunctional. In contrast, if a given material property is achieved by a unique choice of parameters, such optimization for multifunctionality is not possible.

3. 3. Rugged optimization landscapes. If the design problem is convex or otherwise simple, multifunctional solutions might be obtained through analytic methods or numeric techniques such as simultaneous optimization of two different target functions. However, such approaches fail when the optimization landscape is high-dimensional and rugged with numerous local minima. In these circumstances, time-varying design pressure presents a broadly applicable approach that does not require detailed knowledge of the particular design problem.

4. 4. Existing iterative design algorithms for single target behaviors. The algorithms in question can be diverse in nature (e.g., gradient descent methods or evolution strategies such as CMA-ES). The crucial requirement is that the design procedure is an iterative search process over parameter space that can be interrupted before design is complete, so that the design target can be switched.

Self-folding polymer and self-assembly design

A classic design problem that satisfies the previously discussed requirements is that of self-folding polymers. The self-folding polymer design problem targets a specific, folded three-dimensional structure through the choice of polymer sequence.^{Reference Karplus29} Such problems cover natural heteropolymers such as proteins and RNA as well as recent synthetic equivalents such as DNA origami in which a long DNA scaffold strand is directed to fold into specific shapes using short helper strands.^{Reference Rothemund30} For most heteropolymers, multiple sequences can fold into a given structure. Such degeneracy has been investigated in detail for proteins^{Reference England and Shakhnovich31} and RNA.^{Reference Schuster, Fontana, Stadler and Hofacker32,Reference Fontana, Stadler, Bornberg-Bauer, Griesmacher, Hofacker, Tacker, Tarazona, Weinberger and Schuster33} Further, the design problem is known to face a complex rugged landscape in the space of sequences—nearby sequences can vary greatly in their ability to fold into a specific structure.

Thus, multifunctional polymers (i.e., polymers capable of folding into different folded structures with minimal or even no change of sequence) are naturally suited for designing with time-varying pressure. Several research groups have created particular examples of such multifunctional polymers, often by exploiting deep knowledge about the particular kind of polymer used and the particular desired structures. For example, DNA origami usually targets the assembly of one folded structure by combining a long scaffold DNA strand with many short “staple” DNA strands; the staple strands bring specific segments of the scaffold together, thus folding the scaffold into the desired structure. In recent work,^{Reference Dunn, Dannenberg, Ouldridge, Kwiatkowska, Turberfield and Bath34} one scaffold strand was programmed to controllably fold into two distinct structures by changing staples involved. However, the design method did not attempt to minimize the number of staple changes required to switch states. While design algorithms to identify staples for a given structure exist, there are no systematic methods to optimize the simultaneous design of two or more distinct structures. These lessons from biology suggest that such new design algorithms might be unnecessary; simply running existing single-target iterative algorithms with time-varying target structures can produce multifunctional designs.

Similar multitarget polymer folding ideas have been explored earlier in the context of proteins^{Reference Abkevich, Gutin and Shakhnovich35–Reference Fink and Ball37} and RNA.^{Reference Schultes and Bartel38} In a particularly striking demonstration,^{Reference Schultes and Bartel38} two different RNA structures with distinct functions were identified. The two folded structures, shown in Figure 5, occur using two completely unrelated RNA sequences. However, Schultes et al.^{Reference Schultes and Bartel38} exploited the many-to-one nature of the RNA sequence-to-structure map to identify a sequence that could switch between the two folded structures upon a single mutation.

Figure 5.

Examples of adaptable materials that can switch functionality with minimal change in design parameters. (a) A single DNA origami scaffold programmed to kinetically fold into one of two distinct structures. Adapted with permission from Reference Reference Dunn, Dannenberg, Ouldridge, Kwiatkowska, Turberfield and Bath34. © 2015 Springer. (b) Two distinct RNA structures with distinct functions of joining (“ligating”) RNA strands or cleaving an RNA strand. While the natural sequences for these two functions are widely different, Reference Reference Schultes and Bartel38 identified two other sequences, differing by only one mutation, that fold into the two specified functional structures. Adapted with permission from Reference Reference Schultes and Bartel38. © 2000 AAAS. (c) Self-assembling particles with interactions that can kinetically switch assembly between different structures. While self-assembly interactions usually target one specific structure, Murugan et al.^{Reference Zhong, Schwab and Murugan40} showed that some choices of interaction parameters allow multiple targets without interference. (d) Lattice model of protein-inspired allosteric material, evolved to communicate a mechanical signals from one side to the other. When in silico evolution was carried out with ligands that bind to the green and blue binding sites in a time-varying manner, the modular and adaptable solution on the right was obtained. Adapted with permission from Reference Reference Hemery and Rivoire52. © 2015 American Physical Society.

Thus, while multifunctional heteropolymers have been successfully designed in particular cases, few studies have identified systematic design principles. One approach to systematic design used simultaneous optimization^{Reference Abkevich, Gutin and Shakhnovich35,Reference Abkevich, Gutin and Shakhnovich36} of multiple target functions with associated difficulties of a rugged landscape. Time-varying design presents a promising alternative to produce multifunctional heteropolymer sequences by exploiting existing algorithms for single structure design.

Inspired by the ability of biology to reuse components to build different structures, recent works^{Reference Murugan, Zeravcic, Brenner and Leibler39,Reference Zhong, Schwab and Murugan40} have studied the possibility of a “multifarious” soup that reuses components to achieve different structures. However, this approach requires specific interactions between many species of particles, which is currently possible only using DNA hybridization-based interactions.^{Reference Ke, Ong, Shih and Yin41} However, a systematic way of identifying interaction design parameters for multitarget self-assembly—where the assembled structure can be changed through minimal modification of interactions—is still lacking.

Mechanical metamaterials

Recent work on mechanical metamaterials, ranging from origami to elastic networks, has targeted the design of specific deformation modes.^{Reference Bertoldi, Vitelli, Christensen and van Hecke42} The literature on thin sheets and origami has explored how to program different folding modes in sheets with desired geometry.^{Reference Santangelo43–Reference Pinson, Stern, Carruthers Ferrero, Witten, Chen and Murugan45} Other recent works^{Reference Rocks, Pashine, Bischofberger, Goodrich, Liu and Nagel46,Reference Tlusty, Libchaber and Eckmann47} have designed elastic networks with extended deformation modes in response to locally applied forces, mimicking allosteric deformation in proteins. In fact, some of these works^{Reference Yan, Ravasio, Brito and Wyart48} design such deformations through evolutionary or gradient descent processes.

While much metamaterial design has focused on programming a specific functionality into a material, broader applications demand multifunctional elastic materials that can switch between distinct mechanical properties.^{Reference Hawkes, An, Benbernou, Tanaka, Kim, Demaine, Rus and Wood49} Similar multifunctionality has been studied in the biological literature. Raman et al.^{Reference Raman, White and Ranganathan50} and Pumir et al.^{Reference Pumir, Graves, Ranganathan and Shraiman51} determined the structure of a protein before and after mutations and showed how a few mutations can change allosteric modes. In related theoretical work, Hemery et al.^{Reference Hemery and Rivoire52} found that evolution under time-varying selection for different allosteric modes has unique structural fingerprints in a simple lattice model of proteins (e.g., allosteric motion is communicated by a narrow neck whose size depends on the frequency of switching). Such a neck was related to the ability to switch between different deformation modes in response to a few mutations.

Combining ideas from metamaterials design and the evolution of mechanical properties of proteins along the lines previously discussed is likely to offer several benefits. For example, evolution with time-varying requirements can provide a systematic way to create elastic networks whose response to forces can be dramatically changed by a few small local changes of stiffness. Further, this process of designing elastic materials with novel mechanical properties can provide insights into the functioning and evolution of natural proteins themselves.

Shape-sensitive materials

Several recent developments have shown how the ability to modulate the shape of building blocks can control a material’s overall properties. As previously discussed (see Figure 2), changing the shape of constituent particles in a granular material can change the packing fraction over a large range. Similarly, shape can greatly modulate the elastic moduli of packings.^{Reference Miskin and Jaeger10,Reference Athanassiadis, Miskin, Kaplan, Rodenberg, Lee, Merritt, Brown, Amend, Lipson and Jaeger53–Reference Dumont, Houze, Rambach, Salez, Patinet and Damman56} Entropic interactions between shapes, such as depletion interactions or configurational entropy,^{Reference Huntley, Murugan and Brenner57–Reference Young, Personick, Engel, Damasceno, Barnaby, Bleher, Li, Glotzer, Lee and Mirkin59} have been shown to lead to specific interactions that can be exploited for self-assembly and other such applications.

All of these shaped-based approaches involved searching over the high-dimensional space of shape parameterization to achieve desired material properties. Even restricted searches often involve several dozen parameter descriptions of shape.^{Reference Miskin and Jaeger15} Further, shape spaces are often degenerate for most physical properties (e.g., many distinct shapes for constituent particles result in [essentially] the same packing fraction). Finally, any physical property such as packing fraction varies in a rugged manner across the space of shapes, changing dramatically for some shape changes, but not for others. The design of shapes for multifunctional materials satisfies the criteria and seems well suited for design through time-varying goals. Through this approach, one then could create particles that switch their packing propensity between two designated values, ${\phi _A}$ and ${\phi _B},$ simply by morphing (e.g., via swelling and shrinking) between two shapes A and B with minimal change in volume. Alternatively, one might be able to create granular materials that switch between two or more desired bulk elastic properties through small changes in the shape of the constituent particles.