Designing and Analyzing Multistable Mechanisms Using Quadrilateral Boundary Rigid Origami

Origami patterns that are not rigidly foldable, i.e., patterns that require the panels to undergo elastic deformation, often exhibit multistability. Such multistable structures, e.g., twist fold [1], triangulated cylindrical origami [2,3], and pleated hypar [4,5], exhibit rich nonlinear behaviors of the snapping between states. Multistable structures can be used as the building block for the design of mechanical devices and metamaterials, such as switches with tactile feedback [6], vibration isolator [7], mechanical memory [8], cellular materials [9], and meter-scale deployable structures that stay stable in the deployed state [10,11].

In this paper, we propose a novel method for designing origami-based multistable structures, by which we mean (1) to obtain the prescribed overall motion and (2) to control the stiffness of snapping provided by the elastic strain. (1) Most of the existing design approaches are based on the assembly of known multistable origami structures as the building blocks. As such a building block transforms in a predefined manner, it is difficult to design the desired shape deformation. For the versatility of motion design, we adapt quadrangular boundary rigidly foldable origami proposed by Demaine and Ku [12] as the building block and extend the system for bistability design. (2) The stiffness of snapping is controlled by the strain energy potential barrier between stable states (Fig. 1). The strain energy potential and the resulting force can be computed numerically through the angular spring model [1], bar-and-hinge model [4,5], and finite element method analysis [13–15]. However, the inverse problem, i.e., giving the stiffness first and then calculating the structure that realizes it, remains a challenging task. To tackle this design problem, Melancon et al. [11] use the concept of geometric incompatibility to estimate the strain energy barrier. In our study, we take a similar approach to designing geometric incompatibility. We provide a novel incompatibility metric for a quadrangular boundary module.

We solve the design problem of origami-based multistable structures by first prescribing the desired motion using linkage structures with quadrilateral holes, called frames. Frames are relatively easy to design because they form less constrained systems compared to rigid origami without holes; origami without holes tends to be over-constrained. Then, we fill the frames with origami modules, called quadrilateral boundary modules, so that the motions of the modules approximately follow the motion of the frames provided by the linkage. When fitting modules to the linkage, we further create intentional incompatibility between the motions of frames and quadrilateral boundary modules (Fig. 1). This causes strain on the modules, leading to snapping behavior. We propose an optimization-based approach to provide tunable stiffness of snapping by adjusting the amount of incompatibility.

In this paper, we first introduce linkages with quadrilateral holes that can create versatile motions (Sec. 2). Then, we provide the basic geometry of quadrilateral boundary modules (Sec. 3). We also propose a family of modules with sufficiently wide design space. To evaluate equivalent elastic strain due to the incompatibility of frames, we provide the representation model of the motion of the module (Sec. 4). Then, we propose an optimization-based design approach using this representation model to obtain the shapes of quadrilateral modules with desired stiffness (Sec. 5). In Sec. 6, we validate the snapping behavior of the designed modules by comparing them with bar-and-hinge analysis and a loading experiment on physical prototypes. We also build a theoretical model for a two-bit memory based on serially connected modules and experimentally validate the multistability property.

Although our final goal is to make surface structures, we first use linkage structures that leave holes. Figure 2 shows some examples of transformable linkages with quadrilateral frames that can create a variety of motions. A Sarrus linkage can create a linear motion of two parallel panels (Fig. 2(a)). After filling the holes with proper origami modules, it can form a flatly foldable closed tube. Another linkage can create a similar linear extension but with shrinking motion in its radius direction (Fig. 2(b)). Fuller’s Jitterbug structure is an expandable polyhedron composed of corner-connected polygons with six quadrilateral holes. This forms a volume-changing closed polyhedron after filling the holes (Fig. 2(c)). Hoberman’s Switchpitch [16,17] is based on the polyhedral with vertices alternately moving outward and inward. This can be generalized to a family of linkages that transforms between different polyhedra with the same symmetry, e.g., a rhombic dodecahedron linkage can deform into an octahedron or a cube (Fig. 2(d)).

These example structures can further be assembled linearly or tessellated in 3D to form architectured material, where each structure can be designed as a tunable geometric and functional building block. In this paper, we particularly focus on the Sarrus linkage type that creates a linearly expandable bistable structure with tunable stiffness, but our approach can be applied to a variety of linkage structures.

A quadrilateral boundary module is a polyhedral surface homeomorphic to a disk with a quadrilateral boundary, which is not necessarily developable. We consider the kinematics based on the rigid origami model, i.e., the model composed of rigid panels and revolute hinges. A generic quadrilateral boundary module forms a one degree-of-freedom (DOF) mechanism when the module is composed of triangular panels; here, the DOF of triangulated origami is computed as the number of boundary edges minus 3 [18].

The mechanism of quadrilateral boundary module origami has been studied by Demaine and Ku [12]. They focus on the developable module and represent the motion of the module as a curve on two-dimensional parameter space using diagonal lengths. They showed the universality of the design of a module that satisfies multiple boundary configurations. We use a wider family of quadrilateral boundary modules without developability constraints for our design. Also, we provide an improved parameter space that provides a reasonable incompatibility metric.

Quadrilateral boundary modules have been implicitly or explicitly used for the design of approximated bellows structures. As bellows are not realizable with rigid origami models, their transformation requires “cheating” or a small elastic deformation of each panel. Goldberg’s tristable polyhedron is one such mechanism, which can be interpreted as a mechanism formed by two quadrilateral boundary modules. Ma and You designed a collapsible tube by attaching triangulated quadrilateral boundary modules to the corners of a Sarrus linkage [19]; the same composition of linkage and modules for our design. Hoberman designed foldable polyhedral near-bellows with bistability by combining modules with quadrilateral boundaries [10,11]. Tachi proposed cap structures that conform to the motion of the ends of rigid origami tube [20]. Our study is inspired by these studies of approximately achieving prescribed impossible behaviors. We further aim to additionally tune the snapping stiffness by properly defining the incompatibility between the frames and the modules.

We propose a family of quadrilateral boundary modules called mirror-symmetric biconical modules, which is a family that contains the four-panel modules in Ref. [11] and some of the developable modules in Ref. [12]. Refer to Fig. 3. We give a simple open polyline V₀V₁…V_n, called profile polyline, that lies on a plane, mirror plane, and two points O₁, O₂, pivots, symmetric about the plane. A mirror-symmetric biconical module is a quadrilateral boundary module that is a union of two generalized polyhedral cones formed by extruding profile polyline to each of the pivots. A mirror-symmetric biconical module consists of n pairs of triangular panels and has a sequence of n − 1 quadrivalent interior vertices along the profile polyline, which produce the one-DOF mechanism.

We describe the kinematics by finding the profile polyline staying on the mirror plane while the pivots move away from or closer to each other symmetrically about the mirror planes. The initial completely flatly folded state is given when O₁ and O₂ are on the mirror plane. As O₁ moves away by distance h from the plane, each triangle rotates about the edge of the profile polyline. The profile polyline continuously transforms to keep the lifted triangles connected.

The specific motion can be given as follows. Refer to Fig. 4. Let T_min be the triangle with the shortest height h_min, i.e., the distance from the pivot to the edge on the profile polyline. Then, h = h_minsinα > 0 increases as α > 0 increases, where α denotes the rotated angle measured between T_min and the mirror plane. Each triangle T_i with distance h_i is then rotated by $arcsin(hmin/hi)sinα$ to form a one-DOF motion. The distance d takes its maximum when α = 90 deg; we call the corresponding folded state the fully opened state.

The entire deformation behavior of a module from the flat state to the fully opened state is unique. However, after the module hits the fully opened state, the mechanism bifurcates into multiple modes by the flipping of the panels (Fig. 5). Because d = h_minsinα = h_minsin(180 deg − α), T_min at α = 90 deg can either fold in (90 deg, 180 deg) (flips) or fold back in (0, 90 deg) with the same pivot positions d. Each triangle T_i with larger height h_i > h_min cannot flip because it cannot take the rotation angle of 90 deg. If there are multiple triangles T_min with minimum height h_min, then each of T_min can either flip or not, creating $2#ofTmin$ different modes. Although the mode transition by flipping is an interesting phenomenon by itself, it is necessary to prevent it in our design as it can cause unpredictable transitions of states.

The quadrilateral boundary modeled as a four-bar linkage has two DOF, while the rigid origami module inside the quadrilateral boundary is a one-DOF structure. Hence, the deformation of the quadrilateral boundary module can be represented as a one-dimensional path (configuration space) in a two-dimensional coordinate space (parameter space). Natural parameterization of the linkage is to use diagonal lengths as proposed in Ref. [12]; however, the use of diagonal lengths results in a singularity at the flat state, meaning that the infinitesimal (first-order) motion of the frame can be achieved without a first-order length change of diagonals. This is problematic when we use the parameter space to measure the incompatibility between different motions. We provide a representation, such that the distances between configurations can be interpreted as incompatibility metrics with physical meaning.

Figure 6 shows the overview of our parameterization. We first align the frame ABCD such that C is the origin, A is on the z axis, BD is parallel to the x axis, and then scale such that each edge length is 1. Because the edge lengths are equal, every non-degenerate configuration is mirror symmetric with respect to the yz plane and a plane parallel to the xy plane passing through B and D. Then, the entire motion is represented by the motion of point B on a unit sphere, which can be expressed using the coordinates (ϕ, θ) on a sphere, where ϕ ∈ [ − 90 deg, 90 deg] is the azimuthal angle and θ ∈ [0, 90 deg] is the polar angle.

To provide the metrics related to the strain, we further tweak the parameters. We define point P as the nonuniform scaling of B by $1/2$ in the y-coordinate direction, as the reason for scaling is explained in the next section. We call P the representation point of the boundary configuration and the set of P corresponding to the folding motion the representation curve of the quadrilateral boundary module. The representation curve is the configuration space that lies on the ellipsoidal parameter space.

Note that the distortion of the module in actual material is realized by the combination of backlash at the crease and the deformation of panels. In our study, we use an approximated evaluation of elastic strain by assuming that the deformation is concentrated at the vertices. More specifically, we consider rigid origami models of frames and modules connected with zero-length springs attached at the vertices. Under the above-mentioned assumption, we show that the Euclid distances between configurations on the ellipsoidal parameter space can measure the elastic strain.

Consider the quadrilateral boundary ABCD of the module obtained as rigid origami without any strain and the target quadrilateral boundary A′B′C′D′ after distortion. Now, we assume that there exists a zero-length spring with spring constant k between the corresponding two vertices AA′, BB′, CC′, DD′ (Fig. 7). We estimate the strain energy of the distorted module by the elastic energy E, equivalent strain energy, stored in the zero-length spring.

Now, we propose the design scheme of the bistable structure using the representation model. As an example case, we designed a linearly extending bistable structure based on the Sarrus linkage by filling its quadrilateral frames with mirror-symmetric biconical modules. Note that the following procedure can be applied for designs for other types of linkages that are shown in Sec. 2. We let the length of one side of the quadrilateral frame be 10 and the height of the Sarrus linkage be l. Since the fully extended state is a singular state of the original linkage, we avoid the singular state by setting the range in 30 deg ≤ θ ≤ 90 deg, so l moves in 0 ≤ l ≤ 20cos 30 deg ≈ 17.3 (Fig. 8). Here, θ is the polar angle of B of the frame, which is equivalent to the rotation angle of a panel of the Sarras linkage. Therefore, we consider solving the motion between the flat state P₀(ϕ, θ) = (30 deg, 90 deg) and the deployed state P₂(ϕ, θ) = (30 deg, 30 deg).

To solve this design problem, we propose the following optimization scheme. First, we use a mirror-symmetric biconical module with six-vertex polyline and use its parameters as design variables. Next, within the search range of this design variable, we searched for a module that has a representation curve that is close to the target curve as much as possible by optimization. This provided a structure that approximately forms the desired mechanism with a very small strain, called a quasi-mechanism structure. Then, using this quasi-mechanism as the reference, we create bistable structures with tuned stiffness of snapping by designing a module with a representation curve that intersects the target curve at specific points but passes through a point of maximum specified distance from the target curve. The variables and objective functions for both steps are explained below.

We implemented the design scheme using Rhino/Grasshopper to explicitly code the kinematics of biconical modules, create the representation curves, and solve the optimization problem using the evolutionary algorithm of Galapagos [21]. The method proposed is general and can also be applied to the design of various multistability structures using other types of linkages introduced in Sec. 2.

We use the mirror-symmetric biconical module of six vertices, parameterized at the flat state on the xy plane (Fig. 9). We place the pivots on the origin and consider a further symmetry about y axis and use the polar coordinates of two interior vertices (r₁, γ₁) and (r₂, γ₂) as design variables. Here, the boundary vertex is fixed such that r₀ = 10 and γ₀ = 60 deg to have a valid configuration at the flat state. Also, to prevent self-intersections of the polylines, we also set the constraint 0 deg ≤ γ₂ ≤ γ₁ ≤ 60 deg.

To design a bistable structure, we let C_m(γ₁, γ₂, r₁, r₂) pass through two stable points P₀, P₂ on C_t and a point P₁ not on C_t between the two points (Fig. 10, right). The distance from P₁ to the C_t is the estimated equivalent strain causing the potential energy barrier between the stable points and is denoted by $ε1$⁠.

Here, we found that P₁ cannot be arbitrarily chosen. The heuristics for choosing P₁ is to first refer to the shape of the representation curve of the quasi-mechanism module and find the farthest point with maximum equivalent strain. This is around (ϕ, θ) ≈ (30 deg − δ, 50 deg), where δ ≈ 3 deg corresponding to $ϵqMax=1.98%$ in our case. Then, we set P₁ to be reasonably close to this point but further away from the C_t to achieve the prescribed snapping stiffness.

The result of the optimization is shown in Figs. 11 and 12. For a quasi-mechanism structure, we obtained a module where V₀, V₁, V₂ are almost in a straight line (Fig. 11(1a)), thus is approximately realizable only with a four-vertex polyline. The representation curve almost follows that of the frame, where the maximum equivalent strain $εqMax$ of $1.98%$ is obtained at P_max(28.1 deg, 50.8 deg). The objective function G converged to 0.000152, which indicates the average (specifically, root mean square) equivalent strain of $1.23%$⁠. As already explained, P_max is then used as the reference for determining P_1a and P_1b.

For bistable structures, we computed module A with a smaller strain and module B with a larger strain. (A) For module A, we set P₁ to be P_1a(25 deg, 45 deg) which corresponds to the target maximal equivalent strain of $ε1a=4.74%$⁠. This resulted in the maximum equivalent strain for module A is $εqMax=6.22%$ at P_max = (23.5 deg, 45.5 deg) (Fig. 11(2a)). (B) For module B, P_1b(20 deg, 45 deg), we set P₁ to be P_1b(20 deg, 45 deg), which correspond to target maximal strain of $ε1b=9.38%$⁠, respectively, as shown in Fig. 11(3a). This resulted in the maximum equivalent strain for module B is $εqMax=9.35%$ at P_max = (19.4 deg, 41.8 deg). There are problems converged to have γ₂ = 0 with degenerate vertices, so it is realizable only with a five-vertex polyline. The representation curve of the obtained modules almost passes through given P₁ and P₂ (Figs. 11(2b), 11(3b), 11(2c), and 11(3c)); the distances were d_1a = 0.015, d_2a = 0.0094 for module A and d_1b = 0.0023, d_2b = 0.00032 for module B.

Finally, the objective function of modules A and B converged to H = 5.62 × 10⁻⁵ and H = 3.20 × 10⁻⁵ with d_flip = 0.261 and d_flip = 0.258, respectively. The representation curve for module A satisfies $d1/ε1≈0.31<1$⁠, $d2/ε1≈0.20<1$⁠, and $ε1/dflip≈0.18<1$ as desired. Similarly, the representation curve for module B satisfies $d1/ε1≈0.024<1$⁠, $d2/ε1≈0.0034<1$⁠, and $ε1/dflip≈0.36<1$ as desired.

By using the obtained modules, we can design deployable tubular structures by applying the modules to the series of Sarrus linkages (Fig. 12, bottom). The tubular structure with quasi-mechanism modules in Fig. 12(1) (bottom) can be used as mechanical elements such as bellows or linear pneumatic actuators. On the other hand, when n bistable structures are connected serially as in Figs. 12(2) and (3) (bottom), we obtain a multistable structure with 2ⁿ stable positions.

We made a physical model of the quasi-mechanism module and bistable modules A and B with a unit length(r₀) of 35 mm as in Figs. 13(a)–13(d). We used heat-adhesive laminate film of a thickness of 0.25 mm as rigid panels and non-woven fabric (polyester) to create a foldable hinge of close to zero stiffness. To avoid the overlapping of thick panels on the valley fold line, we gave a 0.7 mm gap between the panels.

Figure 14 shows the detailed fabrication process: (1) place a laminate film as the adhesion side on the bottom of the cutting mat and cut with a cutting plotter, (2) attach a temporary sticker sheet onto the laminate film, (3) remove the part that does not become the panel of the module, (4) place a non-woven fabric on the adhesion side of the laminate film and thermally bond with an iron, and (5) remove the sticker sheet. This will complete one-half of the structure. Now we attach two pieces with double-sided tapes at the non-woven fabric surrounding the boundary. Finally, we cut off the unwoven fabric at the corner to reduce the excessive stiffness.

From hands-on interaction with the physical models, we experienced snapping between the flat position and the fully deployed state for all three modules. We also found that when we apply unbalanced force to bistable modules A and B, the structure may be locally stable with one module deployed and two modules folded. This unintended multistability effect is due to the limited stiffness of the quadrangular panels in the Sarrus linkage (Fig. 13(d)).

We first verify the bistability of the mechanism from one unit of the Sarrus linkage with its holes filled with each computed module (Fig. 12, top) by comparing it with the bar-and-hinge mechanism. Consider a linear motion changing the height l ∈ [0, ∼ 62] while applying the load perpendicular to the upper and lower surfaces.

Using the equivalent strain energy computed from the geometric incompatibility and its derivative with respect to the displacement, we can obtain the expected plot of energy–displacement and force–displacement plots. Red curves in Fig. 16 show the estimated curves from the geometric model. Note that the vertical axis for the energy or the force is relative in our model, so they are transformed to fit with the results of bistable module A.

For comparison, we performed structural analysis using a bar-and-hinge model [22]. We used Merlin [23,24] implementing a bar-and-hinge model on matlab to obtain force–displacement plots by applying vertical displacement at the three points on the top triangle while keeping the bottom triangle fixed.

The material parameters of the bar-and-hinge model were as follows: Young’s modulus (E*) is set arbitrarily as we compare only the relative stiffness, while we choose the Poisson ratio and the thickness of the sheet to be 0.3 and 0.25. Also, we let the hinge stiffness of fold lines be effectively zero by setting L* = 100,000 [23]. Blue curves in Fig. 16 show the results from the bar-and-hinge model analysis.

We conducted a cycle of extension and compression using the physical prototypes of bistable modules A and B, which are fabricated in Sec. 6.1. The bottom triangle of modules was attached to the base made of foam board, and the top triangle of the module is attached to the force gauge (ZTA-50N Imada Corp.). The connection part between the force gauge and the upper triangular surface of the module is fabricated using the fused filament fabrication type 3D printer. To eliminate the effect of air pressure, we punched a hole in the rectangle panel of each module (Fig. 21). The speed of the cycle of the force gauge is 50 mm/min.

The outlines of the colored areas in Fig. 15 show the force–displacement curve for a cycle of the extension and compression experiment. Since the physical experiment contains the effect of the plasticity and elasticity of hinges, we cannot directly compare it with the result from the equivalent strain or the bar-and-hinge model. We separate the plasticity and elasticity effects due to the hinges between panels in the following way. As a preprocess, we applied noise removal using Gaussian window filtering of standard deviation σ equal to $0.6%$ of the total displacement of the cycle.

The plasticity or friction of the hinges is observed as the hysteresis curve (colored areas in Fig. 15) of the single unit curve, which would not exist in a conservative system. To separate this effect, we average the load for the forward and backward paths, as shown in the bold-colored curves in Fig. 15. However, since the plot tends to be noisy due to a series of small buckling near the flat state, we replaced the curve near the zero point (0–2.6 mm for the quasi-mechanism module, 0–4.66 mm for the bistable module A, and 0–7.61 mm for the bistable module B) with a linear extrapolation from the half displacement position of the maximum extreme to the last displacement of the buckling range (2.6–15.83 mm for the quasi-mechanism module, 4.66–18.18 mm for the bistable module A, and 7.61–20.67 mm for the bistable module B). The plastic effect leads to the drift of zero points, as shown by the dashed black lines in Fig. 15. We calibrated the zero point by letting the average load (black curves in Fig. 15) become zero.

Refer to Fig. 16 for a comparison between models and the experiment. Note that since we did not specify the material property in the model, the values in the potential energy and force are relative. Therefore, we scaled such that the highest potential energy in module A is equal in the equivalent strain model, bar-and-hinge model, and experiment.

The dots in Fig. 16(a) show the displacement of the peak of each energy result. As designed, we obtained the energy barriers of the quasi-mechanism module being the lowest and the energy barrier of module B being the largest in both the bar-and-hinge model and the experiment. For all three models, we obtained a force–displacement plot in Fig. 16(b) causing the snapping behavior between flat and deployed states. This validates the design method for differentiating the barriers on a qualitative level. However, we also observed that the difference in magnitude of the energy barriers between the three modules is smaller in the bar-and-hinge model and even smaller in the experiment compared to the equivalent strain model. We consider that this is due to the nonlinear effect of the difference between the concentrated strain at the vertices (as modeled in the design model) and the actual strain dispersed to surfaces.

The displacement that takes the peak value and stable points of the potential energy is also similar; however, the displacement of peak values is slightly shifted to the smaller in the bar-and-hinge model and experiment (Fig. 16 (top plots)). The quasi-mechanism module shifts by $1.13%$ (43.75 → 43.05) in the bar-and-hinge result, and $6.95%$ (43.75 → 39.44) in the experiment, of the total stroke, module A shifts by $5.03%$ (48.51 → 45.39) in the bar-and-hinge result and $8.06%$ (48.51 → 43.51) in the experiment, and module B shifts by $8.02%$ (51.54 → 46.57) in the bar-and-hinge result and $7.94%$ (51.54 → 46.62) in the experiment.

The serial composition of bistable modules with different stiffness can take 2ⁿ positions that can be controlled from a sequence of motions on one side; such a behavior can be interpreted as a n-bit mechanical memory [8]. Serially connecting modules A and B (Fig. 17) are expected to be a two-bit memory that exhibits a cycle of transitions between 2² states.

First, we make a theoretical prediction of the two-bit memory hysteresis from the force f and the potential energy E_A and E_B of the bistable modules A and B. We consider a three-dimensional space spanned by displacements l_A, l_B, and the force f as shown in Fig. 18. l_A and l_B indicate the heights of bistable modules A and B, respectively. f indicates the equilibrium force of the corresponding displacement, which is identical to the force of the serial system. The force–displacement plots of modules A and B based on the bar-and-hinge model in Fig. 16(b) are plotted in the l_A–f and l_B–f planes as the red and green curves, respectively, in Fig. 18(a). Here, in order to model the collision of panels at the flat state, we added a vertical curve to the force–displacement plot at l = 0, so that a single module cannot take a negative height.

Now, the force equilibrium for module A is described as the surface obtained by the extrusion of the l_A–f plot to the l_B axis. Similarly, the force equilibrium for module B is described as the extrusion surface of the l_B–f plot along the l_A-axis (Fig. 18(a)). Since the equilibrium path of the serial system satisfies that the load on modules A and B are the same, it can be obtained as the intersection of two extrusion surfaces, as expressed by the black bold curve in Fig. 18. We can plot the relationship between the merged length (l_A + l_B) as the projection to a plane containing the diagonal direction, called the P-plane. By scaling this projected curve by $2$ in the horizontal direction, we can obtain the force equilibrium path of the total displacement for the serial system (Fig. 18(b)).

By observing the total potential energy E_A(l_A, l_B) + E_B(l_A, l_B), we can check the stability of the equilibrium path (Fig. 18(c)). The potential energy takes the local minimum at four stable points [00], [11], [10], and [01]. The l_Al_B-plane is partitioned into four regions (orange dashed lines in Fig. 18(d)) corresponding to the four stable points such that the energy gradient descent path from a state in a region ends up in its corresponding stable state. Figure 18(d) shows the projected intersection curve to the l_Al_B-plane with the contour map of the total energy plot. We can check that the force equilibrium path passes through each stable point and region.

The hysteresis path that appears in the forced displacement loading condition can be drawn by choosing the local minimum point for a varying total displacement l_A + l_B. The path has a snap-through behavior that jumps from an unstable equilibrium point to a stable equilibrium point with the same total displacement. Therefore, the snap-through is visualized by a straight line drawn perpendicular to the P-plane on the energy plot (Fig. 19).

Based on the above discussion, in Fig. 19, we were able to identify the behavior of the forced displacement cycle, which is characterized by the three times snap-through between equilibrium points and hysteresis passing through stable regions [00], [11], [10], [01]. The hysteresis curves computed from bar-and-hinge and equivalent strain models are shown in Figs. 20(f) and 20(g), respectively. Note that the forced displacement path actually hits the stable points [00], [11], [10], but not [01]; as the result, the curve passes through f = 0 at [00], [11], [10], but not at [01]. The estimated behavior of the serial system from equivalent strain energy shows a larger distance between the stable points [01] and [10] and a larger negative force than the bar-and-hinge model.

We connected modules A and B. The bottom triangle of module A was attached to the base made of foam board, and the top triangle of module B is attached to the force gauge (Fig. 20(a)). We started at the flat state (Fig. 20(a)) and stretched until both modules were deployed (Fig. 20(c)). Then the load on the displacement of one cycle returned to the flat state was measured. In this experiment, the speed of the force gauge is 50 mm/min.

We obtained four states and the hysteresis transitions between the states as designed as shown in Figs. 20(a)–20(d). In the deployment process from the flat state (a), module A deploys first (state (b)), and module B follows (state(c)). When this is in compression, module A first compresses (state (d)), then module B compresses (state (a)) and comes back to the flat state. During the entire course of the experiment, the local stability problem due to symmetry breaking described in Sec. 6.1 did not occur. This can be attributed to the fact that the forces applied through the designed connections are evenly distributed.

Figure 20(e) shows the force–displacement plot of an extension and compression cycle path. The result shows the two-bit memory hysteresis and snapping behavior very similar to the expected hysteresis path in Fig. 19. We also record that the three stable states of Figs. 20(a), 20(c), and 20(d) are at the zero force position and the stable state of Fig. 20(b) was not zero force as we expected. However, the mismatch between deployment and compression paths was observed even in the same memory states. This unintended hysteresis can be attributed to the absorption of energy due to the bending of the rectangular panel in Fig. 21(a) or the plasticity and the play of the hinge in Fig. 21(b).

Figures 20(f) and 20(g) compare the hysteresis curves, computed through the bar-and-hinge model and the equivalent strain energy, by fitting the second maximum extreme value to the experimental result. Although the presence of plastic and elastic resistances in the hinges and other elements makes a direct comparison difficult, we observed that the experimental model behaved as the two-bit memory system predicted by our models.

We have proposed a novel design approach for multistable structures by filling the linkage structure with quadrilateral boundary modules so that the equivalent maximum strain takes the designed values. The estimated behavior of this design is compared with the bar-and-hinge analysis and the experiment and is shown to be in qualitative agreement, while we observed that the equivalent strain model tends to overestimate the differences in energy barriers of the models. We also performed experiments on the fabricated two-bit mechanical memory to validate the designed behavior and showed that it performed similarly to the expected hysteresis transitions and path. Since our design example was limited to the linear motion produced by the Sarrus linkage, we plan to apply our method to different linkage and tessellation designs to create different structures and materials in the future.

We thank Chuck Hoberman, Erik Demaine, and Jason Ku for helpful discussions and inspiration. This work is supported by JST PRESTO grant number JPMJPR1927.

There are no conflicts of interest.

The datasets generated and supporting the findings of this article are obtainable from the corresponding author upon reasonable request.