Emergent disorder and mechanical memory in periodic metamaterials

Ordered mechanical systems typically have one or only a few stable rest configurations, and hence are not considered useful for encoding memory. Multistable and history-dependent responses usually emerge from quenched disorder, for example in amorphous solids or crumpled sheets. In contrast, due to geometric frustration, periodic magnetic systems can create their own disorder and espouse an extensive manifold of quasi-degenerate configurations. Inspired by the topological structure of frustrated artificial spin ices, we introduce an approach to design ordered, periodic mechanical metamaterials that exhibit an extensive set of spatially disordered states. While our design exploits the correspondence between frustration in magnetism and incompatibility in meta-mechanics, our mechanical systems encompass continuous degrees of freedom, and thus generalize their magnetic counterparts. We show how such systems exhibit non-Abelian and history-dependent responses, as their state can depend on the order in which external manipulations were applied. We demonstrate how this richness of the dynamics enables to recognize, from a static measurement of the final state, the sequence of operations that an extended system underwent. Thus, multistability and potential to perform computation emerge from geometric frustration in ordered mechanical lattices that create their own disorder.


I. INTRODUCTION
When interactions between different components of a complex system cannot simultaneously be minimized [1,2], the lowest-energy state is a compromise in which some elements of the system are left "unhappy".This frustration can lead to constrained disorder, degeneracy, and multistability [3][4][5].Recent research translated these concepts from magnetic spin systems to engineered soft-matter systems, such as acoustic channels [6], buckled elastic beams [7], and monolayers of colloidal spheres [8,9].Here, frustration in magnetic spin systems is related to incompatibility of soft modes in mechanical systems [10][11][12][13].Mapping of spin systems to metamaterial architectures has provided insight into the mechanical consequences of frustration, such as stress control [14] and domain-wall topology [15,16].In particular, exploring irreversibility and history dependence through these analogies has opened new routes for programmable elastic responses and mechanical memory storage [17][18][19][20][21][22][23][24], and has set the grounds for further advances in mechanical computing [25,26].However, multistability and complex memory formation are typically traits of disordered and amorphous systems [21,27,28], which are harder to predict and control.
On the other hand, periodic mechanical systems are generally not multistable, as their long range elastic interactions tend to resolve frustration with long-range ordered ground states [7,29].Here, the displacement of mechanical degrees of freedom can take intermediate values, leading to an ordered compromise and lifting the degeneracy associated with frustrated spin systems.
In this Article, we introduce ordered mechanical metamaterials that exhibit a large multiplicity of disordered metastable states.This leads to multistability of internal degrees of freedom and to mechanical memory, as configurations can depend on their preparation history.Crucially, our approach relies on spatial separation between the frustrated motifs in the metamaterial, which we achieve via a mapping to vertex-frustrated [30] artificial spin ice (ASI) [31].As a result, longrange cooperative effects are suppressed, and disorder emerges from local rules.Furthermore, the graph of transitions between the metamaterial's states [27] contains irreversible pathways, and its structure leads to non-Abelian behavior.Namely, the system's state depends not only on the external manipulations applied on the metamaterial, but also on their precise sequence.We utilize this to demonstrate how the history of operations acted on the system, and their precise sequence, may be inferred from a static measurement of the system's final state.This takes an important step towards a systematic implementation of computation in mechanical metamaterials [32].
The design of our mechanical metamaterial is inspired by the frustration-based designs of urations often captured by interesting emergent descriptions [35,36].There are strong similarities between ASI and mechanical metamaterials realized by repeated arrangements of simple units endowed with a soft mode: In ASI, nano-islands are arranged along the edges of a lattice, and their magnetization is described by binary arrows pointing toward or from the vertices of the lattice, with certain vertex configurations locally minimizing the energy.Identically, in a mechanical metamaterial, displacements point into or out of the repeating units, and the softest deformation mode of each unit is related to a certain mutual arrangement of these displacements.
The similarity extends to frustration and incompatibility: In vertex-frustrated ASI lattices, not all vertices can simultaneously be in their lowest-energy state [30,37], as demonstrated for the Shakti-lattice ASI in Fig. 1(a).In such lattices, when going along a loop of vertices around any plaquette of the lattice, the spins may not be assigned in a way such that all vertices will be in their energetically-preferred configurations.Equivalently, in combinatorial mechanical metamaterials, one can arrange the units so that they cannot all simultaneously deform according to their soft mode [10,12,13].Vertex frustration in magnetism thus corresponds to incompatibility in mechanical metamaterials.The plethora of frustrated ASI geometries and topologies [30,34,[38][39][40] inspires novel metamaterial designs.We show that in the mechanical system such ASI-inspired design can suppress long-range ordered ground states, and lead to emergent disorder and complex memory formation.

II. THE CHACO MECHANICAL METAMATERIAL
In the Shakti ASI [30,35,38,41,42] In ASI, nanoislands are naturally magnetized.In metamaterials, one can induce an equivalent spontaneous displacement of all edges by prestressing the system [16].Namely, pinning an edge to a distance 2α , with α < 1, imposes a compression, which causes it to buckle by an amount This results in a contraction by a factor α 2 of the two-dimensional lattice.In the prestressed Chaco lattice, the units may not all simultaneously adopt their zero-energy deformations, due to the inherent incompatibility.Nevertheless, in the weak coupling limit, we expect the free nodes to behave in an approximately binary way [16], adopting displacements s i ≈ ±δẑ i from their rest position, where ẑi is the local direction perpendicular to the edge.In this limit, the k 1 bonds are approximately relaxed and most of the strain is concentrated on the k 2 bonds.As we discuss below, the overall compression applied to the lattice breaks the symmetry between stretching and compressing, such that compressing a k 2 bond costs more energy than stretching it.
Experimentally, we create the elastic network shown in Fig. 1

III. MECHANICAL EQUILIBRIUM STATES
The overall compression breaks the Z 2 spin-reversal symmetry between nearest-neighbor degrees of freedom.The mechanical degrees of freedom behave in a spin-like binary way for k 2 k 1 , whereas elastic, continuous deformations are expected otherwise.Consider the interaction energies when ideal displacements ±δ are imposed on pairs of neighboring edges: One in and one out displacement allow the internal k 2 bonds to rotate and remain in their relaxed lengths l 0 = √ 2 and l ∆ 0 = .Two out displacements result in a stretched internal k 2 bond, with lengths l + and l ∆ + , and energies E + and E ∆ + .Two in displacements result in a compressed internal k 2 bond, with lengths l − and l ∆ − , and energies E − and E ∆ − .Due to the overall compression, stretching the internal bonds costs less energy than compressing them further, as shown in Fig. 2.An important consequence of this breaking of the Z 2 symmetry is that excited triangular and square mechanical units have distinct energy hierarchies than the z = 3 and z = 4 Shakti ASI vertices [30].In the ideal binary limit, k 2 k 1 , the energy of a mechanical unit follows from the number of stretched and compressed bonds.There are six energy levels for the square units with energies: 0, 2E + , 2: Z 2 symmetry breaking -Energetics of the ideal nearest-neighbor interactions for the square and triangular mechanical units.Ideal displacements allow the k 1 bonds to be relaxed and place all energy on the internal k 2 bonds.Spin inversion symmetry is broken in the mechanical system: bond extension costs less energy than bond compression.
Similarly to the magnetic Shakti, the most energy-efficient way to resolve the frustration in the mechanical Chaco is to localize stress on half of the triangles.This leaves the remaining triangles and all squares close to their zero-energy states.However, the continuous nature of the mechanical degrees of freedom allows minimizing the energy by long-range pairing of the excited triangles.
Namely, within two adjacent excited triangles, the k 2 bonds do not have to take the ideal values l ∆ + and l ∆ − mentioned above.Instead, due to the nonlinear dependence of energy on displacement, there is an intermediate balance between their deformations which minimizes the energy.In a perfect realization, this leads to an ordered ground state with defect pairing [Fig.3(a)], which has a four-fold degeneracy, since stresses may be placed on either vertical or horizontal pairs of triangles and since spin reversal leaves the energy unchanged.For sufficiently large k 2 /k 1 , the ground state is the only mechanically-stable configuration.However, for small k 2 /k 1 , we also find an extensive manifold of metastable states.Namely, any configuration with all squares in one of their individual zero-energy states becomes metastable.This holds for configurations with stresses localized to exactly half of the triangles, equivalent to any ground state of the Shakti [Fig.3(b)].It also holds for configurations where more than half of the triangles are excited [Fig.3(c)].Although these states have higher energy than the ground state, the excess stress in the triangles does not suffice to overcome the energy barrier for flipping the squares.Thus, since thermal fluctuations are negligible, if the system is in such a disordered metastable state, it will remain there and will not relax to the ordered ground state.See Appendix B for computational details.

IV. MEMORY AND IRREVERSIBLE DYNAMICS
We distinguish between two types of degrees of freedom within the Chaco metamaterialthe square units, and the central edges between pairs of back-to-back triangles.For small k 2 /k 1 , we describe both of them as approximately binary variables, stuck in one of their stable states.
Namely, squares adopt one of their soft configurations, and each central edge buckles to one of two directions.The latter degrees of freedom constitute hysterons [18,43,44], bistable elements whose state is history dependent.Hence, they are a natural representation for memory in the metamaterial.As described above, both states of the squares are always stable.In contrast, the stability of the central edge between two triangles depends on the configuration of the four squares surrounding it.We label the 2 4 = 16 possible states of these squares by 2 × 2 matrices with binary entries, where 0 indicates a displacement into the double triangle and 1 out of it.By symmetry, these 16 states reduce to the seven distinct configurations shown in Fig. 4(a).
We now fix all squares to their zero-energy states, and compute the minimal-energy states of the central edge as function of k 2 /k 1 , for each configuration.The horizontal displacement of the central point is negligible, thus we focus on its vertical displacement y.We find that it may exhibit monostable, metastable or bistable behavior for the different configurations and for different values of k 2 /k 1 , as shown in Fig. 4(b,c): Due to their vertical asymmetry, configurations i − iii are typically monostable, with y > 0: configuration i allows both triangles to adopt zeroenergy states, such that y = δ for any k 2 /k 1 ; configurations ii, iii show a reduced amplitude 0 < y < δ, which decreases as k 2 /k 1 increases.For very small values of k 2 /k 1 , metastable solutions with y < 0 appear for configurations i − iii.Configurations iv − vii are vertically symmetric.As a result, the central edge is bistable for small k 2 /k 1 , while for larger k 2 /k 1 , the central point adopts the value y = 0. From here on, we focus on 0.022 k 2 /k 1 0.057, where configurations i − iii are strictly monostable and configurations iv − vii are bistable, as shown in Fig. 4(c), middle panel.
Our experimental design is similarly aimed to exhibit such behavior for all seven configurations.However, as opposed to the theoretical spring model, in the experimental system, the squares apply torque to the edges of the central beam.In configuration vi the torques are not symmetric under up-down reflection, and consequentially they excite the second bending mode of the central beam, as shown in Fig. 4a.Hence, this configuration is bistable only in the theoretical springs model, while it is monostable in the experimental network of elastic beams.Since the configuration of the squares governs the stability of the central edge, we can precisely control the states of the triangles by manipulating the squares around them.Interestingly, the final state of the central edge depends not only on the configuration of squares, but on the exact order in which they were flipped, and consequently transition pathways are sequence dependent; The transitions between states are best understood using the directed graph [18,27,44] shown in Fig. 5(a), where we separately mark the two possible states of the central edge for the six bistable configurations, resulting in a total of 22 states.All transitions between states with the same direction of the central edge are reversible (gray).An irreversible transition (peach) occurs whenever a bistable configuration is connected to a monostable configuration with the opposite direction of the central edge.
We propose a protocol, shown in simulations and experiments in Fig. 5(b,c), respectively, that demonstrates the aforementioned control: We start from the bistable configuration ( 0 1 0 1 ), with the central edge pointing upwards (top panels).Next, we flip the lower-left square, thus switching to the monostable state ( panels).Finally, we flip the lower-left square back to its initial state, which brings us back to the original bistable configuration ( 0 1 0 1 ).Interestingly, the central edge does not revert back, but remains down (bottom panels).After this irreversible path, if we continue flipping the lower-left square, the system reversibly switches back and forth between the latter two states.

V. NON-ABELIAN RESPONSE
The subsystem of 2 × 2 squares analyzed above suffices to give rise to rich dynamics and pathways within the Chaco metamaterial.The generic motif underlying this richness, we find, is path dependency.The state of the Chaco metamaterial is sensitive to the ordering of external operations, rendering it non-Abelian.Recently, such non-Abelian responses were used to realize both finite state machines [45] and Set-Reset latches [20], which are elementary components of computation.To demonstrate this here, we start from the state ( We now consider a target configuration, for instance the one shown in Fig. 7b, and denote it σ t , the state vector of all double triangles.We ask which flip sequence results in this configuration.
The partial ordering between neighboring flips can be translated to a series of inequalities.Solving

VII. DISCUSSION
The mechanism underlying the non-Abelian response sheds light on the generic emergence of sequence-dependent responses in frustrated mechanical systems [19,20,[47][48][49].The order in which bistable degrees of freedom are manipulated follows a pathway of states, which favors one of the possible final states.In essence, the sequence directs the otherwise spontaneous symmetry breaking of the system at its final state.This framework may be used to cleverly manipulate mechanical systems between multiple functionalities [50][51][52], as different states of the system may exhibit different mechanical responses.This, we envision, may allow adaptable, history dependent mechanical responses via in materia information processing [25,26].
We have focused theoretically and experimentally on the weak interaction limit (measured by the coupling parameter k 2 /k 1 ).Strong interactions between hysterons lead to long-range order and a distinct ordered ground state [7], see Appendix C.However, intermediate interaction strengths may give rise to even richer dynamics [18,21].For example, they can lead to long transients, or complex cycles with longer periods [53][54][55], which can realize an array of computational tasks such as counting [56].
Following the vertex-frustration analogy back to the magnetic realm, non-Abelian protocols may increase the memory capacity in storage devices based on frustrated geometries.Due to the sequence sensitivity, manipulating n binary degrees of freedom may encode more than 2 n unique states, as the system's history may influence additional features that are beyond the states of the manipulated degrees of freedom.Thus the manifold of reachable states may be enriched by frustrated interactions.
Finally, the Chaco mechanical metamaterial inherits from the Shakti ASI a topological structure: the possible allocation of its incompatibilities can be mapped both in the free fermion point of a six vertex model [41] and thus into a dimer cover model [42].These height models can be thought of as instances of so-called classical topological order [57] for which boundary conditions strongly constrain the manifold in the bulk, and in the case of Chaco metamaterial, that would imply an encoding via the boundaries.Moreover, using the same dimer mapping and the design strategy by triangle rotation pioneered in ref [12] for the kagome lattice, an extensive number of Chaco-based metamaterials can be produced.
To realize the frustrated Chaco lattice, we design a network of elastic beams, tracing the geometry of the Hookean springs in our theoretical model.Using a 3D printer (Prusa MK3 i3) we manufacture molds for the network, which are then used to cast samples made from silicone rubber (Mold Max™) with a Shore hardness of 30 A. The geometric dimensions of the rubber samples are shown in Fig. 9.The basic beam length, constituting half the edge of each square or triangular unit in the Chaco lattice is equal to = 15mm.To induce prestress in the samples, we constrain the elastic network to dimensions smaller than its rest configuration, using a set of pins, 2 ± 0.01mm in diameter.
Each pin acts as a rotation axis fixed in space.The spacing between neighboring pins is set to induce a uniform compression factor of α = 0.92 ± 0.01.To limit the deformation of the central beam between two triangles to its first buckling mode, we remove from the triangles the elastic beams which trace the k 2 springs.Instead, the region between squares meeting at a vertex of the lattice serve as a coupling mechanism between the deformation of adjacent beams in each triangle.The k 2 springs are not removed from the squares, however, to ensure the squares do not exhibit additional metastable states.The control sequence protocols are performed manually.The beams corresponding to the k 2 springs allow easy manipulation of the squares between their two stable states.We force the squares to flip by pushing one of their outer beams.See supplemental movie [46].Experiments are documented using a digital camera (Sony alpha3).In Figs. 4, 5 and 6, the background is digitally removed for clarity.ization E ∆ /E 0 → 1/2 as k 2 /k 1 → 0. This behavior is also seen for the Shakti-lattice ASI ground state [Fig.3(b)], which implies that the energy is also localized to the triangles, and further implies that the gap between the Shakti ASI ground state and the true mechanical ground state closes as k 2 /k 1 → 0. For larger k 2 /k 1 the triangles and the squares in the Shakti ASI ground state have higher energy compared to the mechanical ground state.The energy of the squares in the Shakti ground state increases smoothly until k 2 /k 1 ≈ 0.15, where sharp jumps start to appear, indicating that some squares flip so that excited triangles can become paired.
For relaxed squares in random orientations [Fig.3(c)], the energy is also localized to the triangles, but the presence of higher-order excitations on the triangles causes the energy stored to be nearly twice that of the triangle energies in the mechanical ground state, giving E ∆ /E 0 → 1 as k 2 /k 1 → 0. Finally, the state with random edge displacements has large energy stored in both the squares and the triangles.As k 2 /k 1 increases, edges become unstable and flip, taking the configuration towards the ordered ground state.The three metastable states shown all reach the ground state for k 2 /k 1 > 0.5.In larger systems, the true ground state is not always fully reached because several domains with different ground-state orientations may form.show the same behaviors for each class.The important major difference that appears in the larger system is some of the energy curves do not collapse to merge with the mechanical ground state curves as k 2 /k 1 increases.Inspection of the real space evolution of the configurations for the larger system in Fig. 11 shows that the failure to collapse completely to the ground-state curves is a result of the presence of multiple competing ground-state domains.

FIG. 1 :
FIG. 1: Design principle -(a) The Shakti lattice artificial spin ice (ASI) exhibits vertex frustration: A conflict (red arrows) results when trying to place all vertices around a rectangular plaquette in their lowest-energy states (as individually shown in b).(c) The Chaco lattice mechanical metamaterial inherits this vertex frustration of the Shakti ASI: The five units surrounding any corner in the metamaterial may not simultaneously deform according to their softest mode (as individually shown in d).Thick lines indicate k 1 springs forming the edges of the square and triangular units.Thin lines indicate k 2 springs responsible for the interactions between neighboring edges within each unit.Spontaneous edge displacements are generated by pinning the corners of all units (white dots) to distances smaller than the relaxed edge length.(e) Experimental realization of the Chaco metamaterial, and (f) the softest deformations of its units.In the triangular units, the different stiffnesses at the corners connecting adjacent edges give rise to the softest modes obtained by the k 2 springs in the theoretical model.
(e), top.Similarly to the spring model, a thin enough elastic beam of length 2 compressed by a factor α tends to buckle via its first deformation mode.We prestress the rubber network by constraining the corners of all square and triangular units to an array of metal pins connected to a rigid substrate, and spaced such that the network is uniformly compressed by factor α [Fig.1(e),bottom].The geometry of the units leads to mechanical soft modes [Fig.1(f)], which reproduce those of the theoretical springs model [Fig.1(d)].See Appendix A for experimental details.

FIG. 3 :
FIG. 3: Multiple stable configurations -Square and triangle unit energies for sample configurations with increasing energy in a 4 × 4 Chaco lattice with periodic boundary conditions, all for k 2 /k 1 = 0.033 and α = 0.9: (a) Ordered ground state with energy localized on vertical pairs of triangles.(b) Complex stress distribution of a metastable state derived from a magnetic Shakti ground state, in which each stressed triangle is adjacent to a relaxed triangle.(c) Metastable configuration with all squares in relaxed random orientations and stress localized to more than half of the triangles.

FIG. 4 :
FIG. 4: Stability of double triangles -(a) The seven distinct configurations of a pair of back-to-back triangles and their surrounding four squares within the Chaco lattice, with 1 (0) indicating displacement out (in) of the triangle.We divide them by the symmetry of forces applied to the central beam; Configurations i-iii are monostable.Configurations iv,v,vii are bistable.Configuration vi is bistable in the theoretical springs model, while experimentally its central beam exhibits a single stable state with a second buckling mode.(b) Stable vertical positions of the central point between the two triangles as a function of k 2 /k 1 for α = 0.9 and assuming that the four external edges are fixed with ideal displacements s i = ±δẑ i .Inset indicates the vertical displacement y.(c) Potential energy E(y) as function of central beam displacement for configurations ii (left) and vi (right) representing the two classes.At high k 2 /k 1 the symmetric configurations become monostable (bottom right).At low k 2 /k 1 the asymmetric configurations develop metastability (top left).In between is our operational regime where all configurations follow the behavior shown in panel (a).

0 1 1 1 FIG. 5 :
FIG. 5: Local memory and irreversibility -(a) Transition graph between the possible states.Bistable configurations are shown twice to include the possible deflections of the central edge.Reversible transitions are plotted in grey.Irreversible transitions, which flip the central edge, are plotted in peach.(b,c) Irreversible control sequence for the internal bistable edge in simulations (b) and experiment (c).The starting configuration (top) is bistable with the internal edge up (cyan arrow).Flipping the lower-left square (middle) creates a monostable configuration, so the internal edge flips down (yellow arrow).The lower-left square is returned to its original bistable state (bottom), but the internal edge remains pointing down.The pathway corresponding to panels (b,c) is highlighted in black on the transition graph (a).

FIG. 6 :
FIG.6: Non-Abelian response -starting from state ( 0 1 1 0 ), consider flipping the bottom squares marked A and B as two operations acting on the system.The response of the Chaco metamaterial to these operations depends on the order in which they are applied, as indicated by the displacement of the central edge.Black dots mark the square that has just flipped.

FIG. 7 :FIG. 8 :
FIG. 7: Sequence Recognition -(a) Initial configuration for the sequence recognition protocol, where the squares marked A − D are flipped; (b) Target configuration after all squares have been flipped.The state of the double triangles, denoted σ 0 translates to a set of inequalities for the flipping order; (c) distance from the target d = i (σ i − σ t i ) 2 along any flip sequence.The sequence BCDA which satisfies all inequalities is the only one reaching the target, and the distance metric allows sequence recognition.

Figure 10 FIG. 10 :
Figure10plots the relative energy curves for two system sizes and also shows three separate

FIG. 11 :
FIG.11: Real space evolution of the metamaterial configurations for a system of 8 × 8 squares, for the three classes of initial condition, unpaired defects corresponding to a magnetic Shakti ground state, random initial choice of squares in their minimum energy orientations, and random choice of the edge displacements.As k 2 /k 1 increases, edges inside the system start to flip to bring the configuration closer to the mechanical ground state.However, for the k 2 /k 1 values shown here, there are still several competing ground-state domains, separated by domain walls.