A Computational Design Synthesis Method for the Generation of Rigid Origami Crease Patterns

Origami has received considerable attention in recent decades when mathematicians, engineers, and artists recognized the benefits of origami. These benefits are numerous: origami is scale-independent [1], allows for a compact or flat stowage in the unfolded state as well as complex three-dimensional motion during deployment [2]. The deployment can be achieved by a low number of degrees-of-freedom (DOF), which minimizes the amount of resources required to actuate the folding motion and enables a reliable control [3]. The facilitated actuation and the complex motion then enable programmable structures that can change mechanical properties, shape, and function on demand [4,5]. In addition, origami can be produced in the flat state using additive manufacturing [6], which enhances realization possibilities and reduces assembly time.

Due to its benefits, the principle of origami offers potential for various scientific fields that pose different engineering design tasks. The transformation of origami into engineering applications includes an entire spectrum of more abstract to more direct implementations. The abstract end of the spectrum corresponds to so-called origami-inspired products [7] that fold and exhibit origami-like geometry but otherwise show little resemblance to origami. More closely related to traditional origami are origami-adapted mechanisms [8] that are based on origami crease patterns but use nonpaperlike materials and accommodate for finite thickness. Finally, origami-applied systems [7] use paperlike material and exhibit little to no alteration of the underlying crease pattern. Together, these categories constitute the entire spectrum of origami-based design [9,10], which finds applications in mathematics [11], material science [12–15], DNA [16] and biomedical research [17–21], mechanical engineering [22], robotics [23–26], consumer goods [27], architecture [28], and space [29–31].

The widespread applicability of origami is accompanied by the adaptation and introduction of various computational tools [10,32–38] that facilitate the application of origami in engineering design tasks. However, despite the trend toward computer-aided processes, today most crease patterns used in technical applications are selected from a handful of well-known origami principles [2]. Computational algorithms capable of generating novel crease patterns either target artistic origami [39,40] or shape matching [41] without guaranteeing rigid foldability, are limited to quadrilateral creased paper [42] using only one pattern [33] and generation of origami tessellations with discrete sector angles [43]. Furthermore, most of the methods do not embed direct design knowledge on how to generate rigidly foldable crease patterns. Two examples for the latter case are works by Fuchi et al. [44] and Gillman et al. [45] that utilize ground structures that only allow for sector angle configurations that include angles of 45 deg or multiples thereof, which offers limited design freedom. This only leads to regular patterns, and thus resembles an undirected search for feasible solutions rather than the purposeful design of novel crease patterns that satisfy given engineering design tasks. In both cases, a ground structure is used as a genotype, in which each element can be switched on or off while optimizing using genetic algorithms (GA). Further approaches that rely on GA and meta-heuristics to generate solutions are works [46] and [47], respectively. Although those approaches follow encodings similar to Refs. [44,45], they are not limited to a fixed ground structure, but the question of how to perform a search that guarantees rigid foldability of generated candidate solutions is still unresolved. In an effort to improve the performance of GA and topology optimization in Ref. [45], the subsequent work [48] uses Bayesian optimization to design origami for a specific task. Rather than focusing on generating patterns, the focus is set on an efficient response evaluation since nonlinear material properties are considered. Again, a predefined ground structure is used for which the optimization process tunes the locations of vertices to satisfy the objective function.

The lack of computational methods for the generative design of crease patterns for engineering applications arises from a multitude of geometric complexities intrinsic to origami. First, many engineering applications require an origami to fold rigidly for the incorporation of rigid materials or electronics. Second, each internal vertex in a crease pattern exhibits two rigid body modes (RBMs) that allow for both “upward” and “downward” motion [49,50], which results in an exponential growth of the number of possible modes as the number of internal vertices increases. Third, the relation between kinematic determinacy, DOF, and pattern symmetry is still largely unexplored in origami. Fourth, real-world applications do not allow for self-intersection for which there still exist no intrinsic conditions. Finally, converting infinitely thin patterns into realizations with finitely thick materials imposes problems on crease pattern and hinge design [51].

In a recent work [52], whose findings will be briefly explained in Sec. 2, the authors addressed some of the above complexities by presenting the Principle of Three Units (PTU) that leads to an analytical, intrinsic condition for the rigid foldability of single degree-n vertices, where n is the number of crease lines incident to a vertex. In addition to this condition, the PTU yields an analytical kinematic model for single vertices to which the RBMs are inputs rather than the results of simulation [53]. The PTU [52] further offers guidelines for the generation of crease patterns whose kinematic motion can be explicitly determined. The condition for the rigid foldability of a single vertex can then be extended to each vertex in the generated crease pattern to achieve global rigid foldability if the underlying crease pattern is acyclic.

Based on these findings, in this paper we introduce a Computational Design Synthesis (CDS) method [54] for the generative design of novel crease patterns to develop origami concepts for engineering applications. We focus on the generation of simple crease patterns that fit the scope of engineering applications, and use the term “concept” to clarify that the method addresses the abovementioned geometric complexities except the problem of finite thickness. This leads to a method that generates rigidly foldable, infinitely thin, and intersection-free crease patterns that satisfy a given design task.

Section 2 outlines the PTU and the corresponding findings, after which the proposed CDS method is presented in Sec. 3. The method for the generative design of novel origami concepts is then applied to two engineering tasks that involve the design of rigid origami grippers and robotic arms are given in Sec. 4. Subsequently, the results are presented, the method is discussed in Sec. 5, and the contributions and findings of the paper are summarized in Sec. 6.

This section describes the PTU and the condition for the rigid foldability of degree-n vertices, after which the corresponding kinematic model is summarized. The section concludes with the implications for crease pattern generation.

The PTU [52] predicates on the fact that every single degree-n vertex requires n − 3 inputs and 3 outputs to fold in a kinematically determinable manner [55]. The inputs are the dihedral angles that drive the motion of the vertex, here called the driving angles, which need to be prescribed in order to fold the vertex. The outputs are the three remaining dihedral angles, here called the unknown dihedral angles, whose behavior is determined by the driving angles and the size of the sector angles around the vertex. As an example, Fig. 1(a) depicts a degree-8 vertex with five driving angles ρ₁₋₄ and ρ₆ assigned to an arbitrary set of crease lines, and the corresponding unknown dihedral angles ρ₅, ρ₇, and ρ₈.

Independent of the chosen set of driving angles and the degree of the vertex (for n ≥ 4), virtually cutting the vertex at the locations of the unknown dihedral angles reveals three parts [55] called unitsu₁, u₂, and u₃ (Fig. 1(b)) [52]. What remains within these units are the sector and driving angles determined by the user, allowing for the entire kinematic behavior of each individual unit to be expressed analytically. Let a unit u_i be denoted with an ordered set comprising the sector and driving angles within the unit, $ui=(αui,ρui)$⁠, where the set of driving angles can be empty if a unit consists of a single sector. For the example in Fig. 1, this notation results in u₁ = ((α₁₋₅), (ρ₁₋₄)), u₂ = ((α₆, α₇), (ρ₆)), and u₃ = (α₈).

Every crease pattern is perceived as a network of vertices (nodes) that transform n − 3 inputs into 3 outputs [55], turning an origami into a directed graph where the inputs and the outputs are represented as incoming and outgoing edges, respectively. For an internal vertex v_i of degree n_i, its indegree (number of incoming edges) $deg−(vi)=ni−3$ linearly scales with n_i and is thus unrestricted as long as n_i ≥ 4. However, its outdegree $deg+(vi)$ (number of outgoing edges) must always satisfy the condition $deg+(vi)=3$⁠. Thus, any number of incoming edges can be added to a vertex as long as the vertex has three outgoing edges. This is the first condition that needs to be satisfied while generating crease pattern graphs. The second condition stems from the kinematic model of the PTU that requires a graph to be acyclic in order to model its motion.

The proposed CDS method (Fig. 2) represents crease patterns as graphs and uses a graph grammar system to generate all possible crease pattern topologies (“Topology” in Fig. 2). The inputs to the method are an initial graph G₀, which serves as a starting point for the graph grammar generation of origami concepts, and a maximum number of internal vertices N_max contained within these concepts. Both G₀ and N_max are user defined. The generation uses two different rule types, rule r₁: extend vertex, and rule r₂: combine vertices, within a matching process to identify the vertices to which the two rules can be applied. Simultaneously, the graph grammar automatically introduces the variables and constraints of the generated graph topologies for the optimization in subsequent steps. The generated topologies are then checked for redundancy and eliminated based on isomorphism. Further elimination criteria can be supplied manually.

Subsequently, the method steps into a loop (“Geometry” in Fig. 2) in which a graph is first embedded within three dimensions by associating it with the vertex coordinates and the relevant variables to represent its geometry. Then, each graph is associated with all different RBM assignments to its internal vertices, corresponding to the enumeration of all design alternatives. Using the formulation of a design task provided by the user that involves an objective function and other optimization-related inputs, the evaluation step first optimizes each alternative and, if an optimized alternative meets the design criteria supplied with the design task, subjects it to an intersection check. A design alternative that fails this intersection check is discarded, otherwise it is added to the list of feasible origami concepts. In the following subsections (Secs. 3.1–3.4.2), we detail the method according to the steps in Fig. 2.

A crease pattern is represented as a directed, labeled graph G = (V, E, L_V, L_E), where V is a nonempty set of vertices v_i and E is a nonempty set of directed edges e_i,j from vertices v_i to v_j such that i ≠ j. Vertices v_i are called the predecessors of v_j and v_j are called the successors of v_i. L_V is a nonempty set of vertex labels and L_E is a nonempty set of edge labels. Let Σ_V be a map Σ_V:V → L_V that labels each vertex v_i to its ordered predecessors P_i and a type T_i, so that each vertex label L_i in L_V is an ordered pair defined as L_i = (P_i, T_i). If a vertex v_i has no predecessors, $Pi=(∅)$⁠, then v_i is a source vertex. Otherwise, P_i is populated with all predecessors of v_i ordered in the counter-clockwise direction with respect to the crease pattern graph in the neighborhood of v_i. The type T_i of vertex v_i adopts the symbol χ when the vertex can be extended with new outgoing edges, or it adopts the symbol $Ti=∅$ when the vertex is either a source or has already been extended. Let Σ_E be a map Σ_E: E → L_E that assigns an edge label L_i,j in L_E to each directed edge e_i,j, where $Li,j=(FLi,j,FRi,j)$ is an ordered pair of the facets f located on the left and on the right side of e_i,j, respectively. While the sets of vertices and edges define the topology of a graph, the labels are introduced to control the generation of new graphs.

To clarify the notation, Fig. 3 shows the simplest possible initial graph G₀. This graph consists of just two vertices v₁ and v₂ connected by a directed edge e_1,2. The vertex label $L1=((∅),∅)$ denotes that the vertex v₁ is a source and cannot be extended, respectively. The vertex label L₂ = ((v₁), χ) denotes that v₁ is the single predecessor of v₂ and that v₂ can be extended, respectively. Finally, the two components $FL1,2=f1$ and $FR1,2=f2$ of the edge label L_1,2 signify that the edge e_1,2 is adjacent to the (yet invisible) facet f₁ on its left and f₂ on its right side.

In addition to the initialization associated with the representation of the topology, the user has to relate the initial graph to its engineering purpose and thus to its geometry. This relation is achieved by allocating in-plane coordinates $xi$ to all vertices v_i and defining a driving angle ρ_i,j for every incoming edge incident to a vertex v_j that can be extended. As an example, for the graph in Fig. 3 this could be defined as $x1=(0,0)$⁠, $x2=(1,0)$⁠, and ρ_1,2 = t where t represents a linear driving angle that goes from zero to some value t_max ∈ [−π, π]. In short, Fig. 3 thus represents an initial graph with a single extendable vertex (v₂), with a geometry (vertex coordinates), and with a single DOF (driving angle ρ_1,2 = t).

In general, the definition of the initial graph is an input to the method and can involve different topologies, including multiple DOF. However, the two guidelines given above also restrict the definition of possible initial graphs. First, all initial graphs must be acyclic. Second, the outdegree of all initial vertices must be smaller than 3; the method models the kinematics of a vertex while extending a vertex, and it does not extend vertices that already satisfy $deg+(vi)=3$⁠.

The origami graph grammar GG is defined by the triple $GG=(G0,R,∅)$ where G₀ is the initial graph, $R=(r1,r2)$ is the set of rules containing r₁ and r₂, and $∅$ is the terminal symbol with respect to the type T_i of a vertex v_i that prevents rule r₁ from being applied to a vertex with $Ti=∅$⁠.

The rule set $R=(r1,r2)$ is designed to conform to the two guidelines for the generation of rigidly foldable crease pattern graphs. The first rule r₁ extends a vertex by new outgoing edges and incident vertices and ensures that each extended vertex has three outgoing edges after the application. The second rule r₂ combines two vertices and their incoming edges to enable the generation of higher order vertices (degree-5, 6, etc.) while guaranteeing that the generated graphs stay acyclic. In addition to the graph transformations, the rules r₁ and r₂ embed a graph G in the plane, model extended vertices kinematically, and build or adjust the sets of the optimization variables and constraints. To automate this process, the method initializes empty sets for the in-plane coordinates x, the units u, the dihedral angles ρ, the global rotation matrices R, the optimization variables Φ, and the optimization constraints ψ. Moreover, the method initializes a set of three-dimensional vertex coordinates $X={Xi}$ where $Xi$ are the in-plane coordinates of all vertices v_i contained within the initial graph G₀, which are appended with a zero that represents the z-coordinate. Then, these sets are automatically filled and adjusted when the rules r₁ and r₂ are applied.

Rule r₁ is a production of the type r₁: LHS₁ → RHS₁, where the LHS₁ is a single vertex v_a to which the rule is applied and where the RHS₁ is the same vertex v_a with a number of new outgoing edges and incident vertices. Depending on the initial graph G₀, the outdegree $deg+(va)$ on the LHS₁ can be 0, 1, or 2, and r₁ should accordingly generate $3−deg+(va)$ new edges and vertices. Thus, rule r₁ is parametric with respect to the number of successors of v_a on the LHS₁, which results in six possible rule application scenarios r_1,1 − r_1,6. For clarity, only the scenario r_1,1 (Fig. 4) is described in detail, whereas all other rule application scenarios r_1,2 − r_1,6 are listed in a more concise form in the Supplementary Material (see Supplemental Figs. 1 and 2 available in the Supplemental Materials on the ASME Digital Collection, and Supplemental Tables 1 and 2 available on the ASME Digital Collection). Figure 4 depicts the case r_1,1 in which v_a has no successors on the LHS_1,1, where LHS_1,1 and RHS_1,1 are highlighted in red and the graph surrounding v_a in gray as a reference. This surrounding graph consists of any number of predecessors $vpm$ for 1…m ordered in the counter-clockwise order, where the respective edges are divided by the sector angles $αp1$ to $αpm−1$⁠.

A match $M1,1$ of the LHS_1,1 is found if v_a is extendable as defined by its type T_a = χ and if $deg+(va)=0$ (which corresponds to this specific scenario).

On the RHS_1,1, rule r_1,1 generates three new successor vertices $vs1−3$ numbered in counter-clockwise order that are connected to v_a with directed edges $ea,s1−3$⁠. All generated vertices have identical vertex labels $Ls1−3=((va),χ)$⁠, making them applicable for a successive rule application of r₁. While the edge labels adjacent to existing facets adopt the edge labels of the corresponding edges (⁠$FLp1,a$ and $FRpm,a$ in Fig. 4 right), two new facets are generated and denoted with $FLa,s1$ and $FLa,s2$⁠. Hence, the edge labels result in $La,s1=(FLa,s1,FRpm,a)La,s2=$$(FLa,s2,FLa,s1)$⁠, and $La,s3=(FLp1,a,FLa,s2)$⁠. All newly generated indices s₁₋₃ are in accordance with the enumeration of the nodes in G prior to the rule application. The application of r₁ to a graph G, $G⇒r1G′$⁠, results in the graph $G′=(V′=V⊎{vs1−3},E′=E⊎{ea,s1−3}$⁠, $L′V=LV⊎{Ls1−3},L′E=LE⊎{La,s1−3})$ where $⊎$ represents the disjoint union of sets [57]. In addition, the type of the vertex v_a is set to $Ta′=∅$ on the RHS_1,1 to prevent any further extension of v_a.

In addition to the sector angles and edge lengths, all scenarios of r₁ introduce an optimization variable M_a that stands for the RBM of the vertex v_a. The set of optimization variables is then expressed as $Φ′=Φ⊎{αs1−3,la,s1−3,Ma}$⁠.

If v_a exhibits only a single predecessor, the sum over $αpj$ is neglected. Then, the units in Eq. (6) are added to the set of units, $u′=u⊎{u1a−3a}$⁠. By setting these units and the RBM M_a into f_φ, the method determines the dihedral angles $φ1a−3a$ and associates them with the generated edges, $ρa,s1−3=φ1a−3a$⁠, which are then added to the set of dihedral angles $ρ′=ρ⊎{ρa,s1−3}$⁠.

The expression of units and dihedral angles enables the method to kinematically model each extended vertex individually. However, the locations of vertices need to be expressed in three-dimensional space to represent an actual origami, which requires the method to transfer the individual (local) kinematics to a global kinematic behavior. To do so, the method first finds the local rotation matrices $R1a−3a$ of v_a by applying $fR$⁠. Then, these local rotation matrices are turned into global rotation matrices, $Ra,si=Rp1,aRz(π)Ria$⁠, which are added to the set of global rotation matrices $R′=R∪{Ra,s1−3}$⁠. A special case arises when v_a belongs to the initial graph since the global rotation matrix $Rp1,a$ has not been determined by an application of r₁. In this case, the method sets $Rp1,a=Rz(αp)$ where α_p is the angle between the edge $ep1,a$ and the x-axis. As described in Sec. 2, the first crease line of a vertex is locally associated with the x-axis, with which the method expresses the three-dimensional coordinates $Xsi$ as $Xsi=Xa+la,si*Ra,si⋅(1,0,0)T$⁠. These coordinates are subsequently added to the set of three-dimensional vertex coordinates, $X′=X∪{Xs1−3}$⁠.

The second type of constraint introduced by r₁ corresponds to the condition for rigid foldability in Eq. (1) that states U_max(t) ≤ U_med(t) + U_min(t). The difficulty for the application of this condition arises from the fact that the unit angles change their size over time depending on the parameter t, which effectively means that Eq. (1) has to be satisfied for the whole folding motion.

The constraints in Eqs. (7)–(9) are then added to the set of optimization constraints ψ′ after each application of r₁.

The second rule r₂ is defined as a production of the type r₂:LHS₂ → RHS₂ in Fig. 5, where LHS₂ and RHS₂ are both shown in red, the surrounding graph is shown in gray, and adjustments to the constraint system are shown in blue.

Rule r₂ accepts two vertices $vc1$ and $vc2$ that share an edge label on the LHS₂ and combines them into a single vertex $vc1$ on the RHS₂. The incoming edges formerly incident to $vc2$ are reassigned to $vc1$ and a new constraint is introduced between the reassigned edges and the edges incident to $vc1$⁠. Then, r₂ adjusts the set of optimization variables and constraints by replacing all instances of the lengths and sector angles corresponding to the reassigned edges, such as $lp2,c2$ and $αc2$ in Fig. 5 left, respectively, with the parts colored in blue in Fig. 5 right. The vertex v_b signifies any vertex that was previously coupled to the sector angles of the reassigned edge (⁠$αc2$⁠).

A match $M2$ of the LHS₂ in G involves four conditions. Both vertices $vc1$ and $vc2$ have to be extendable, $Tc1=Tc2=χ$⁠. This condition guarantees the generation of acyclic graphs since no vertex can be combined with any of its predecessors.

The second condition requires that both vertices $vc1$ and $vc2$ belong to the same facet. This can be checked by comparing the edge labels of the incoming edges incident to $vc1$ and $vc2$⁠, which is shown in Fig. 5 on the LHS₂ where $FRp1,c1=FLp2,c2$ has to be satisfied. Note that while Fig. 5 illustrates the case in which both vertices $vc1$ and $vc2$ each exhibit one predecessor, multiple predecessors to both vertices are possible. In this case, the method compares the edge labels of the first and last predecessors of both vertices. Identifying these first and last predecessors is straightforward since the predecessors of a vertex are always ordered in the counter-clockwise direction as explained for rule r₁.

To prevent the combination of vertices that share the same predecessor, such as $vs1−3$ in Fig. 4, the third condition states that $vc1$ and $vc2$ cannot exhibit the same predecessors, $Pc1∩Pc2={∅}$⁠.

Since r₂ deletes a vertex, the fourth condition states that $vc1$ and $vc2$ cannot both be contained in the initial graph, which can be assessed by checking if their in-plane coordinates are dependent on any optimization variables. If one vertex is contained within the initial graph, it becomes the “dominant” vertex that is not deleted, and if neither belongs to the initial graph then either vertex can be eliminated. The following graph transformation describes the generic case in which $vc1$ is dominant.

Once r₂ is applied, $G⇒r2G″$⁠, the RHS₂ includes $vc1$ while $vc2$ is subtracted from the set of vertices, $V″=V∖{vc2}$⁠. The edge formerly incident to $vc2$ is subtracted from the set of edges and a new edge is introduced from the predecessors of $vc2$ to $vc1$⁠, $E″=(E∖{ep2,c2})⊎{ep2,c1}$⁠. The predecessors of both vertices are united in the correct order such that the predecessors of the vertex contributing the right side of the edge label (⁠$FRp1,c1$ in Fig. 5) go first. In addition, vertex $vc1$ on the RHS₂ can still be extended, resulting in the vertex label $Lc1″=((vp1,vp2),χ)$⁠. The edge label $Lc2$ is then subtracted from the set of vertex labels, $LV″=LV∖{Lc2}$⁠. The edge label formerly corresponding to the edge incident to $vc2$ is transferred identically to the reassigned edge, $Lp2,c1″=Lp2,c2$⁠. Finally, the edge label of the edge incident to $vc2$ is subtracted from the set of edge labels, $LE″=LE∖{Lp2,c2}$⁠. If $vc2$ is incident to multiple incoming edges on the LHS₂, the above procedure is applied to all edges and edge labels.

The in-plane location and thus the coordinates $xj$ of a vertex v_j are always determined by a single sector angle α_j and a single edge length $lpi,j$⁠. Since r₂ deletes the vertex $vc2$⁠, the corresponding sector angle and edge lengths of $vc2$ are eliminated from the set of optimization variables, $Φ″=Φ∖{αc2,lp2,c2}$⁠.

Equation (10) is then added to the set of optimization constraints ψ.

When $vc2$ has multiple predecessors on the LHS₂, this replacement procedure needs to be analogously performed for all reassigned edges. This also applies when the sector angles and edge lengths are already expressed in terms of the in-plane coordinates, i.e., when r₂ has already been applied to $vc2$⁠.

The rule set $R=(r1,r2)$ is designed to comply with the guidelines for the generation of rigidly foldable and acyclic crease pattern graphs and is compact given the small number of rules required for the generation of such graphs [52]. Together, the rules r₁ and r₂ are able to generate existing crease patterns such as slender origami [59], origami strings [60], or the gripper in Ref. [61]. To automate the generation step (Fig. 2), the method applies the two rules whenever a respective LHS match $M1$ or $M2$ is found and enumerates all possible rule application sequences that arise from the initial graph G₀ and the maximum number of internal vertices N_max.

However, distinct rule application sequences do not guarantee distinct graphs [62], which is why the method checks the generated crease pattern graphs for redundancy and filters them before they are subjected to the guidance step. For generic design tasks, the method filters graphs based on isomorphism while filters specific to a design task have to be introduced manually if required. The filtering of isomorphic graphs is described here, while two additional filters specific to the design tasks of the gripper and the robotic arm are explained in Sec. 4.2.

For the rule system presented, two graphs are isomorphic if they exhibit an edge-preserving vertex bijection. If this condition is satisfied, the set of vertex labels L_V is guaranteed to coincide, whereas the set of edge labels L_E does not play a role in terms of graph morphology (i.e., there is no difference in kinematics when facets are numbered differently). Isomorphism checks are implemented in many standard programs, and the method uses the function IsomorphicGraphQ integrated in Mathematica 11.

Both the automated generation of crease pattern graphs and the isomorphism check are listed in the pseudocode in Supplemental Table 3 in the Supplementary Material available on the ASME Digital Collection. In short, the method automatically enumerates all possible crease pattern graphs that arise from the initial graph and N_max by applying both rules r₁ and r₂ whenever a LHS match is found. As shown in Fig. 2, after the generation step all graphs G are forwarded to the guidance step.

The RBMs M contained in the set of the optimization variables Φ play a special role since every M only adopts the discrete states “up,” M_↑, and “down,” M_↓. Each internal (or extended) vertex in a graph G contributes one variable M, which results in 2^N different design alternatives for a crease pattern with N internal vertices. Due to the lack of more knowledge about the search space associated with the RBMs, the method enumerates all possible RBM assignments and subjects each design alternative to the evaluation.

As illustrated in Fig. 2, the evaluation includes both the optimization of a crease pattern as well as an intersection check. While all design alternatives are optimized, only the optimized alternatives that meet the design criteria are subjected to the intersection check.

Much of what is required for the optimization is intrinsic to the rule system. By applying the rules, the method automatically generates the set of optimization variables Φ and constraints ψ of each graph G, both of which are expressed analytically. What remains to be defined by the user to run the optimization is a design task involving an objective function Ω and the numerical values for the variable boundaries l_min, l_max, α_min, and α_max. The objective function can include and should be dependent on the dihedral angles ρ_i,j or the three-dimensional vertex coordinates $Xi$⁠. In practice, a user defines one or multiple states t = t_opt ∈ [ − π, π] at which the dihedral angles ρ_i,j(t_opt) or the vertex coordinates $Xi(topt)$ should equal certain values or locations, respectively. The definition of the objective function for the exemplary design tasks is given in Sec. 4.3.

Depending on the objective function Ω, the proposed method adjusts the set of optimization variables as a last step before the optimization. For an origami crease pattern, the lengths of the edges that are incident to degree-1 vertices on the border of the paper do not influence the kinematics of the crease pattern. Thus, all variables corresponding to the edge lengths l_i,j that are present neither in the objective function nor in the set of constraints ψ are discarded from the set of optimization variables Φ, resulting in an adjusted set Φ_opt. For the subsequent representation of the origami, the discarded edge lengths in the set of three-dimensional vertex coordinates are substituted with the value of the minimum edge length, X = X(l_i,j ← l_min).

If an optimized design alternative satisfies $Ω*≤ϑ$⁠, where $ϑ$ is the target criterion provided by the user in the design task, the design alternative is subjected to an intersection check that requires the method to first transform the set of optimized variables $Φopt*$ into an origami with facets. To do so, the optimized variables are first substituted into the set of three-dimensional vertex coordinates $X*=X(Φopt←Φopt*)$ that is then only dependent on t. To find the vertices that constitute a facet, the method iterates through all edges and for each edge assigns both incident vertices to the facets contained in the respective edge label. This procedure results in sets of vertices that represent the facets, and the method assigns to each of these sets a polygon. By associating each vertex v_i with its three-dimensional coordinates $Xi*$⁠, a design alternative is represented as a set of polygons that move with respect to t. Subsequently, the design alternative is subjected to the intersection check described in Ref. [64]. Since the intersection check is purely numerical, the method performs the check 16 times throughout the folding motion from t = 0 to t = t_max in even intervals. The specific number of checks is determined empirically to achieve an appropriate trade-off between a low time-consumption and a reliable intersection check.

In this section, the proposed method is applied to two design tasks that include a gripper and a robotic arm. The principle of rigid origami has been used for gripping tasks before, such as in the Oriceps [23] or in Ref. [61], or trajectory following motions for surgical applications [20]. However, all of these existing origami-inspired designs were designed by hand. To display the usefulness of the proposed method, this work additionally introduces an obstacle that lies between the crease pattern and the point object that the final origami should be able to grip (Fig. 6(a)). In the robotic arm task, 11 points are approximated so that the tip of the arm follows a given trajectory (Fig. 6(b)), equivalent to the path generation synthesis of spatial mechanisms [65]. The remainder of the section first defines the input for both design tasks, and then in addition to the filter based on isomorphism, presents two additional filters to check for redundant graphs within the generation step. The section concludes with the optimization objectives definition, and the results produced by the method.

The design tasks of the gripper and the robotic arm are illustrated in Figs. 6(a) and 6(b), respectively. The initial graph G₀ (identical for both tasks) with vertices v₁ and v₂ end edge e_1,2 is shown in black and corresponds to the graph described in Fig. 3 with the same vertex labels $L1=((∅),∅)$ and L₂ = ((v₁), χ), edge label L_1,2 = (f₁, f₂), as well as in-plane coordinates $x1=(0,0)$ and $x2=(1,0)$⁠. The driving angle of the edge e_1,2 in both cases is linear, ρ_1,2 = t, where t goes from zero to t_max = π/2, applied mirror-symmetrically with respect to the xz-plane. The actuation thus corresponds to a single DOF, and the maximum number of internal vertices for both design tasks is set to N_max = 3. Note that graphs with one or two internal vertices are possible, N_max is just the upper limit.

The point object P_G to grip in Fig. 6(a) is located at P_G = (0, 0, 1) and the obstacle is represented by a cylinder with radius 0.25 whose axis is coincident with the line segment going from (0, − 1, 0.5) to (0, 1, 0.5). This setup allows for the generation of only one side of the gripper that emerges from the initial graph shown in black, after which the arm is rotated, once optimized, by $180deg$ around the z-axis, which is illustrated in Fig. 6(a) by the graph in gray corresponding to $v2′$⁠. Figure 6(b) shows the trajectory to approximate and the 11 points $PTi$ located at (1, 0, 0) + (2cos t_i, 0, 2 sin t_i) where t_i goes from 0 to π/2 in steps of π/20.

To facilitate the description of a crease pattern graph, here the following notation is introduced for the rule application sequences: when a vertex v_i is extended by r₁, i is added to the sequence, and when vertices v_i and v_j are combined by r₂, (i, j) is added to the sequence. As an example, extending vertices v₂ and v₃ in the given order and then combining vertices v₄ and v₈ yields a sequence denoted as 2, 3, (4, 8).

Since both design tasks illustrated in Fig. 6 are mirror symmetric with respect to the xz-plane, a graph can be filtered if it is mirror symmetric to another graph with respect to the x-axis. An example for the two symmetric graphs 2, 3 and 2, 5 is given in Fig. 7(a). These two graphs will result in the same (symmetric) optimized configuration, which renders one of the graphs redundant. In addition, any rule applied to, e.g., 2, 5 will generate a graph that is symmetric to a graph that can be generated by applying the corresponding rule to 2, 3. Thus, in both design tasks, the method filters all symmetric graphs and their descendants that result from applying more rules to symmetric graphs.

The second additional filter for both problems results from the task that the origami concept must perform. One of the generated vertices contained in the graphs has to approximate the point to grip P_G in the gripper task or the trajectory points $PTi$ in the robotic arm task. For simplicity, such a vertex is called the gripping vertex for both tasks. Extending a vertex by r₁ only adds value to the resulting graph, not its descendants, if one of its successors becomes the gripping vertex, which is why only the successors of the vertex extended last can perform the gripping task. Figure 7(b) shows the example of a crease pattern graph 2, 4, 3 that contains unnecessary vertices and edges, here called a semantically invalid graph [66]. Since v₃ is the vertex extended last, v₉, v₁₀, or v₁₁ are the gripping vertex candidates. In this case, extending v₄ is semantically invalid since all of its successors v₆, v₇, and v₈ can be cut from the graph without changing the performance of the origami gripper. In contrast to symmetric graphs, however, the descendants of semantically invalid graphs can lead to useful crease pattern graphs, which is why their descendants are not filtered.

An optimization of the robotic arm is considered successful if the objective value $Ω*$ of a design alternative is smaller than or equal to $ϑ=3*10−1$⁠. The variable boundaries for both tasks are l_min = 0.1, l_max = 1.5, where α_min = π/18 and α_min = 17π/18. The length boundaries prevent the generation of heavily distorted faces, whereas the sector angle boundaries contribute prevent sickle-like patterns or trivial folds for α = π, respectively.

With the initial graph G₀ in Fig. 6, N_max = 3, and rules r₁ and r₂, the graph grammar GG generates a total of 291 crease pattern graphs for both design tasks. After filtering isomorphic, symmetric, and semantically invalid graphs, 52 of the original 291 graphs remain. The design space of all possible crease pattern topologies and rule application sequences for both design tasks is visualized in the search tree in Fig. 8. The filled nodes illustrate the 52 meaningful graphs generated, whereas both the symmetric (sym.) as well as the semantically invalid (inv.) graphs are illustrated using white nodes. Each level of the search tree contains the graphs that exhibit the same number of applied rules starting from one application at the top to seven applications at the bottom.

From the 52 distinct graphs emerge 1170 distinct design alternatives (and as many objective function evaluations) that comprised the different assignments of RBMs and gripping vertices. While the topology generation of graphs is instantaneous, the optimization of a single design and the intersection check of a successful design takes from a few seconds to minutes depending on the complexity of the graph. The graphs corresponding to the first three levels in Fig. 8 each take between 0.5 and 30 s to optimize, graphs generated by four rule applications take about 10–60 s, and graphs with more rule applications show an increasing tendency where the maximum time recorded is 3 min for a single optimization. The intersection check averages at around 10 s for each design alternative. The total runtime involving the enumeration of the design space, the optimization within it, as well as the intersection check, takes 35 h for the gripper task and 40 h for the robotic arm task on an Intel i7 processor with 8 GB RAM.

With respect to the gripper, 836 of the 1170 design alternatives are able to grip the target point P_G, and thereof, 148 alternatives do not self-intersect. When the cylinder is present, the number of feasible origami concepts reduces to 36. The gripping motion of three such feasible origami concepts is depicted in Fig. 9 at discrete folding states from left to right. Figure 9(a) shows a gripper concept with a thin stem and big, arrow-like arms whose motion runs along the axis of the cylinder. The gripper in Fig. 9(b) occupies little space because most of its surface in the flat state lies beneath the cylinder, and its gripping motion is angled at approximately $45deg$ to the axis of the cylinder. The gripper in Fig. 9(c) is larger than the two other concepts, grips perpendicularly to the cylinder, and contains a degree-5 vertex. The details involving the optimized vertex locations for the designs in Fig. 9 can be found in the Supplementary Material (see Supplemental Fig. 3 available in the Supplemental Materials on the ASME Digital Collection and Supplemental Tables 4 and 5 available on the ASME Digital Collection).

With respect to the robotic arm task, 56 of the 1170 design alternatives are able to successfully approximate the trajectory points $PTi$⁠, and thereof, 29 alternatives do not self-intersect. The folding motion of three design alternatives is depicted in Fig. 10 from left to right together with the trajectory points $PTi$⁠. The three design alternatives exhibit an objective value of $Ω*=2*10−2$ (a), $Ω*=3*10−2$ (b), and $Ω*=4*10−2$ (c). The details involving the optimized vertex locations for the designs in Fig. 10 can be found in the Supplementary Material (see Supplemental Fig. 3 available in the Supplemental Materials on the ASME Digital Collection and Supplemental Tables 4 and 5 available on the ASME Digital Collection).

In comparison to related numerical approaches [67] that take multiple hours for the optimization of a single design alternative, the graph-based representation offered in this work enables the optimization of a design alternative in seconds to minutes. This is primarily due to the fact that the presented method incorporates the design knowledge required to conceptualize origami-adapted structures, involving the condition for rigid foldability, the kinematic model that determines the motion of all vertices depending on their RBM assignments, as well as the corresponding rule sets. The formalized knowledge allows to formulate the search space of crease pattern topologies, and subsequently exhaustively searches that space to generate rigidly foldable solution candidates. Such a generic approach to search space generation and exploration is not found in related work, e.g., [37,38,45], even including methods [40,41] that reside to graph-theoretical approaches for resolving some of the complexity of origami crease pattern generation. Due to the generality of the presented method, the developed concepts of the origami-adapted structures can satisfy three-dimensional design tasks, fold rigidly with the predefined number of DOF, and are free of self-intersection.

Figure 11(a) plots the distribution of generated design alternatives over the number of rule applications for the gripper design task. The distribution of the 836 designs that are able to grip P_G peaks at (A) with 311 designs and four rule applications, in contrast to the total distribution of the 1180 designs having a peak with 432 designs and five rule applications at (B). Figure 8 shows that at four rule applications, all graphs except one are generated by three applications of r₁ and one application of r₂, while only r₂ can be applied thereafter. Rule r₂ introduces constraints into the system without adding variables, which is why the number of designs that can grip P_G steadily decreases from four applications onward. For the same reason, the percentage of successful gripping designs per total number of designs is highest at three rule applications (C divided by D), where five of seven graphs are generated purely by r₁ (Fig. 8, third level of the search tree). The distributions for the 148 designs without self-intersection and the 36 designs without any intersection peak at points (E) and (F) with 79 and 14 designs, respectively, showing a behavior similar to the distribution of the 836 designs that are able to grip P_G.

When at first one could have argued not to apply r₂ since it only introduces constraints, the increase in feasible designs from three to four rule applications show that a certain number of applications of r₂ can drive the optimized designs out of intersecting configurations. Without specific constraints that prevent intersection within the optimization, r₂ thus provides more variability in the crease pattern graphs that satisfy the design task. However, if r₂ is applied too many times, the designs are either unable to grip P_G or lead to intersection, which leads to zero feasible designs at six (G) and seven rule applications for the gripper task.

Figure 11(b) refers to the robotic arm task and shows the distribution of the 56 successful designs that are able to approximate the given trajectory and the 29 designs that in addition do not self-intersect. The peaks appear at four rule applications (H and J), which again shows the usefulness of r₂. In contrast to the results of the gripper, for the robotic arm task there are feasible concepts generated by six rule applications (K), which stems from the problem formulation (Fig. 10). Although the given trajectory is a perfect quarter circle that seems straightforward to approximate, the maximum crease line length l_max = 1.5 prevents a simple solution, allowing for complex concepts to be generated. The target criteria for both the gripper and the robotic arm task listed in Sec. 4.3 are determined empirically and would slightly change the plots in Fig. 11 if chosen differently. However, since the method uses a global optimization, the only influence of the target criteria on the results are the number of feasible design alternatives, leaving a choice for the user between precision and variability of the design alternatives.

Applied to both the gripper and the robotic arm design task, the enumeration of the design space yields 52 distinct crease pattern graph topologies depicted in Fig. 8. While the search tree is narrow at the top where r₁ predominates, the midsection of the search tree becomes broader since more extended vertices offer more possibilities for the application of r₂. The number of possible rule applications then progressively decreases with each line because of the maximum number of extended vertices N_max = 3, which is why the search tree narrows again at the bottom. However, even with a small increase in N_max to 8 vertices, for example, the maximal number of RBMs for only one particular topology is 2⁸ making the enumeration of all crease patterns for their respective RBMs impractical. By adding to the search space all variants of crease pattern graphs that can be generated for certain limiting N_max, the enumeration influences the scalability of the method significantly and blocks more complex engineering applications. For example, for N_max = 8 the design space yields 3123 distinct graphs and 147,078 distinct design alternatives that needed to be optimized in the subsequent step. To improve the scalability of the method in the future, a more informed rule application should be performed. In both design tasks, there is no statistical difference in the performance of design alternatives with respect to the assignment of RBMs; both modes M_↑ and M_↓ are equally represented across the range. Considering the search method, a transfer from the design space enumeration to more efficient search methods, such as branch-and-bound algorithms [68], is required to expand the application of the method to larger design spaces.

Another interesting path for future work involves the expansion of the presented graph grammar. The rule set $R=(r1,r2)$ enables the generative design of a variety of novel crease patterns such as the ones shown in Figs. 9 and 10, as well as existing crease patterns such as slender origami [59], origami strings [60], or the gripper in Ref. [61]. However, many known crease patterns make use of specific sector angle configurations and mode assignments that fold rigidly although they cannot be obtained by using the graph grammar in this work. To generate such crease patterns, the rule system presented could be expanded by an adjusted version of rule r₂, here denoted as $r2′$⁠. Rule $r2′$ would also combine two vertices, but the vertex $vc1′$ on the RHS₂ in Fig. 5 would lie on the line segment between the predecessors $vp1$ and $vp2$ and would become a sink vertex (⁠$deg+(vc1′)=0$⁠) of the type $Tc1′=∅$⁠. Figure 12 illustrates how a Miura-ori pattern could be generated by the expansion of the presented rule system with rule $r2′$⁠.

The crease pattern graph depicted in Fig. 12(a) can be generated by the initial graph G₀ in Fig. 3 with the rule application sequence 2, 3, 4, 8. By applying $r2′$ to vertices v₉ and v₁₄, resulting in the rule application sequence 2, 3, 4, 8, (9, 14)′ and the crease pattern graph shown in Fig. 12(b) that depicts a Miura-ori unit cell. The generation of the pattern in Fig. 12(b) would require a condition for the dihedral angles, $ρ8,9′=ρ4,9′$⁠, and thus demand a condition for the sector angle configuration of the entire crease pattern. Such a rigidly foldable crease pattern can only be achieved by introducing some sort of symmetry [69] or by complying with the fold angle multipliers [70]. Introducing such conditions within the constraint system generated by the graph grammar is complex and requires future work.

Finally, the concepts developed in this paper are of zero thickness and applying them in real-life scenarios will require adaptations toward finitely thick materials. Adapting zero thickness concepts works best for crease patterns that exhibit only degree-4 vertices and uniform sector angle configurations. To achieve the latter, the variable boundaries α_min and α_max in Eqs. (7) and (10) can be defined to lie within a small interval. However, introducing such intervals might deplete the number of feasible design concepts, and the user of the method has to estimate the trade-off between the number and the realizability of these solutions. To achieve the former, users can select the feasible designs generated only by applications of rule r₁ or a minimum number of applications of rule r₂. Nonetheless, adapting the feasible concepts to finite thickness still requires future work, especially if the adaptations were to be integrated into the method itself. A range of possible adaptation techniques is given in Ref. [51].

In this paper, we introduced a generative method for the generation of novel origami concepts for engineering design tasks. The method includes an origami graph grammar that generates rigidly foldable crease patterns as well as the constraints for both the boundaries of the optimization variables based on the conditions for rigid foldability. Although the method generates rather simple patterns, the enumeration of all possible design alternatives based on the RBMs to the internal vertices shows a surprising richness of the underlying search space resulting in a considerable number of design concepts that have been generated. Finally, the paper demonstrates that the PTU enables the efficient optimization of design alternatives even including the self-intersection check. The future work will involve a development of a more efficient PTU-based generative method that does not rely on the exhaustive search and is able to generate origami patterns of increased complexity, involving origami strips, tessellations, and nondevelopable patterns.

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. The authors attest that all data for this study are included in the paper.

f =

facet

l =

crease line length

n =

degree of a vertex

r =

graph grammar rule

t =

time parameter

u =

set of units

v =

vertex

x =

set of in-plane vertex coordinates

$d$ =

normalized direction vector

$x$ =

in-plane vertex coordinates

E =

set of edges

G =

graph

M =

generic rigid body mode

N =

number of internal vertices

R =

set of global rotation matrices

V =

set of vertices

X =

set of three-dimensional vertex coordinates

$X$ =

three-dimensional vertex coordinates

$G$ =

set of graphs

$L$ =

language of a grammar

$M$ =

left-hand-side match for graph grammar rules

$R$ =

set of graph grammar rules

e_i,j =

edge from vertex v_i to vertex v_j

u_i =

unit, i = 1, 2, 3

$fR$ =

function determining local rotation matrices

f_φ =

function determining unknown dihedral angles

L_i =

vertex label of vertex v_i

L_i,j =

edge label of edge e_i,j

L_E =

set of edge labels

L_V =

set of vertex labels

P_i =

predecessors of a vertex v_i

T_i =

type of a vertex

U_i =

unit angle, i = 1, 2, 3

$Ri$ =

local rotation matrix, i = 1, 2, 3

$Ri,j$ =

global rotation matrix of the edge e_i,j

$FLi,j,FRi,j$ =

facets on the left and right of edge e_i,j, respectively

M_↑, M_↓ =

rigid body modes up and down, respectively

P_G, P_T =

point to grip and points of trajectory, respectively

α =

sector angle

ρ =

set of dihedral angles

ρ_i,j =

dihedral angles

$ϑ$ =

target criterion

Σ_E, Σ_V =

maps for edge and vertex labels

φ_i =

unknown dihedral angle, i = 1, 2, 3

Φ =

set of optimization variables

χ =

symbol for type T meaning vertex can be extended

ψ =

set of optimization constraints

Ω =

objective function

$∅$ =

empty and terminal symbol

$⊎$ =

disjoint union of sets