From the elastica compass to the elastica catapult: an essay on the mechanics of soft robot arm

An elastic rod is clamped at one end and has a dead load attached to the other (free) end. The rod is then slowly rotated using the clamp. When the load is smaller than the buckling value, the rod describes a continuous set of quasi-static forms and its end traces a (smooth, convex and simple) closed curve, which would be a circle if the rod were rigid. The closed curve is analytically determined through the integration of the Euler’s elastica, so that for sufficiently small loads the mechanical system behaves as an ‘elastica compass’. For loads higher than that of buckling, the elastica reaches a configuration from which a snap-back instability occurs, realizing a sort of ‘elastica catapult’. The whole quasi-static evolution leading to the critical configuration for snapping is calculated through the elastica and the subsequent dynamic motion simulated using two numerical procedures, one ad hoc developed and another based on a finite-element scheme. The theoretical results are then validated on a specially designed and built apparatus. An obvious application of the present model would be in the development of soft robotic limbs, but the results are also of interest for the optimization analysis in pole vaulting.


Introduction
The design of innovative devices for advanced applications is being driven by the need for compliant mechanisms, which are usually inspired by nature [1,2] and is part of a transition from traditional robotics to soft robotics [3][4][5]. Compliant   An elastic rod of length l, with bending stiffness B and linear mass density γ , has attached a lumped mass m at one end and is constrained at the other end by a slowly rotating clamp, inclined at an angle α (increasing function of time t) with respect to the direction of the gravity, so that a force P = mg is applied to an end of the rod. The rotation of the rod's axis with respect to the undeformed (straight) configuration is measured through the angle θ (s, t), with s being the curvilinear coordinate, s ∈ [0, l]. Thex −ŷ and x − y reference systems are reported, both centred at the clamp point, the former attached to the rotating clamp, whereas the latter fixed as the loading direction. The polar coordinates r − ϕ defining the rod's end position are also reported. (Online version in colour.) development and use of nonlinear mechanical models such as the Kirchhoff rod [6] and Euler's elastica, which allow the description of large deflections in elastic bars and the modelling of snake locomotion [7,8], as well as object manipulation [9][10][11], useful in robotic assistance during surgery [12,13] and for physical rehabilitation [14].
In this article, a basic mechanical model for a soft robot arm is addressed through new theoretical, numerical and experimental developments. In particular, the deformable mechanical system sketched in figure 1 is considered, in which an elastic rod is clamped at one end and subject to a dead load at the other. The load is provided by the weight of a mass predominantly higher than that of the rod. The clamp rotates slowly, so that starting from a configuration in which the rod is subject to purely tensile axial load, the system quasi-statically evolves in a number of elastic forms at varying clamp angle. When the load is inferior to that corresponding to buckling of the straight and uniformly compressed configuration, a whole quasi-static 360 • rotation of the clamp is possible and the edge of the rod describes a (smooth, convex and simple) closed curve, which is theoretically solved using Euler's elastica. In the case of a rigid rod, this curve would be a circle, so that the mechanical system behaves as an 'elastica compass', thus tracing the curve described by the elastica. However, when the load is higher than that which would lead to buckling for the straight rod, an unstable configuration is quasi-statically reached, at which the rod suffers a snapback instability and dynamically approaches another configuration, so that the system behaves as an 'elastica catapult'. 1 The description of the quasi-static path of the system and the determination of the unstable configuration is solved in an analytical form by means of elliptic functions through an extension of results pioneered by Wang [15], who employed a numerical integration procedure and provided asymptotic estimates valid in some particular cases (very stiff and very compliant rods, and nearly vertical equilibrium configurations). In addition to the quasi-static solution, the dynamics of the snap instability is addressed numerically. The set-up of a numerical technique is a complex problem, which was analysed from several points of view, but not still completely solved [16][17][18][19][20][21]. To this purpose, two approaches are presented, one is a standard use of a finite-element software (Abaqus), whereas the other is developed as a perfection of a technique introduced for pneumatic soft robot arms [22]. The latter approach, in which the elastic rod is reduced to a nonlinear spring governed by the elastica, is elegant, but the kinematics is limited to the first deformation mode and an axial deformation and viscous damping have to be added to prevent spurious numerical instabilities, issues which may be circumvented through the finite-element approach.
Finally, the experimental validation of the elastic system was performed using a mechanical set-up specifically designed and realized at the 'Instabilities Lab' of the University of Trento (http://www.ing.unitn.it/dims/ssmg/). Experimental results (also available as a movie in the electronic supplementing material) fully validate the theoretical modelling, thus confirming that the elastica allows for solutions useful in the kinematics of a soft robot arm. The performance of the robot arm is also assessed in terms of (i) the maximum and minimum distances that can be reached without encountering loss of stability of the configuration and (ii) the maximum energy release that can be achieved when the system behaves as a catapult. These results open the way to a rational design of deformable robot arms and, as a side development, may find also application in the analysis of the pole vault dynamics and the optimization of athletes' performance [23,24].

Formulation
An inextensible planar rod with bending stiffness B, length l, linear mass density γ , and straight in its undeformed configuration, has a lumped mass m attached at one end, whereas the other end is constrained by a clamp having inclination α with respect to gravity direction (figure 1). Denoting with g the gravitational acceleration, the rod is then loaded by the weight P owing to the lumped mass m, so that P = mg, and the rod distributed weight γ g (the latter neglected in the quasi-static analysis). The clamp angle α smoothly and slowly increases in time t, so that a quasi-statically rotating clamp is realized. For simplicity of presentation, the dependence on time t is omitted in the notation in the following of this section. The rotation of the rod's axis with respect to the undeformed (straight) configuration is denoted by θ(s), function of the curvilinear coordinate s ∈ [0, l], with s = 0 singling out the position of the clamp (where θ(0) = 0) and s = l of the loaded rod's end. With respect to the undeformed straight configuration, 'frozen' at the inclination angle α, the coordinatesx(s) andŷ(s) measure the position of the rod's axis in the rotating system, and, owing to the inextensibility condition, are connected to the rotation field through the following differential relations [25]x where a prime denotes the derivative along the curvilinear coordinate s. The position can be described through the coordinates x(s) and y(s) (the former orthogonal and the latter parallel, but with opposite direction, to the gravity), which are connected with the positionsx(s) andŷ(s) through the following relationships During the rotation, the position of the clamp (s = 0) is considered fixed and taken as the origin of the reference systems, Neglecting the rotational inertia of the lumped mass and of the rod, the Lagrangian functional for the system is given by (2.5) where N x (s) and N y (s) are the internal forces aligned in parallel with the x-and y-directions (which play the role of Lagrangian multipliers), T is the kinetic energy of the system (with a dot denoting the time derivative) and V is the sum of the elastic energy stored within the rod and the negative of the work done by the dead load P and that by the rod distributed weight γ g, The governing equations can be derived through application of the least action principle on the functional with t 0 and t * being arbitrary initial and final instants of the analysed time interval. Reference is made to a perturbation in the fields {x(s), y(s), θ (s)} by means of the small parameter and of the field variations {x var (s), y var (s), θ var (s)} satisfying the time and space conditions x var (s) = y var (s) = θ var (s) = 0, for t = t 0 and t = t * and Therefore, the minimum of the functional A, equation (2.8), is obtained by imposing the following condition and the boundary conditions and

Quasi-static response
The weight of the rod γ gl is considered here negligible when compared with that of the lumped mass attached at the end of the rod P. Moreover, when quasi-static conditions prevail, the acceleration of the rod's axis can be neglected,ẍ(s) =ÿ(s) = 0, so that from the boundary conditions (2.12) 2 and (2.12) 3 , and the differential equations (2.11) 2 and (2.11) 3 , the rod's internal forces result to be constant along the rod N x (s) = 0 and N y (s) = −P, (3.1) so that the governing equation for the rotation field (2.11) 1 reduces to the elastica [25,26] Bθ (s) − P sin[θ (s) + α] = 0, (3.2) to be complemented with the boundary conditions, equation (2.12) 1 and θ(0) = 0. Introducing the symbol λ 2 = P/B and the auxiliary rotation ψ(s) = θ(s) + α − π (measuring the inclination of the rod tangent with respect to the y axis, see figure 1), the elastica (3.2) and the boundary conditions can be rewritten as ψ (s) + λ 2 sin ψ(s) = 0, ψ(0) = α − π and ψ (l) = 0, (3.3) so that, through the following change of variables the integration of the differential problem (3.3) leads to the following expression providing the relation between the applied load P and the rotation θ(l) = θ l at the loaded end of the rod In equation (3.5), n = 1, 2, 3 . . . is an integer representing the nth mode solution, whereas K(k) and K(σ 0 , k) are, respectively, the complete and the incomplete elliptic integral of the first kind of modulus k In the following, the analysis is mainly addressed to the deformation mode n = 1, because only the equilibrium configurations related to this mode are stable [25,27]. However (when existing), the unstable equilibrium configurations related to deformation mode with n = 1 will also be considered in the case P/P cr > 4 to provide the whole equilibrium paths of the system. Once the rotation θ l is computed from the nonlinear equation (3.5) for a given load P and a clamp angle α, the rotation field θ (s) is obtained from integration of equation (3.3) as the solution of the equation where Sn(·, k) represents the sinus of the Jacobian amplitude am(·, k) of modulus k. With reference to the buckling load of the purely compressed clamped rod (α = π ), namely P cr = π 2 B/(4 l 2 ), the solution of the nonlinear equation (3.5) displays the two following different behaviours: -when P ≤ P cr , a unique value of the rotation at the loaded end θ l corresponds to a unique value of the clamp inclination α; -when P > P cr , more than one solution for the rotation at the loaded end θ l may exist when the clamp inclination α falls within the interval (2π − α s , α s ), with α s ∈ [π , 2π ].
These two behaviours are highlighted in figure 2, where the rotation (with respect to the vertical direction) θ l + α of the loaded rod's end is reported as a function of the clamp inclination α, solution of the nonlinear equation (3.5). Six values of the ratio P/P cr have been considered, namely {0.5, 1, 3} and {5, 8, 9.161} in the left and right parts of figure 2, respectively. Uniqueness of the end rotation θ l as a function of α is displayed only when the load P does not exceed the critical load P cr (case P = 0.5P cr ), whereas more than one equilibrium configuration may be found when P > P cr for some set of values for the clamp angle α. For instance, in the case P = 3P cr ( left), three equilibrium configurations (related to n = 1) are displayed for the clamp angle α within the interval (2π − α s , α s ). Note that when the equilibrium configuration is not unique, only two of the equilibrium configurations are stable (represented as continuous lines in the figure), whereas the others are unstable (represented as discontinuous lines). The limit case P = P cr (figure 2, left) is also reported, for which a vertical tangent is displayed at α = π , which defines the transition between the two behaviours.
For completeness, it is observed that when P ≤ q 2 P cr (q ∈ N), the nonlinear equation (3.5) may admit solutions for nth mode, with n < q. For instance, when P ∈ [9, 16]P cr , a total of five equilibrium configurations may exist for the same value of α (figure 2, right, case P = 9.161P cr ).
The uniqueness/non-uniqueness of the quasi-static solution defines two qualitatively different mechanical responses for the analysed elastic system. Indeed, considering a monotonic increase of the clamp angle α from 0 to 2π : -when P ≤ P cr , the rotation θ l changes continuously, so that the end of the rod describes a (smooth, convex and simple) closed continuous curve. In this condition, the system behaves as an elastica compass; -when P > P cr , the rotation θ l reaches a critical value (corresponding to the snap clamp inclination α s ∈ [π , 2π ]), for which a further increase in the clamp angle necessarily yields a jump in the rotation θ l . Such a jump involves a release of elastic energy and a dynamic snap to another, non-adjacent configuration. In this condition, the system behaves as an elastica catapult.
The clamp angle α s for which the snap-back instability occurs has been numerically evaluated and is reported as a function of the load ratio P/P cr in figure 3. It can be noted that α s is always greater than π (limit value attained when P coincides with the buckling load, P = P cr ) and is an increasing function of the applied load, so that as the applied load increases the snap-back occurs 'later'. Results from the dynamic analyses (provided in §7) and from the experimental tests on the developed physical prototype (described in §8) are also reported in the figure. The snap condition related to the load value P si , defining the limit condition of self-intersection (see §4), is highlighted. The agreement of the analytical predictions with the results obtained from both the

The elastica compass and the elastica catapult
A further insight into the mechanics of the elastica compass and catapult requires the development of the rod's kinematics. From the rotation field (3.7), the integration of equation (2.3) provides the position of the rod's axis, that in the x − y reference system becomes (see [28,29]) and Cn denotes the Jacobi cosine amplitude function Cn(·, k) = cos(am(·, k)), whereas E is the incomplete elliptic integral of the second kind, Considering a fixed weight P, the quasi-static evolution of the deformed configuration can be represented by varying the clamp angle α using the kinematical description (4.1) (figure 4). It is found that, when P = P si ≈ 9.161P cr , the deformed configuration displays a self-contact point with the clamp at the verge of the snap-back, α s ≈ 1.784π . Therefore, the load P si defines the lowest value of the load needed to achieve self-intersection of the elastic rod during rotation of the clamp. Depending on the out-of-plane geometry of the rod, two behaviours can be attained in the case P > P si : (i) if the geometry permits self-intersection, the present solution holds and selfintersecting elastica are displayed, whereas (ii) if the geometry does not permit self-intersection, a contact point is formed within the configuration of the rod at increasing the clamp rotation. The response of the system in the case P > P si will be the subject of future investigation.
The trajectory travelled by the loaded end (playing the role of the pencil lead of the 'elastica compass') can be traced by evaluating the coordinates (4.1) at the loaded end, s = l, at varying clamp inclination α ∈ [0, 2π ]. The quasi-static trajectories are reported in figure 5 for different  Figure 4. A sequence of quasi-static configurations of the elastica catapult (P/P cr > 1) at increasing the clamp angle α for different values of P/P cr ending with snap in the equilibrium configuration. The self-intersection of the rod occurring for P ≥ P si is displayed. In particular, the limit case of self-intersection is shown in the central sequence (P = P si ≈ 9.161P cr ), whereas self-intersection is shown in the lower sequence (P = 12P cr ). (Online version in colour.)  values of the ratio P/P cr . It can be observed that the trajectories have the shape of (smooth, convex and simple) closed curves in the case P < P cr . Furthermore, unstable positions for the rod's end are reported as a discontinuous line in the case P > P cr (as dashed and dotted line for the first and second mode configurations, respectively), so that the snap-back instability is initiated at the point where the continuous line ends. The position of the loaded end can also be described in a polar reference system through the radius r = x(l) 2 + y(l) 2 and the angle ϕ = 3π/2 − arctan[y(l)/x(l)], which result The polar coordinates r (made dimensionless through division of the rod's length l) and ϕ which describe the rod's end trajectory are reported in figure 6 at varying the clamp angle α for different values of the ratio P/P cr .  From figures 5 and 6, it can be observed that -the behaviour of the usual, namely undeformable, compass is recovered in the limit of vanishing P/P cr , for which the elastic rod behaves as a rigid bar, r(α) = l and ϕ(α) = α; -owing to the inextensibility assumption, the loops drawn by the elastica compass always lie inside the circle of radius l, therefore, r(α) ≤ l; -for dead loads smaller than the buckling load, P < P cr , the loops drawn by the elastica compass are nearly circular despite the large difference between ϕ and α. This is due to the fact that the maximum percentage decrease in the radius length is about 6%; -for dead loads larger than the buckling load, P > P cr , the polar coordinate ϕ is limited by the upper bound ϕ max (P/P cr ) = max α ϕ (α, P/P cr ), described by the dashed curve reported in figure 7a. Defining ϕ s as the polar angle at the verge of the snap-back instability, namely ϕ s (P/P cr ) = ϕ(α s (P/P cr ), P/P cr ) (reported as continuous curve), it is observed that ϕ s (P/P cr ) ≤ ϕ max (P/P cr ), where the equality holds for P/P cr ≤ 8.3.
In figure 7b, three soft robot arms with the same bending stiffness B, subject to the same weight The reaction moment M(0) is reported in figure 8, showing a change in sign at the snap-back, as a change in sign for the rod's curvature occurs. This feature has been exploited to detect the snap inclination α s from the results obtained with the numerical and experimental investigations explained below.

Robot's arm performance and design
The performances of the soft robot arm are investigated in terms of extremum distances that the loaded end can attain, which represent fundamental quantities in the design of soft robot arm, to achieve targeted positions for the hanged weight.
The maximum horizontal distance d x , the maximum height d y , and the minimum radial distance r min reached by the hanged weight during the clamp rotation and before the possible snap-back instability, are introduced as where α ∈ [0, π ] for the elastica compass (P < P cr ) and α ∈ [0, α s ] for the elastica catapult (P > P cr ). Considering fixed both the length l and the stiffness B and exploiting the kinematical description (4.1), the distances d x , d y , r min have been evaluated at varying load P and are reported in figure 9a, together with the respective clamp rotations α x , α y and α r for which these distances are attained, figure 9 (right). The distance d h (also called 'longest horizontal reach' by [30]) is defined as the horizontal distance attained by the weight when its vertical coordinate vanishes (namely the weight and the clamp are at the same height): and has been evaluated and reported in figure 9, together with the respective angle α h . From figure 9, the following conclusions can be drawn at varying the load P and considering constant both the length l and the stiffness B.
-The maximum horizontal distance d x attains its maximum value in the rigid limit (P/P cr = 0), for which d x = l, and is a decreasing function of the load P. The distance d x is always attained for a clamp angle α x ∈ [π/2, π ]. -In the case of the elastica compass, the maximum height corresponds to the robot arm length, d y = l, independently of the load P ≤ P cr and is attained for α = π . In the case of the elastica catapult, the maximum height d y decreases at increasing load P and is attained at the verge of the snap-back instability, α y = α s .  Figure 9. Maximum horizontal distance d x , maximum height d y , and minimum radial distance r min (a) and corresponding angles α x , α y and α r (b) plotted as functions of the load ratio P/P cr . The distance d h , equation (5.2), called 'longest horizontal reach' [30], the related angle α h , and the snap angle curve α s are also reported. (Online version in colour.) always remains below the clamp before the occurrence of the snap-instability, d y < 0. For such a load range, the longest reach d h corresponds only to unstable configurations (and therefore cannot be attained).
-The minimum radial distance attains its minimum value r min = 0.104 l for the load P/P cr = 7.985. Moreover, when P/P cr > 5.665, the minimum distance is achieved for a clamp rotation at the verge of instability (α r = α s ). -The distance d h is always smaller than or equal to the horizontal distance d x , namely d h ≤ d x , where the equality holds only in the case of vanishing load, P = 0. Moreover, the distance d h is defined only when P/P cr < 7.464, because, otherwise, the configurations are unstable at null vertical coordinate, y(l) = 0.
It is worth remarking that the curves reported in figure 9 are plotted in a dimensionless plane, in which both axes are affected by a change in the rod length l. Therefore, playing with this length would permit maximizing of the physical distances effectively reached by the loaded end when the hanged weight P is kept constant (as well as the bending stiffness B). Such maximum distances can be evaluated by seeking the load ratios P/P cr for which the following function attains a maximum P P cr d j (P/P cr ) The numerical implementation of this procedure leads to the following load ratios for which the soft arm displays the maximum possible distances P P cr d max  The maximum distances d y and d h and related clamp inclination α y and α h are reported as spots in figure 9. The optimization of the soft robot arm in terms of kinematics is completed by considering the maximum radial distance r max (ϕ), defined as the maximization of the radial distance r for a given polar angle ϕ (keeping fixed the applied load P and the bending stiffness B). Following equation (5.3), the load ratio P/P cr for which the maximum radial distance r max is achieved can be found by seeking the maximum of the function P P cr r(P/P cr , ϕ) l . (5.9) Restricting the attention to the polar coordinate ϕ ranging between π/4 and π , the optimal load ratio P/P cr has been reported in figure 10a as the result of the maximization of equation (5.3). Under this loading condition, the maximum radial distance r max (normalized through division by the length l) and the clamp inclination α (for which such a distance is attained) are evaluated at varying the angle ϕ in figure 10 (respectively, b and c). Furthermore, the arm length l and the maximum radial distance r max (made dimensionless through division by the constant parameters P and B) are provided in figure 10d. It can be observed that -all the curves display a jump in their first derivative for the value ϕ = 0.794π . For polar coordinates ϕ higher than this value, the load ratio P/P cr maximizing the radial distance corresponds to an angle equal to the maximum possible for that loading condition, ϕ = ϕ max (figure 10a); -the maximum radial distance is never observed for loads P smaller than the critical one P cr ; -the results reported in figure 10d represent an extension of the recent result by Wang [30] and Batista [28], referred to the case of an hanged load located along the horizontal axis x (at ϕ = π/2) and expressed by the distance d h , see equation (5.7) 2 . On the other hand, in the limit of ϕ = π , the value of the radial distance r max approaches the value of the length l, see equation (5.7) 1 . -The curve reported in figure 10d defines the values of the optimal lengths l for which the maximum radial distance is attained. In other words, lengthening does not always realize an increase in the achieved distances so that, using words similar to Wang [30], shortening of rods may provide longer distances, a key concept in the design of soft robot arms.
An applicative example of kinematic performance design towards the achievement of the maximum radial distance r max , for given load P and stiffness B, is reported in figure 7b for three soft robot arms differing only in their length. The deformed configurations, subject to the condition of loaded end lying along the polar coordinate ϕ = π/4, are shown for the three systems. It is observed that the maximum distance is attained when the soft arm has the optimal length l = 4.005 √ B/P, which value is provided by the curve reported in figure 10d for ϕ = π/4. Finally, with reference to the strength, the design of the arm cross section is ruled by the maximum bending moment M max experienced along the elastic arm, s ∈ [0, l], at varying the clamp inclination

M(s).
For the considered load range, P/P cr ∈ [0, 16], it is numerically found that the maximum bending moment always occurs at the clamp (s = 0) for the clamp inclination α = α x , so that (5.10)

Energy release
In the case of the elastic catapult, when the clamp inclination approaches the value α s , the snapback instability induces a dynamic motion in the elastic system. Such a motion is due to the energy release of the system, provided in terms of the release of both elastic energy stored in the rod and the potential energy of the hanged load. With reference to the two stable equilibrium configurations associated to the clamp rotation α ∈ [2π − α s , α s ], the energy difference E can be computed as  Because of the symmetry of the paths in figure 2, the following relation for the energy difference holds so that E(P, α = π ) = 0, (6.4) and therefore the clamp angle α = π represents the Maxwell line for the system [25,31,32]. The energy difference E(P, α) is reported in figure 11a at varying the load ratio P/P cr and the clamp angle α ∈ [π , α s ]. It is observed that the maximum energy difference, max α E(P, α) is attained at the verge of the snap instability, namely, at the clamp angle α = α s . Therefore, the energy release E R occurring at the snap-back, evaluated from the quasi-static solution, corresponds to the maximum energy difference possible for the elastic system, E R (P) = E(P, α s ) (6.5) which is reported as a function of the applied load P in figure 11b.

Snap-back dynamics
The dynamics of the rotating cantilever rod, occurring as a consequence of the instability when the clamp is quasi-statically rotating, is here investigated through two models: (i) a simplified model (called 'massless rod model') where the rod's inertia is neglected (because considered small when compared to the inertia of the point mass m) and regularizations are introduced through viscous damping owing to air drag acting on the mass and axial compliance of the clamp constraint, and (ii) a finite-element model (performed in Abaqus v. 14.3), in which the inertia of the rod is fully accounted for and Rayleigh internal damping is present. A slightly different version of model (i), which will be presented below, was proposed in the modelling of a tube-type manipulator arm by Snyder & Wilson [22]. It is shown that both models provide useful insights into the dynamics of the mechanical system after snap-back instability. However, while model (ii) is capable of describing the complete dynamic evolution of the system up to a wide range of loading parameters, the model (i) can predict the full dynamics only until deformation modes higher than the first do not enter the solution (in other words, the massless rod model works correctly for P/P cr smaller than ≈2).  rod's free end. Therefore, the rotation field of the rod's axis is governed by where the load R(t) is the force resultant acting at the free end and β(t) is the inclination angle of the force resultant with respect to the (straight) undeformed configuration ( figure 12), In equation (7.2), R x (t) and R y (t) are the resultant components of R(t) along the x and y axes which are given as the sum of the respective translational inertia, the (linear) viscous damping force, and the dead load where c is the damping coefficient owing to air drag. To model the axial compliance of the system (neglected in the quasi-static calculations), an axial spring of stiffness k c is introduced in the clamp, so that the kinematical constraint (2.3) is now replaced by Introducing the reference time T = ml 3 /B and the following dimensionless quantities 3) with the differential equations (7.4) leads to the following nondimensionalized differential system, governing the snap-back dynamics and Note that the elastica, equation (7.1), is the same as for the quasi-static problem (3.2), but with the two functions of time R(t) and β(t) replacing now P and α, respectively. Therefore, the integration in space for the dynamic problem, equation (7.1), can be computed similarly to the quasi-static case ( §3) restricted to the first deformation mode, n = 1. This integration yields the following relation between the normalized resultant μ(τ ), the rotation of the rod at the free end θ (l, τ ), and the resultant inclination β(τ ) On the other hand, the normalized displacements at the free end, ξ (τ ) and η(τ ), follow from the displacement field for the quasi-static problem, equation (4.1), and from the axial compliance of the clamp, equation (7.5), as and The substitution of the expressions (7.8) and (7.10) for μ(τ ), ξ (τ ), η(τ ) in the differential equations (7.7) leads to a nonlinear second-order differential implicit system for the evolution of the free end rotation θ(l, τ ) and for the resultant inclination β(τ ). For given initial conditions and rotation law α(τ ), the integration of the differential system (7.7) can finally be numerically performed.
In the present analysis, initial rest condition for the system and a linear time evolution law for the clamp rotation are considered, so that α(t) = α 0 + ωt, (7.11) where α 0 is the initial clamp rotation and ω is the clamp angular velocity, and which can be rewritten in terms of dimensionless time as α(τ ) = α 0 + Ωτ , (7.12) with Ω being the clamp angular velocity referred to the dimensionless time, Ω = ωT. To avoid a trivial solution, the evolution of the system was analysed starting at τ = 0 from an 'almost' undeformed state for the rod, namely the initial clamp inclination was set to be α 0 = 10 −5 . The rest condition for the system at τ = 0 implies that initially the inclination β(0) coincides with the clamp inclination α(0), so that the force resultant R(0) momentarily coincides with the dead load The four conditions (7.13) imply four initial conditions for the resultant inclination β(τ ) and the free end rotation θ (l, τ ) β(0) = α 0 , θ (l, 0) = θ l0 ,β(0) =β 0 andθ(l, 0) =θ l0 , (7.14) where the values of θ l0 ,β 0 andθ l0 can be computed using equations (7.8) and (7.10). The evolution of the resultant inclination β(τ ) and the free end rotation θ(l, τ ) are obtained from the numerical integration of the nonlinear differential implicit equations (7.7) through the function NDSolve in Mathematica c (v. 10). Considering the clamp angular velocity ω = 0.014 rad s −1 and assuming a small normalized damping coefficient υ = 0.069 and high normalized axial stiffness κ c = 310.647, the numerical integration has been performed at varying the normalized mass parameter, namely for ζ = {1.157, 1.987, 2.818, 3.648, 4.478} corresponding to P/P cr = {0.469, 0.805, 1.142, 1.479, 1.815} and which are representative of some of the experimental set-ups considered in §8.
The massless model is able to capture reasonably well the dynamics of the system for values of the normalized mass parameter ζ lower than 5, corresponding to P/P cr < 2. For higher values of ζ , the numerical integration fails to converge soon after the snap-back instability occurs, namely, the free end experiences a very fast and large oscillation in its acceleration, and the Mathematica solver reveals that the differential system becomes stiff. The encountered difficulty in the numerical treatment is related to the fact that rod configuration is imposed to assume the first quasi-static mode during the dynamics, a condition which becomes not representative at values P/P cr > 2. Indeed, it is experimentally shown in §8 that for this loading condition the dynamics of the rod is characterized by both transverse and longitudinal oscillation, whereas the present simplified model is reliable only when transverse oscillations prevail in the dynamic response.

(b) Finite-element analysis
A finite-element analysis of the system was performed in Abaqus with the purpose to provide a full description of the dynamics of the rod for every loading condition, thus overcoming the limitations found with the simplified model. The rod is modelled through 100 linear elastic planar beam elements, where the first element has the external edge constrained by the rotating clamp while the lumped mass m is attached on the external edge of the final element. Introducing energy dissipation through Rayleigh damping, large displacement analyses were performed through the following two steps: Step 1. Static-the gravitational force is applied in order to establish the deformed initial configuration. Here a quasi-static deformation is produced by the weight of both the elastic rod and the lumped mass (although the effect of the former is much smaller than the effect induced by the latter); Step 2. Dynamic implicit-a slow rotation to the clamped end is imposed as a boundary condition. The analysis takes into account the inertial forces generated during the clamp rotation.
Simulations were run considering geometry, material, inertia and clamp angular velocity representative of the experimental set-up described in the next section. Rayleigh damping was set through the mass-proportional damping coefficient A = 10 −2 s −1 and the stiffness-proportional damping coefficient B = 5 × 10 −3 s. Several low rotation velocities for the clamp have been investigated (ω = {0.014, 0.06, 0.14} rad s −1 ), basically showing no influence of small ω on the dynamics of the system.

Experiments versus modelling for the elastica compass and elastica catapult
An experimental set-up was designed and manufactured at the Instabilities Lab of the University of Trento (http://ssmg.unitn.it) for the analysis of the rotating clamped elastic rod (figure 13). A rack and pinion actuator (figure 13c) was used to transform the linear motion of an electromechanical testing machine (ELE Tritest 50 from ELE International Ltd) into the rotation of the clamp at required angular speeds (different low clamp velocities were tested in the experiments, without noting substantial differences, so that only results for ω = 0.014 rad s −1 are reported). The angular position of the rotating clamp was measured through a contactless rotary position sensor (NRH280DP/180/360). The moment transmitted to the clamp was measured with a lever system connected to a load cell (Leane DBBSN, RC 10 kg), as evident from figure 13c.
During the tests, the modulus of the acceleration of the lumped mass, attached at the end of the beam a(l, t) = ẍ(l, t) 2 +ÿ(l, t) 2 was obtained mounting two miniaturized mono-axial accelerometers (352A24, PCB piezotronics, 0.8g) perpendicular to each other. All the data were acquired with a NI compactRIO system interfaced with Labview 2013 (from National Instruments).
Two rods were used, both made up of a solid polycarbonate strip (white 2099 Makrolon UV from Bayer, E = 2350 MPa, Poisson's ratio ν = 0.37 and volumetric mass density ρ = 1180 kg m −3 ). One rod, 345 mm long, has a 25 × 3 mm cross section and the other, 695 mm long, has a 25 × 2.85 mm cross section. The masses attached at the end of the two rods were chosen to produce the following values of P/P cr = {0.469, 0.805, 1.142, 1.479, 1.815    function of the clamp rotation α for the case of the elastica compass (P/P cr = 0.805, where a quasistatic path is followed, figure 14) and in the case of the elastica catapult (P/P cr = 1.142 and P/P cr = 7.795, where the dynamics prevails after snap-back, figure 15a,b, respectively). In the case of the elastica compass the acceleration a(l, t) was found to be negligible (the maximum value was about 10 −2 g), whereas in the case of the elastica catapult the acceleration a(l, t) (reported in figure 15c,d) was found to have a peak of about 0.6g and 55g for P/P cr = 1.142 and P/P cr = 7.795, respectively. These peaks in the acceleration were found to occur just after the snap-back instability. Photos taken at different rotation angles are superimposed in figure 16 for two cases in which the system behaves as an elastica compass (figure 16a, P/P cr = 0.805) and an elastica catapult (figure 16b, P/P cr = 7.795). In the latter case, the photo labelled 12 corresponds to the last photo taken before snap-back instability (α = 1.744π ), whereas photos 13 and 14 are blurred because they were taken (using a NEX-5N Sony camera) during post-snap dynamics. Note that the shutter speed was set on 1/20 s, which provides a measure of the velocity of the rod during snap-back. In figure 16, the position of the attached lumped mass is compared with the theoretically predicted rod's end trajectory (obtained from the quasi-static analysis and reported as yellow line), showing a very nice agreement.
Movies of the experiments, taken by both high-speed camera Sony PXW-FS5K (240 fps) and a Sony Handycam HDR-XR550, are available as electronic supplementary material.

Conclusion
In the simplest set-up for a flexible robot arm, an elastic rod is subject to a prescribed slow rotation at one end and to a concentrated mass in a gravity field at the other. Solving this system through the elastica has evidenced a discriminating role for the applied load. When the load is smaller than the buckling load (for the straight configuration), the loaded end describes a closed continuous curve (not far from an ellipse), whereas when it is higher the continuous displacement of the loaded end terminates in an unstable geometrical form, from which a snap-back instability leads to a non-adjacent configuration. In the former case, the system behaves as an 'elastica compass' (so that the loaded end of the arm would describe a circle if the rod would be rigid), whereas in the latter case, the system behaves as an 'elastica catapult' (so that the mass would be thrown away if detachment were allowed). The dynamic motion after the snap instability has been analysed with a standard finite-element program and an ad hoc developed software. The quasi-static and dynamical behaviours have been validated with a specifically designed experimental set-up. Results show that the basic problem of soft robot arm addressed in this article, and systematically explored in terms of mechanical performances, can be analytically and numerically described with a great accuracy and therefore is ready for exploitation in real devices. Authors' contributions. All the authors designed the research, co-wrote the paper, carried out the mathematical analysis, carried out the experimental tests and gave final approval for publication.