## Abstract

Elastic gridshell is a class of net-like structure formed by an ensemble of elastically deforming rods coupled through joints, such that the structure can cover large areas with low self-weight and allow for a variety of aesthetic configurations. Gridshells, also known as X-shells or Cosserat Nets, are a planar grid of elastic rods in its undeformed configuration. The end points of the rods are constrained and positioned on a closed curve—the final boundary—to actuate the structure into a 3D shape. Here, we report a discrete differential geometry-based numerical framework to study the geometrically nonlinear deformation of gridshell structures, accounting for non-trivial bending-twisting coupling at the joints. The form-finding problem of obtaining the undeformed planar configuration given the target convex 3D topology is then investigated. For the forward (2D to 3D) physically based simulation, we decompose the gridshell structure into multiple one-dimensional elastic rods and simulate their deformation by the well-established discrete elastic rods (DER) algorithm. A simple penalty energy between rods and linkages is used to simulate the coupling between two rods at the joints. For the inverse problem associated with form-finding (3D to 2D), we introduce a contact-based algorithm between the elastic gridshell and a rigid 3D surface, where the rigid surface describes the target shape of the gridshell upon actuation. This technique removes the need of several forward simulations associated with conventional optimization algorithms and provides a direct solution to the inverse problem. Several examples—hemispherical cap, paraboloid, and hemi-ellipsoid—are used to show the effectiveness of the inverse design process.

## 1 Introduction

Traditional three-dimensional shell structures can resist external loads through their inherent shapes; however, if regular holes are made in the shell, with the removed material concentrated into the remaining strips, a structurally flexible gridshell can be achieved. Several spectacular architectures, e.g., Helsinki Zoo’s observatory tower and Centre Pompidou Metz, were manufactured with a network of one-dimensional beams, such structures serve both aesthetical and functional purposes in the civil engineering community [1]. Besides the construction of buildings in civil engineering [24], abundant applications in mechanical systems, e.g., micro/nanostructures [57], stretchable electronics [8,9], and bio-inspired patterns [10,11], employ gridshell as a major structural component in their design step to achieve specific functionalities. While the gridshells studied by Baek et al. [12,13] had joints that were free to rotate and twist, recent work by Panetta et al. [14] constrained the bending and twisting at the joints. This leads to non-trivial twisting and bending coupling between two rods at the joints, which can improve the robustness of the structure and increase the design space of the architectural shapes [3]. Computationally efficient numerical simulation tools for this class of structures can allow simulation-guided design and eliminate the need for painstaking trial-and-error prototyping.

In the computational mechanics community, modeling and simulation of thin elastic objects, e.g., rods and shells, are of sufficient general interest because of the preponderance of geometrically nonlinear deformation. Finite element method has been the most commonly used method in structural analysis over the past few decades [15,16]. Recently, discrete differential geometry (DDG)-based methods [17] are becoming popular in the computer graphics community to simulate the thin elastic structures, e.g., hair and clothes, due to the computational efficiency and the robustness in handling geometric nonlinearity, collision, and contact. Previous DDG-based methods have shown surprisingly successful performance in simulating slender structures, e.g., rods [1822], beams [23], beam networks [24], ribbons [25,26], and plates/shells [26,2730]. Gridshell, on the other hand, usually represents a curved surface comprised of multiple 1D elastic rods and differs from the traditional 1D rods or 2D shells. This leaves room for new numerical methods for accurate and efficient simulation of gridshells. Baek et al. first proposed a method based on discrete elastic rods (DER) to investigate the buckling instability and form-finding of gridshells [12] and found excellent agreement between experiments and simulations. A stiff spring is used in that framework to simulate the joint between two rods and the spring force is treated using an explicit approach. The joint between two intersecting rods is free to twist as well as rotate such that the twisting and bending coupling between two rods are not taken into account. This numerical framework was later used to study the elastic rigidity of hemispherical gridshells [13]. Numerical methods to capture the bending and twisting coupling at the joints either use a penalty energy between the neighboring material frames of two rods system [31,32], or a geometric constraint-based energy functional [14]. Finite element-based numerical methods to simulate this class of structures have also been introduced [33,34].

An even more intriguing feature of elastic gridshell is its form-finding process. Figure 1(a) shows three examples of convex 3D gridshells, whereas Fig. 1(b) describes the actuation process. In Fig. 1(b), the undeformed gridshell is planar and the extremities of the elastic rods fall on a closed curve, $G0$. In order to actuate the gridshell, the end points of the rods are constrained to fall on a second closed curve, $G$. The form-finding problem, i.e., the inverse problem in this case, calls for computation of $G0$ given the target 3D shape and the final boundary, $G$. This transformation between the 2D planar structures and the complex 3D topologies by using the geometry and structural instability is of interest [5,6,11,34,35] and might lead many applications in mechanical systems [3648]. Prior works on mechanically guided assembly of 3D structures range from macroscopic origami-inspired structures [49,50] to microscopic buckling of elastic ribbons attached to a pre-stretched substrate [6]. While a number of studies investigated the forward dynamics, we focus on a computationally efficient method to solve the inverse design problem of finding the initial planar shape with a given 3D target configuration. Reference [51] considered Chebyshev net theory to map a group of rods onto a given surface, e.g., human face, to design wire mesh. Prior works on the inverse problem include analytical solution to a pair of ordinary differential equations (ODEs) on the basis of Gauss equation [12], or numerical optimization coupled with physics-based simulations [11,14,30]. Recently, a genetic algorithm-based method [30] and an optimization-based simulation framework [14] have been introduced to study the form-finding problem in elastic gridshells; however, these methods require running the physics-based simulation numerous times in order to find the optimal solution, especially when a good initial guess is not available. As an example, the form-finding problem of a hemispherical gridshell in Fig. 1(a) may take approximately 102 generations with 5 × 102 individuals in a population in genetic algorithm, corresponding to 5 × 104 forward simulations. The proposed method reduces this problem to a single forward simulation.

Fig. 1
Fig. 1
Close modal

Here, we develop a numerical method for the inverse form-finding problem of gridshells. Different from above analytic and optimization methods that typically require numerous “forward” simulations to predict the deformation of the gridshell under various boundary conditions imposed by $G$ and $G0$, this method implements a mechanics-based forward simulation of a gridshell draping around the target rigid shape under gravity. This single forward simulation can offer an excellent solution to the inverse problem. The simulation relies on a DER-based numerical framework, where both the rods and the joints are represented by the discrete elastic rod model. Discrete equations of motions, based on the balance of elastic and external forces, are solved to update the structural configuration with time. The main contribution of this paper is a numerical method for form-finding of convex gridshells based on contact [52]. In Fig. 1, we show several 3D configurations of convex gridshell structures as well as their corresponding initial planar boundaries constructed by the contact-based method described in this paper. The boundary of the 2D undeformed shape, $G0$, can be almost exactly obtained by draping the elastic gridshell under gravity over the rigid 3D target surface. This calls for simulation of contact between the gridshell and the rigid surface and is handled via the modified mass method [27,53]. Discrete simulations are naturally suited to handle contact, which underlines the need for DDG-based methods in the study of form-finding of gridshells. The initial planar pattern of grid can be easily obtained by only running the physically based simulation once, which can significantly reduce the computational time when solving the form-finding problem.

Our paper is organized as follows. In Sec. 2, we discuss the DER-based numerical framework of gridshell simulation, with a focus on geometric decomposition and bending-twisting coupling at rotational joints. Then, in Sec. 3, we introduce a modified mass-based contact algorithm for the form-finding problem associated with gridshells. Finally, concluding remarks and avenues for future research are presented in Sec. 4.

## 2 Numerical Method

In this section, we discuss the forward physically based simulation of gridshells. Gridshell is a type of structure that comprises multiple one-dimensional rods connected through joints. These joints may twist and rotate [14]. Here, in Sec. 2.1, we first briefly review the DER method for a single elastic rod, then extend this method to a numerical framework for the simulation of gridshells in Sec. 2.2. Finally, Sec. 2.3 presents two simple cases to demonstrate the twisting and bending coupling between two rods at the joint.

### 2.1 Discrete Elastic Rods Method.

In the discrete setting of DER, shown schematically in Fig. 2(a), the rod centerline is discretized into N nodes: [x1, x2, …, xN], and (N − 1) edges: [e1, e2, …, eN−1], with ei = xi+1xi, where i = 1, …, N − 1. Each edge, ei, has an orthonormal adapted reference frame ${di(1),di(2),ti}$ and a material frame ${mi(1),mi(2),ti}$; both the frames share the tangent ti = ei/‖ei‖ as one of the directors. The reference frame is updated at each time-step through parallel transport in time, and, referring to Fig. 2(b), the material frame can be obtained from a scalar twist angle θi, see Ref. [21] for a detailed exposition of the DER algorithm. Nodal positions (total 3N) and the twist angles (total N − 1) constitute the (4N − 1)-sized degrees-of-freedom (DOF) vector, q = [x1, θ1, x2, …, xN−1, θN−1, xN], of the discrete rod. Based on this kinematic representation, in the remainder of this section, we discuss the formulation of elastic energies, elastic forces, and the time stepping procedure of the rod solver.

Fig. 2
Fig. 2
Close modal
An elastic rod is modeled as a mass-spring system, with a lumped mass (and angular mass) at each node (and edge), and associated discrete stretching, bending, and twisting energies. For a rod with Young’s modulus E, shear modulus G, and isotropic circular cross section, the elastic energies—stretching, bending, and twisting—are given by [18,19]
$Es=12∑i=1N−1EA(ϵi)2‖e¯i‖$
(1a)
$Eb=12∑i=1NEIΔli(κi−κi¯)2$
(1b)
$Et=12∑i=1NGJΔli(τi)2$
(1c)
where A is the area of the cross section, I is the area moment of inertia, J is the polar moment of inertia, εi is the stretching strain associated with the ith edge, $‖e¯i‖$ is its undeformed length, $κi=[κi(1),κi(2)]$ is the bending curvature at the ith node ($κ¯i$ is the curvature in the undeformed configuration), τi = θi − θi−1 + mi is the twist at the ith node (mi is the twist of the reference frame, details can be found in Refs. [19,21]), and Δli = (‖ei‖ + ‖ei+1‖)/2 is its Voronoi length. The strain measures, i.e., εi, $κi$, and τi, can be expressed in terms of q (specifically, xi−1, θi−1, xi, θi, xi+1). The case of non-circular cross-section can be included in the above formulation with minor changes [14,18,19].
At each degree-of-freedom, qj (jth element of the vector q), the elastic forces (associated with nodal positions) and elastic moments (associated with the twist angles) are given by
$Fjint=−∂∂qj(Es+Eb+Et)$
(2)
where j is an integer between 1 to 4N − 1 and $Fjint$ is the jth element of the (4N − 1)-size elastic force vector.
In the time stepping scheme of the DER method, implicit Euler integration is used to solve the following 4N − 1 equations of motion and update the DOF vector q as well as its velocity (time derivative of DOF) $v=q˙$ from time-step tk to tk+1 = tk + h (h is the time-step size):
$MΔq(tk+1)−hMv(tk)−h2[Fint(tk+1)+Fext(tk+1)]=0$
(3a)
$q(tk+1)=q(tk)+Δq(tk+1)$
(3b)
$v(tk+1)=1hΔq(tk+1)$
(3c)
where Fext is the external force vector (e.g., gravity, damping, or reaction force due to contact) and M is the diagonal mass matrix comprised of the lumped masses [54]. The Jacobian associated with Eq. (3a) necessary for Newton’s iteration can be expressed as
$Jij=miδij−h2(∂Fiint∂qj+∂Fiext∂qj)$
(4)
where i and j are integers between 1 to 4N − 1, the mass (or angular mass) associated with the ith DOF is mi, and both the elastic internal force, $Fiint$, and external force, $Fiext$, are evaluated at t = tk+1. If the gradient of the external forces cannot be evaluated, $∂Fiext/∂qj$ can be ignored so that these forces are handled explicitly. Importantly, the Jacobian $J$ is a banded matrix and the time complexity of this algorithm is O(N) (the computational time linearly scales with the number of nodes [19]). This computational efficiency has motivated its application in the animation industry (e.g., hair simulation for movies) as well as its adoption in mechanical engineering [5557].

### 2.2 Discrete Elastic Gridshells.

Now we formulate the discrete numerical framework for simulation of gridshells. We consider a basic element of the gridshell in Fig. 3(a), where two rods intersect each other at a shared point. The position of the jth node on the ith rod within the gridshell system is denoted as xi,j. The twist angle of the jth edge on the same rod is θi,j. In Fig. 3(a), the two nodes, x1,3 and x2,3, from two different rods overlap at the joint. A straightforward method to enforce the coincidence of two nodes at the joint is a linear spring-like energy of the form
$Ec=12C‖x1,3−x2,3‖2$
(5)
where C is the Lagrange multiplier. Its negative gradient, −(∂Ec/∂q), is included as an external force in Eq. (3a). The Hessian of this energy can be trivially computed to aid the Newton’s method in the solution of the equations of motion. An alternative approach to enforce this condition will be discussed next.
Fig. 3
Fig. 3
Close modal

In addition to the coincidence of two nodes, there is a non-trivial coupling between the twisting and bending modes at the joints, e.g., twisting rod 1 in Fig. 3(a) can cause rod 2 to rotate. Here, we consider the pin-joints with a specific constraint for rotations at the contact area. To account for this coupling at the joints, we decompose the basic gridshell element into four elastic rods in Fig. 3(b): the first two are the physical rods denoted as rod 1 and 2; the other two are linker rods with 3 nodes to model the joints. Hereafter, we use subscripts to denote quantities associated with the physical rods, e.g.,x1,1 is the first node on the first rod, and superscripts when associated with linker rods, e.g., x1,1 is the first node on Linker 1. Each rod can be simulated by the conventional DER method. A penalty energy can be used to account for the coupling between twisting and rotating at the joints. For the first linker and the first physical rod, the penalty energy is

$Ec′=12C1‖x1,2−x1,1‖2+12C1‖x1,3−x1,2‖2+12C2‖θ1,2−θ1,1‖2$
(6)
where C1 and C2 represent the stiffness of the joint against the rotation and twist coupling. A similar penalty energy exists between the first linker and the second physical rod. At sufficiently high values of C1 and C2, the rods at the joint cannot twist or rotate with respect to one another. We use C1 = C2 = 106EI in the current numerical investigation after a convergent study [30]. The external force and Jacobian associated with these energies can again be trivially computed. We should keep in mind that, when we decompose the basic element of the gridshell structure into two rods and two linkers, the mass and stiffness of the rods at the joint should not be double counted, e.g. the lumped mass at the joint node should be divided by four and then used as the mass associated with x1,3, x2,3, x1,2, and x2,2.

In our numerical implementation, at every time-step, the equations of motions for the physical and linker rods are independently solved. This allows us to take advantage of the banded nature of the Jacobian matrix. The penalty forces in Eqs. (5)(6) are then calculated and included as external force in the next time-step, i.e., the penalty forces are treated explicitly. An alternative to this approach of solving a number of smaller systems and subsequently bringing them together is to solve a large system, consisting of all the physical and linker rods, with an implicit treatment of the penalty forces. The large system would no longer have a banded Jacobian matrix since the Hessian matrix of the penalty energies would occupy non-banded entries within the Jacobian. A second alternative is to forego the use of penalty energies and treat the overlapping nodes (e.g. x1,3, x2,3, x1,2, and x2,2) and edges (e.g. θ1,2 and θ1,1) with the same degrees-of-freedom. For example, instead of using 3 × 4 degrees-of-freedom for the overlapping nodes in Fig. 3(b), we can introduce three degrees-of-freedom, xjoint, for the joint node and apply the sum of forces from all the four nodes onto the newly introduced single node. A simulation code developed for one method can be easily re-purposed to employ a different method. While solving extremely large systems, correct choice of the time integration scheme may depend on the computer memory as well as the degree of parallelism. A detailed comparison among the explicit method (used in the current study), implicit method, and the mapping method can be found in Ref. [30].

### 2.3 Demonstration of Bending and Twisting Coupling.

We use two simple demonstrations to show the coupling between two rods at the joint. In Fig. 4(a), we show the response of the basic element of a gridshell when one rod is twisted. The first twist angle of one rod, θ1,1, referring to Fig. 4(a), is rotated with a prescribed angular velocity, ω = 10 rpm, such that θ1,1(t) = ωt. Due to the two linkers at the joint, the centerline of the other rod rotates about the former one with the same angular velocity ω. This demonstrates the twisting coupling between two rods.

Fig. 4
Fig. 4
Close modal
We now turn to bending coupling in Fig. 4(b) and specify the position of the first node of the first rod, x1,1 ≡ [x1,1 (t), y1,1 (t), z1,1 (t)], highlighted in Fig. 4(b), where
$x1,1(t)=12Lcosωt$
(7a)
$y1,1(t)=12Lsinωt$
(7b)
$z1,1(t)=0$
(7c)
L is the total length of the undeformed rod, and ω = 10 rpm. The location of the middle node on the first rod that falls on the joint is kept fixed with time to avoid rigid body motion. Due to these fixed degrees-of-freedom, the first rod rotates about the z-axis with a prescribed angular velocity ω. Also, because of the coupling of bending between two rods at the joint, the second rod also rotates about the z-axis at the same angular velocity.

## 3 Form-Finding

We now focus on the inverse problem of obtaining the undeformed 2D shape, given the target 3D shape. This method relies on draping the gridshell around a rigid body with the target geometry. In this section, first, the discrete gridshell algorithm is coupled with the modified mass method [27,53] to simulate the draping process; then, the procedure to obtain the initial boundary, $G0$ in Fig. 1(a), is described, accompanied by a number of examples.

### 3.1 Modified Mass Method.

In Fig. 5(a), the target rigid shape is described by the function z = f(x, y). The position of a node, xi(tk), in a discrete gridshell structure at time t = tk, approaches the target rigid surface, z = f(x, y). If the rigid surface is not accounted for, the position of this node at the next time-step is, say, xi′(tk+1) ≡ [x0, y0, z0]. In the time marching scheme of the simulation, if this node falls under the target surface so that
$z0
(8)
a correction is required to move xi′(tk+1) onto the target surface along the the surface normal vector
$pn(x0,y0)=[−∂f(x0,y0)∂x0,−∂f(x0,y0)∂y0,1]$
(9)
and the necessary displacement vector for correction is
$δi=‖δiz‖cosψpn‖pn‖$
(10)
where $δiz=[f(x0,y0)−z0]nxy$, ‖pn‖ is the magnitude of the vector pn, and ψ is the angle between nxy (the normal vector to the xy plane), and pn (the surface normal vector), as shown in Fig. 5(b).
Fig. 5
Fig. 5
Close modal
In order to enforce this displacement, the equations of motion for the nodal positions in Eq. (3) have to be slightly modified. For the three degrees-of-freedom at the ith node, xi = [xi, yi, zi], the updated form of Eq. (3a) is
$Δx˙i(tk+1)−hmiSi(tk+1)[Fiint(tk+1)+Fiext(tk+1)]−Δvipre(tk+1)=0$
(11)
where mi is the lumped mass; $Δvipre$ is the prescribed change in velocity that can be obtained from the prescribed displacement, $δi$; $Fiint$ is the three-element elastic force vector on the ith node; $Fiext$ is the external force vector; and the modified mass matrix is
$Si={IiffreeDOFofithnode=3(I−pnpnT‖pn‖2)iffreeDOFofithnode=2(I−pnpnT‖pn‖2−qnqnT‖qn‖2)iffreeDOFofithnode=10iffreeDOFofithnode=0$
(12)
where $I$ is the 3 × 3 identity matrix; pn is the constrained direction when free DOF=2; and pn and qn are the constrained directions when free DOF=1. Note that when a node is free, $Δvipre=0$, and Eq. (12) can thus be simplified to Eq. (3). If the node is fully constrained ($Si=0$), Eq. (12) reduces to $Δx˙i(tk+1)−Δvipre(tk+1)=0$ and the change in position (as well as the velocity) is enforced to take the prescribed value. In our case, we only constrain the node along surface normal, pn, such that the number of free DOFs of the ith node is 2. In our numerical implementation, we employ inelastic collision between the i-th node and the target 3D surface, i.e. once the node is in contact with the target surface, we manually set its velocity to zero at the end of the current time-step.

Every time-step in simulation accounting for the contact with a rigid surface may require integration of the equations of motion twice. The first solve is the predictor step that determines if any node fell under the target surface. The optional second solve is the corrector step that is only necessary if any node was detected to fall through the rigid surface. A pseudo code of contact-based form-finding simulation is provided in Algorithm 1.

#### Modified mass-based contact simulation of elastic gridshell

Algorithm 1

Require:m — Number of rods

Require:n — Number of linkers

Require:$q(i)(0)$, where $i∈[1,m]$ — Initial DOFs of rods

Require:$q(i)(0)$, where $i∈[1,n]$ — Initial DOFs of linkers

Require:h — time step size

Require:T — total simulation time

Ensure:$q(i)(T)$, where $i∈[1,m]$ — Final DOFs of rods

Ensure:$q(i)(T)$, where $i∈[1,n]$ — Final DOFs of linkers

$k←0$

$tk←0.0$

while$tk≤T$do

$tk+1=tk+h$

$k←k+1$

for$i=1$ to $i=m$do

$solved←0$

while solved = 0 do

Update DOFs in i-th rod ($q(i)$)

$solved←1$

for each node on i-th rod do

$[x0,y0,z0]←$$x,y,z$ coordinates of the node

if$z0then

Constrain the normal direction of $xi,j$

$solved←0$

end if

end for

end while

end for

fori = 1 to i = ndo

Update DOFs in i-th linker (q(i)) based on Eq. (3b)

end for

end while

return q(i) (T), where $i∈[1,m]$.

return q(i) (T), where $i∈[1,n]$.

### 3.2 Initial Boundary From the Draping Method.

Three target shapes—hemisphere, paraboloid, and hemi-ellipsoid shapes—are used as examples to demonstrate the form-finding process of gridshell structures. The physical parameters are as follows: rod length s = 1.2 m, rod radius r0 = 1 mm (and, therefore, second moment of inertia, $I=πr04/4$, polar moment of inertia, $J=πr04/2$, and cross-sectional area $A=πr02$), Young’s modulus E = 1 MPa, shear modulus G = E/3 (assuming incompressible material), material density ρ = 1.0 g/cm3, distance between two parallel rods Δs = 3 cm, discrete edge length $‖e¯‖=5mm$, and distance between planar X-shell and 3D surface top (choosing arbitrarily) is H = 5 cm. As long as the rod can be assumed to be soft enough and inextensible, these parameters do not significantly influence the actuated shape of the gridshell [12]. The geometries of the target surfaces are given by
$Γ1(x,y,z):x2Rh2+y2Rh2+z2Rh2=1,with0⩽x2+y2⩽Rh2andz⩾0$
(13a)
$Γ2(x,y,z):x2Rp2+y2Rp2+zHp=1,with0⩽x2+y2⩽Rp2andz⩾0$
(13b)
$Γ3(x,y,z):x2ae2+y2be2+z2ce2=1,with0⩽x2ae2+y2be2⩽1andz⩾0$
(13c)
where Rh = 0.2 m (for hemisphere); Rp = 0.2 m and Hp = 0.12 m (for paraboloid); ae = 0.2 m, be = 0.15 m, and ce = 0.12 m (for hemi-ellipsoid).

In Figs. 6(a1)6(c1), the undeformed planar gridshells are located above the 3D rigid surfaces described by Eqs. (13). The elastic rods are symmetrically distributed about the x and y-axes in case of the hemisphere (17 × 17 grid) and the paraboloid (15 × 15 grid); for the hemi-ellipsoid, on the other hand, there are 11 rods along the x-axis and 15 along the y-axis. Note that the rod number for each case is determined by the size of the desired shapes, i.e. we want at least one node on each rod to contact the target surface. The planar gridshells are dropped under a gravity-type load that is large enough to drape the structure around the target rigid surface. In Fig. 6, gravitational acceleration of g = 9.81 m/s2 was sufficient. Figures 6(a2)6(c2) shows the deformed shapes of the gridshells. Parts of the gridshell are in contact with the rigid surface (located above the xy plane) and the other parts remain suspended under gravity below the xy plane. The suspended parts (i.e. nodes that fall below the minimum z-coordinate of the target rigid shape) are trimmed in Figs. 6(a3)6(c3) to obtain the new extremities (first and last nodes) on each elastic rod. This describes the final boundary $G$ of the form-finding problem (also see Fig. 1). In Figs. 6(a4)6(d4), the extremities upon trimming are mapped back to the initial planar gridshell, i.e., the planar shape is also trimmed to get rid of the suspended portions. This gives the initial boundary, $G0$, of the gridshell. Then, the target 3D pattern described in Eqs. (13) can be obtained by moving the nodes on the extremities of the rods from the initial footprint, $G0$, to the final boundary, $G$.

Fig. 6
Fig. 6
Close modal

The analytical solution to the initial boundary in case of a hemisphere [12] is also shown in Fig. 6(a4). For the cases of paraboloid and hemi-ellipsoid, the analytical solutions are not easy to derive and, therefore, we compare the planar boundaries obtained from the draping process and the ones found by genetic algorithm-based optimization [30] in Fig. 6(b4) and 6(c4). The good match indicates the correctness and the validity of the proposed method. Even when the solution from the process outlined in Fig. 6 is not accurate enough, it provides an excellent initial guess for conventional optimization algorithms.

For a physical understanding of this method, we consider the balance of forces. Each node in the simulation is balanced by three forces: (1) gravity, (2) contact force from the target rigid surface, and (3) elastic forces (primarily bending). This competition of forces yields a deformed shape that conforms to the target surface. On the other hand, in the “pop-up” fabrication process [12] of gridshell where the nodes on $G0$ are moved to $G$, gravity and contact forces are replaced by forces acting on the extreme nodes (located on the boundary) by an external agent. Our results show that, surprisingly, the deformed shape remains almost the same despite substitution of gravity and contact with boundary conditions on a handful of nodes. Next, to demonstrate that the numerical method is robust against initial grid spacing, in Figs. 7(a) and 7(b), we show the hemispherical gridshell with different grid spacings. Here, the distance between two parallel rods are Δs = 4 cm (for Fig. 7(a)) and Δs = 5 cm (for Fig. 7(b)); the rod number is changed to 13 × 13 and 11 × 11, respectively, to ensure that each rod comes in contact with the target surface. As shown in Figs. 7(a2) and 7(b2), the initial planar grids match well with the analytical solution in both of these two cases.

Fig. 7
Fig. 7
Close modal

### 3.3 Computational Time.

Next, we highlight the computational efficiency of the presented contact-based numerical simulation of elastic gridshells. In Fig. 8(a), we plot the computational time as a function of time-step size, h, for three different cases. Here, the number of DOFs in all scenarios is fixed at ∼8000. The total simulation times are approximately 50 s (for hemisphere), 30 s (for paraboloid), and 20 s (for hemi-ellipsoid), separately. Then, referring to Fig. 8(b), we show the reliance of the computational time on the number of DOFs. The time-step size in this figure is set to be h = 0.1 s. Unsurprisingly, computational time dramatically increases as the number of DOFs increases. The simulations are performed on a single thread of Intel Core i7-6600U Processor @ 3.4 GHz. Overall, reasonable predictions can be obtained within one minute.

Fig. 8
Fig. 8
Close modal

### 3.4 Limitations.

Via comparison with Refs. [12,30], the effectiveness of our proposed direct-contact method can be indirectly validated. The proposed method relies on the approximation that replacing the reaction forces at the boundary of the gridshell with a gravity-like load on all the nodes and constraints from the rigid target surface results in the same deformed shape of the gridshell. In the future, limitations of this approximation can be analyzed using theoretical mechanics. We also limited ourselves to convex surfaces with analytical solutions. Future directions of work can include surfaces with negative curvature and concave shapes; our contact-based form-finding method fails to handle negative curvature, e.g. hyperbolic surface. In this case, the final footprints are no longer planar. In addition, the proposed frictionless contact-based method would fail to fully drape around non-convex geometry. One possible approach in these cases is to use the proposed method to obtain an initial guess and then use this guess to start a more comprehensive inverse design process, e.g., genetic algorithm [30] and generative adversarial networks [58]. Another potential direction is to design multiple simple gridshells and then sew the solutions together to achieve more complex (potentially concave) gridshells. This concept was alluded to in Ref. [12].

The amount of twist along the rods in the examples explored in this manuscript is negligible. We did not find any significant effect of the coupling between bending and twisting on multiple rods at the joints. In the future, the role of the mechanics of the joints on the overall shape of the gridshell can be explored.

## 4 Conclusions

We introduced a numerical framework for the simulation of gridshells and solved the form-finding problem directly, without any numerical optimization. For the forward physical simulation, we first decomposed the gridshell as well as its joints into multiple elastic rods, such that each component can be treated using the well-established DER method. For the inverse problem of form-finding, we formulated a modified version of the discrete gridshell simulation algorithm by coupling it with the modified mass method to account for the contact between an elastic gridshell and the target rigid 3D surface. We showed that the gridshell, upon draping around the target shape, can be simply trimmed to directly get the initial planar boundary. A good match between the analytical solution and the contact-based result in case of a hemispherical target shape indicates the potential use of our method in form finding problems. Here, we limited ourselves within the convex surfaces with analytical solutions. The shape construction for arbitrary surfaces may need to introduce the frictional contact between the stretchable gridshells and target surfaces. We hope that our results and methodology will instigate future work on buckling induced mechanically guided assembly in physical systems (e.g., pop-up actuation of a planar grid to a target shape) from macro scale (e.g., domes in architecture) to micro-scale (e.g., controlled buckling of slender rods for stretchable electronics).

## Acknowledgment

We acknowledge support from the National Science Foundation (Award # IIS-1925360) and the Henry Samueli School of Engineering and Applied Science, University of California, Los Angeles.

## Conflict of Interest

There are no conflicts of interest.

## Data Availability Statement

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.

## References

1.
Ghiyasinasab
,
M.
,
Lehoux
,
N.
, and
Ménard
,
S.
,
2017
, “
Production Phases and Market for Timber Gridshell Structures: A State-of-the-Art Review
,”
BioResources
,
12
(
4
), pp.
9538
9555
.
2.
Quagliaroli
,
M.
, and
Malerba
,
P.
,
2013
, “
Flexible Bridge Decks Suspended by Cable Nets. a Constrained Form Finding Approach
,”
Int. J. Solids. Struct.
,
50
(
14–15
), pp.
2340
2352
. 10.1016/j.ijsolstr.2013.03.009
3.
Tayeb
,
F.
,
Caron
,
J.-F.
,
Baverel
,
O.
, and
Du Peloux
,
L.
,
2013
, “
Stability and Robustness of a 300 M2 Composite Gridshell Structure
,”
Constr. Build. Mater.
,
49
, pp.
926
938
. 10.1016/j.conbuildmat.2013.04.036
4.
Lefevre
,
B.
,
Douthe
,
C.
, and
Baverel
,
O.
,
2015
, “
Buckling of Elastic Gridshells
,”
J. Int. Assoc. Shell Spatial Struct.
,
56
(
3
), pp.
153
171
.
5.
Li
,
Y.
,
Song
,
J.
,
Fang
,
B.
, and
Zhang
,
J.
,
2011
, “
Surface Effects on the Postbuckling of Nanowires
,”
J. Phys. D: Appl. Phys.
,
44
(
42
), p.
425304
. 10.1088/0022-3727/44/42/425304
6.
Xu
,
S.
,
Yan
,
Z.
,
Jang
,
K.-I.
,
Huang
,
W.
,
Fu
,
H.
,
Kim
,
J.
,
Wei
,
Z.
,
Flavin
,
M.
,
McCracken
,
J.
,
Wang
,
R.
, and
,
A.
,
2015
, “
Assembly of Micro/Nanomaterials Into Complex, Three-Dimensional Architectures by Compressive Buckling
,”
Science
,
347
(
6218
), pp.
154
159
. 10.1126/science.1260960
7.
Zhao
,
H.
,
Li
,
K.
,
Han
,
M.
,
Zhu
,
F.
,
,
A.
,
Guo
,
P.
,
Xie
,
Z.
,
Park
,
Y.
,
Chen
,
L.
,
Wang
,
X.
, and
Luan
,
H.
,
2019
, “
Buckling and Twisting of Advanced Materials Into Morphable 3D Mesostructures
,”
Proc. Natl. Acad. Sci. U. S. A.
,
116
(
27
), p.
201901193
.
8.
Fan
,
J. A.
,
Yeo
,
W. -H.
,
Su
,
Y.
,
Hattori
,
Y.
,
Lee
,
W.
,
Jung
,
S. -Y.
,
Zhang
,
Y.
,
Liu
,
Z.
,
Cheng
,
H.
,
Falgout
,
L.
, and
Bajema
,
M
,
2014
, “
Fractal Design Concepts for Stretchable Electronics
,”
Nat. Commun.
,
5
, p.
3266
. 10.1038/ncomms4266
9.
Jang
,
K.-I.
,
Li
,
K.
,
Chung
,
H. U.
,
Xu
,
S.
,
Jung
,
H. N.
,
Yang
,
Y.
,
Kwak
,
J. W.
,
Jung
,
H. H.
,
Song
,
J.
,
Yang
,
C.
, and
Wang
,
A.
,
2017
, “
Self-Assembled Three Dimensional Network Designs for Soft Electronics
,”
Nat. Commun.
,
8
, p.
15894
. 10.1038/ncomms15894
10.
Jang
,
K.-I.
,
Chung
,
H. U.
,
Xu
,
S.
,
Lee
,
C. H.
,
Luan
,
H.
,
Jeong
,
J.
,
Cheng
,
H.
,
Kim
,
G.-T.
,
Han
,
S. Y.
,
Lee
,
J. W.
, and
Kim
,
J.
,
2015
, “
Soft Network Composite Materials With Deterministic and Bio-Inspired Designs
,”
Nat. Commun.
,
6
, p.
6566
. 10.1038/ncomms7566
11.
Xu
,
Z.
,
Fan
,
Z.
,
Fu
,
H.
,
Liu
,
Y.
,
Zi
,
Y.
,
Huang
,
Y.
, and
Zhang
,
Y.
,
2019
, “
Optimization-Based Approach for the Inverse Design of Ribbon-Shaped Three-Dimensional Structures Assembled Through Compressive Buckling
,”
Phys. Rev. Appl.
,
11
(
5
), p.
054053
. 10.1103/PhysRevApplied.11.054053
12.
Baek
,
C.
,
Sageman-Furnas
,
A. O.
,
Jawed
,
M. K.
, and
Reis
,
P. M.
,
2018
, “
Form Finding in Elastic Gridshells
,”
Proc. Natl. Acad. Sci. U. S. A.
,
115
(
1
), pp.
75
80
. 10.1073/pnas.1713841115
13.
Baek
,
C.
, and
Reis
,
P. M.
,
2019
, “
Rigidity of Hemispherical Elastic Gridshells Under Point Load Indentation
,”
J. Mech. Phys. Solids.
,
124
, pp.
411
426
. 10.1016/j.jmps.2018.11.002
14.
Panetta
,
J.
,
Konaković-Luković
,
M.
,
Isvoranu
,
F.
,
Bouleau
,
E.
, and
Pauly
,
M.
,
2019
, “
X-shells: A New Class of Deployable Beam Structures
,”
ACM Trans. Graph. (TOG)
,
38
(
4
), p.
83
. 10.1145/3306346.3323040
15.
Wood
,
R. D.
, and
Zienkiewicz
,
O.
,
1977
, “
Geometrically Nonlinear Finite Element Analysis of Beams, Frames, Arches and Axisymmetric Shells
,”
Comput. Struct.
,
7
(
6
), pp.
725
735
. 10.1016/0045-7949(77)90027-X
16.
Liu
,
Y.
,
Wang
,
X.
,
Xu
,
Y.
,
Xue
,
Z.
,
Zhang
,
Y.
,
Ning
,
X.
,
Cheng
,
X.
,
Xue
,
Y.
,
Lu
,
D.
,
Zhang
,
Q.
, and
Zhang
,
F.
,
2019
, “
Harnessing the Interface Mechanics of Hard Films and Soft Substrates for 3d Assembly by Controlled Buckling
,”
Proc. Natl. Acad. Sci. U. S. A.
,
116
(
31
), pp.
15368
15377
. 10.1073/pnas.1907732116
17.
Grinspun
,
E.
,
Desbrun
,
M.
,
Polthier
,
K.
,
Schröder
,
P.
, and
Stern
,
A.
,
2006
, “
Discrete Differential Geometry: An Applied Introduction
,”
ACM SIGGRAPH Course
,
7
, pp.
1
139
.
18.
Bergou
,
M.
,
Wardetzky
,
M.
,
Robinson
,
S.
,
Audoly
,
B.
, and
Grinspun
,
E.
,
2008
, “
Discrete Elastic Rods
,”
ACM Trans. Graph. (TOG)
,
27
(
3
), p.
63
. 10.1145/1360612.1360662
19.
Bergou
,
M.
,
Audoly
,
B.
,
Vouga
,
E.
,
Wardetzky
,
M.
, and
Grinspun
,
E.
,
2010
, “
,”
ACM Trans. Graph. (TOG)
,
29
(
4
), p.
116
. 10.1145/1778765.1778853
20.
Audoly
,
B.
, and
Pomeau
,
Y.
,
2010
,
Elasticity and Geometry: From Hair Curls to the Non-linear Response of Shells
,
Oxford University Press
,
Oxford
.
21.
Jawed
,
M. K.
,
Novelia
,
A.
, and
O’Reilly
,
O. M.
,
2018
,
A Primer on the Kinematics of Discrete Elastic Rods
,
Springer
,
New York City
.
22.
Li
,
X.
,
Huang
,
W.
, and
Jawed
,
M. K.
,
2020
, “
Discrete Elasto-Plastic Rods
,”
Extreme Mech. Lett.
,
38
, p.
100767
. 10.1016/j.eml.2020.100767
23.
Li
,
X.
,
Huang
,
W.
, and
Jawed
,
M. K.
,
2019
, “
A Discrete Differential Geometry-Based Approach to Numerical Simulation of Timoshenko Beam
,”
Extreme Mech. Lett.
,
35
, p.
100622
.
24.
Lestringant
,
C.
, and
Kochmann
,
D. M.
,
2020
, “
Modeling of Flexible Beam Networks and Morphing Structures by Geometrically Exact Discrete Beams
,”
ASME J. Appl. Mech.
,
87
(
8
), p.
081006
. 10.1115/1.4046895
25.
Shen
,
Z.
,
Huang
,
J.
,
Chen
,
W.
, and
Bao
,
H.
,
2015
, “Geometrically Exact Simulation of Inextensible Ribbon,”
Computer Graphics Forum
, Vol.
34
,
Helwig
Hauser
, ed.,
Wiley Online Library
,
New Jersey
, pp.
145
154
.
26.
Huang
,
W.
,
Wang
,
Y.
,
Li
,
X.
, and
Jawed
,
M. K.
,
2020
, “
Shear Induced Supercritical Pitchfork Bifurcation of Pre-Buckled Bands, From Narrow Strips to Wide Plates
,”
J. Mech. Phys. Solids.
,
145
, p.
104168
. 10.1016/j.jmps.2020.104168
27.
Baraff
,
D.
, and
Witkin
,
A.
,
1998
, “
Large Steps in Cloth Simulation
,”
Proceedings of the 25th Annual Conference on Computer Graphics and Interactive Techniques
,
Orlando, FL
, ACM, pp.
43
54
.
28.
Grinspun
,
E.
,
Hirani
,
A. N.
,
Desbrun
,
M.
, and
Schröder
,
P.
,
2003
, “
Discrete Shells
,”
Proceedings of the 2003 ACM SIGGRAPH/Eurographics Symposium on Computer Animation
,
San Diego, CA
, pp.
62
67
.
29.
Batty
,
C.
,
Uribe
,
A.
,
Audoly
,
B.
, and
Grinspun
,
E.
,
2012
, “
Discrete Viscous Sheets
,”
ACM Trans. Graph. (TOG)
,
31
(
4
), p.
113
. 10.1145/2185520.2185609
30.
Qin
,
L.
,
Huang
,
W.
,
Du
,
Y.
,
Zheng
,
L.
, and
Jawed
,
M. K.
,
2020
, “
Genetic Algorithm-Based Inverse Design of Elastic Gridshells
,”
Struct. Multidisci. Optim.
,
62
(
5
), pp.
2691
2707
. 10.1007/s00158-020-02639-8
31.
Spillmann
,
J.
, and
Teschner
,
M.
,
2009
, “
Cosserat Nets
,”
IEEE Trans. Visual. Comput. Graph.
,
15
(
2
), pp.
325
338
. 10.1109/TVCG.2008.102
32.
Pérez
,
J.
,
Thomaszewski
,
B.
,
Coros
,
S.
,
Bickel
,
B.
,
Canabal
,
J. A.
,
Sumner
,
R.
, and
,
M. A.
,
2015
, “
Design and Fabrication of Flexible Rod Meshes
,”
ACM Trans. Graph. (TOG)
,
34
(
4
), p.
138
.
33.
Fu
,
H.
,
Nan
,
K.
,
Froeter
,
P.
,
Huang
,
W.
,
Liu
,
Y.
,
Wang
,
Y.
,
Wang
,
J.
,
Yan
,
Z.
,
Luan
,
H.
,
Guo
,
X.
, and
Zhang
,
Y.
,
2017
, “
Mechanically-guided Deterministic Assembly of 3D Mesostructures Assisted by Residual Stresses
,”
Small
,
13
(
24
), p.
1700151
. 10.1002/smll.201700151
34.
Yan
,
Z.
,
Han
,
M.
,
Yang
,
Y.
,
Nan
,
K.
,
Luan
,
H.
,
Luo
,
Y.
,
Zhang
,
Y.
,
Huang
,
Y.
, and
Rogers
,
J. A.
,
2017
, “
Deterministic Assembly of 3D Mesostructures in Advanced Materials Via Compressive Buckling: A Short Review of Recent Progress
,”
Extreme Mech. Lett.
,
11
, pp.
96
104
. 10.1016/j.eml.2016.12.006
35.
Zhang
,
Y.
,
Zhang
,
F.
,
Yan
,
Z.
,
Ma
,
Q.
,
Li
,
X.
,
Huang
,
Y.
, and
Rogers
,
J. A.
,
2017
, “
Printing, Folding and Assembly Methods for Forming 3D Mesostructures in Advanced Materials
,”
Nat. Rev. Mater.
,
2
(
4
), p.
17019
. 10.1038/natrevmats.2017.19
36.
Rogers
,
J.
,
Huang
,
Y.
,
Schmidt
,
O. G.
, and
Gracias
,
D. H.
,
2016
, “
Origami MEMS and NEMS
,”
Mrs Bull.
,
41
(
2
), pp.
123
129
. 10.1557/mrs.2016.2
37.
Kim
,
D.-H.
,
Lu
,
N.
,
Ghaffari
,
R.
,
Kim
,
Y.-S.
,
Lee
,
S. P.
,
Xu
,
L.
,
Wu
,
J.
,
Kim
,
R.-H.
,
Song
,
J.
,
Liu
,
Z.
, and
Viventi
,
J.
,
2011
, “
Materials for Multifunctional Balloon Catheters With Capabilities in Cardiac Electrophysiological Mapping and Ablation Therapy
,”
Nat. Mater.
,
10
(
4
), p.
316
. 10.1038/nmat2971
38.
Dai
,
X.
,
Zhou
,
W.
,
Gao
,
T.
,
Liu
,
J.
, and
Lieber
,
C. M.
,
2016
, “
Three-Dimensional Mapping and Regulation of Action Potential Propagation in Nanoelectronics-Innervated Tissues
,”
Nat. Nanotechnol.
,
11
(
9
), p.
776
. 10.1038/nnano.2016.96
39.
Yu
,
X.
,
Wang
,
H.
,
Ning
,
X.
,
Sun
,
R.
,
,
H.
,
Salomao
,
M.
,
Silva
,
A. C.
,
Yu
,
Y.
,
Tian
,
L.
,
Koh
,
A.
, and
Lee
,
C. M.
,
2018
, “
Needle-shaped Ultrathin Piezoelectric Microsystem for Guided Tissue Targeting Via Mechanical Sensing
,”
Nat. Biomed. Engin.
,
2
(
3
), p.
165
. 10.1038/s41551-018-0201-6
40.
Feiner
,
R.
,
Engel
,
L.
,
Fleischer
,
S.
,
Malki
,
M.
,
Gal
,
I.
,
Shapira
,
A.
,
Shacham-Diamand
,
Y.
, and
Dvir
,
T.
,
2016
, “
Engineered Hybrid Cardiac Patches With Multifunctional Electronics for Online Monitoring and Regulation of Tissue Function
,”
Nat. Mater.
,
15
(
6
), p.
679
. 10.1038/nmat4590
41.
Liu
,
Y.
,
Pharr
,
M.
, and
Salvatore
,
G. A.
,
2017
, “
Lab-on-Skin: a Review of Flexible and Stretchable Electronics for Wearable Health Monitoring
,”
ACS. Nano.
,
11
(
10
), pp.
9614
9635
. 10.1021/acsnano.7b04898
42.
Yao
,
S.
,
Swetha
,
P.
, and
Zhu
,
Y.
,
2018
, “
Nanomaterial-Enabled Wearable Sensors for Healthcare
,”
,
7
(
1
), p.
1700889
.
43.
Song
,
Z.
,
Ma
,
T.
,
Tang
,
R.
,
Cheng
,
Q.
,
Wang
,
X.
,
Krishnaraju
,
D.
,
Panat
,
R.
,
Chan
,
C. K.
,
Yu
,
H.
, and
Jiang
,
H.
,
2014
, “
Origami Lithium-ion Batteries
,”
Nat. Commun.
,
5
, p.
3140
. 10.1038/ncomms4140
44.
Tang
,
R.
,
Huang
,
H.
,
Tu
,
H.
,
Liang
,
H.
,
Liang
,
M.
,
Song
,
Z.
,
Xu
,
Y.
,
Jiang
,
H.
, and
Yu
,
H.
,
2014
, “
Origami-enabled Deformable Silicon Solar Cells
,”
Appl. Phys. Lett.
,
104
(
8
), p.
083501
. 10.1063/1.4866145
45.
Xiao
,
X.
,
Zhou
,
W.
,
Kim
,
Y.
,
Ryu
,
I.
,
Gu
,
M.
,
Wang
,
C.
,
Liu
,
G.
,
Liu
,
Z.
, and
Gao
,
H.
,
2015
, “
Regulated Breathing Effect of Silicon Negative Electrode for Dramatically Enhanced Performance of Li-ion Battery
,”
,
25
(
9
), pp.
1426
1433
46.
Alameh
,
Z.
,
Deng
,
Q.
, and
Sharma
,
P.
,
2018
, “
Emergent Magnetoelectricity in Soft Materials, Instability, and Wireless Energy Harvesting
,”
Soft. Matter.
,
14
(
28
), pp.
5856
5868
. 10.1039/C8SM00587G
47.
Lee
,
J.-H.
,
Koh
,
C. Y.
,
Singer
,
J. P.
,
Jeon
,
S.-J.
,
Maldovan
,
M.
,
Stein
,
O.
, and
Thomas
,
E. L.
,
2014
, “
25th Anniversary Article: Ordered Polymer Structures for the Engineering of Photons and Phonons
,”
,
26
(
4
), pp.
532
569
48.
Valentine
,
J.
,
Zhang
,
S.
,
Zentgraf
,
T.
,
Ulin-Avila
,
E.
,
Genov
,
D. A.
,
Bartal
,
G.
, and
Zhang
,
X.
,
2008
, “
Three-Dimensional Optical Metamaterial With a Negative Refractive Index
,”
Nature
,
455
(
7211
), p.
376
. 10.1038/nature07247
49.
Filipov
,
E. T.
,
Tachi
,
T.
, and
Paulino
,
G. H.
,
2015
, “
Origami Tubes Assembled Into Stiff, Yet Reconfigurable Structures and Metamaterials
,”
Proc. Natl. Acad. Sci. U. S. A.
,
112
(
40
), pp.
12321
12326
. 10.1073/pnas.1509465112
50.
Reis
,
P. M.
,
Jiménez
,
F. L.
, and
Marthelot
,
J.
,
2015
, “
Transforming Architectures Inspired by Origami
,”
Proc. Natl. Acad. Sci. U. S. A.
,
112
(
40
), pp.
12234
12235
. 10.1073/pnas.1516974112
51.
Garg
,
A.
,
Sageman-Furnas
,
A. O.
,
Deng
,
B.
,
Yue
,
Y.
,
Grinspun
,
E.
,
Pauly
,
M.
, and
Wardetzky
,
M.
,
2014
, “
Wire Mesh Design
,”
ACM Trans. Graphics
,
33
(
4
), pp.
66
. 10.1145/2601097.2601106
52.
Bouhaya
,
L.
,
Baverel
,
O.
, and
Caron
,
J.-F.
,
2010
, “
Mapping Two-way Continuous Elastic Grid on An Imposed Surface: Application to Grid Shells
,”
Symposium of the International Association for Shell and Spatial Structures. (50th. 2009. Valencia) Evolution and Trends in Design, Analysis and Construction of Shell and Spatial Structures: Proceedings, Editorial Universitat Politècnica de València
,
Valencia
.
53.
Huang
,
W.
,
Huang
,
X.
,
Majidi
,
C.
, and
Jawed
,
M. K.
,
2020
, “
Dynamic Simulation of Articulated Soft Robots
,”
Nat. Commun.
,
11
(
1
), pp.
1
9
.
54.
Huang
,
W.
, and
Jawed
,
M. K.
,
2019
, “
Newmark-Beta Method in Discrete Elastic Rods Algorithm to Avoid Energy Dissipation
,”
ASME J. Appl. Mech.
,
86
(
8
), p.
084501
. 10.1115/1.4043793
55.
Jawed
,
M. K.
,
Da
,
F.
,
Joo
,
J.
,
Grinspun
,
E.
, and
Reis
,
P. M.
,
2014
, “
Coiling of Elastic Rods on Rigid Substrates
,”
Proc. Natl. Acad. Sci. U. S. A.
,
111
(
41
), pp.
14663
14668
. 10.1073/pnas.1409118111
56.
Jawed
,
M. K.
,
Khouri
,
N.
,
Da
,
F.
,
Grinspun
,
E.
, and
Reis
,
P. M.
,
2015
, “
Propulsion and Instability of a Flexible Helical Rod Rotating in a Viscous Fluid
,”
Phys. Rev. Lett.
,
115
(
16
), p.
168101
. 10.1103/PhysRevLett.115.168101
57.
Jawed
,
M. K.
,
,
N.
,
Parks
,
D.
, and
Reis
,
P.
,
2018
, “
Patterns of Carbon Nanotubes by Flow-Directed Deposition on Substrates With Architectured Topographies
,”
Nano Lett.
,
2
(
3
), p.
165
.
58.
Liu
,
Z.
,
Zhu
,
D.
,
Rodrigues
,
S. P.
,
Lee
,
K.-T.
, and
Cai
,
W.
,
2018
, “
Generative Model for the Inverse Design of Metasurfaces
,”
Nano. Lett.
,
18
(
10
), pp.
6570
6576
. 10.1021/acs.nanolett.8b03171