- Chair of Automatic Control Engineering (LSR), Department of Electrical and Computer Engineering, Technical University of Munich (TUM), Munich, Germany

The potential of large elastic deformations in control applications, e.g., robotic manipulation, is not yet fully exploited, especially in dynamic contexts. Mainly because essential geometrically exact continuum models are necessary to express these arbitrarily large deformation dynamics, they typically result in a set of nonlinear, coupled, partial differential equations that are unsuited for control applications. Due to this lack of appropriate models, current approaches that try to exploit elastic properties are limited to either small deflection assumptions or quasistatic considerations only. To promote further exploration of this exciting research field of large elastic deflection control, we propose a geometrically exact, but yet concise a beam model for a planar, shear-, and torsion-free case without elongation. The model is derived by reducing the general geometrically exact the 3D Simo–Reissner beam model to this special case, where the assumption of inextensibility allows expressing the couple of planar Cartesian parameters in terms of the curve tangent angle of the beam center line alone. We further elaborate on how the necessary coupling between position-related boundary conditions (i.e., clamped and hinged ends) and the tangent angle parametrization of the beam model can be incorporated in a finite element method formulation and verify all derived expressions by comparison to analytic initial value solutions and an energy analysis of a dynamic simulation result. The presented beam model opens the possibility of designing online feedback control structures for accessing the full potential that elasticity in planar beam dynamics has to offer.

## 1 Introduction

Dynamic robotic manipulation of highly deformable objects is still a rarely considered field in literature due to a lack of appropriate dynamic models. Elastic objects in robotic manipulation are usually either considered only quasistatically (Bretl and McCarthy, 2014) or dynamically under the assumption of small deflections which in return results in applications where elasticity is often treated as an undesired property that needs to be avoided or compensated (Tavasoli, 2015). There are some approaches that suggest exploiting elasticity instead, e.g., in terms of safety (Bicchi and Tonietti, 2004; Haddadin et al., 2009), payload estimation (Malzahn et al., 2015), or energy storage (Tantanawat and Kota, 2007; Haddadin et al., 2011). Some authors exploit elasticity also for quasistatic manipulation (McCarragher, 2000) or even to design new strategies for dynamic manipulation, e.g., Pekarovskiy et al. (2014) or Pham and Pham (2017). Nonetheless, we propose that elastic dynamics have even more potential that can be exploited when large deflections are considered.

Models that allow taking advantage of such arbitrarily large deformations need to fulfill two main criteria. They have to

1. describe the dynamics *geometrically exact*, i.e., independently of the magnitude of deformation

2. be concise and simple enough to allow the design of online feedback control

However, such models are not yet available in literature. In this study, we start filling this gap for the case of a planar elastic beam undergoing pure bending, by proposing a single-dimensional dynamic continuum model.

### 1.1 Related Work

Large deflection dynamics are typically treated with one of two approaches—discretized approximations using a finite element method (FEM) or continuum models expressed as partial differential equations (PDEs). The FEM-based descriptions, for which recent literature such as the study of Duriez (2013) reaches real-time capable control, are limited to quasistatic deformable structures, as the computational cost of FEM descriptions for true dynamic large deformation online feedback control is still out of reach. Whereas, the geometrically exact continuum models result in multivariate and highly coupled nonlinear PDE systems, and control approaches are thus limited to oscillation damping (Ito, 2001; Hegarty and Taylor, 2012). Unlike existing control methods for linear beam models (Krstic et al., 2006a; Krstic et al., 2006b), the literature body on theory of nonlinear PDE systems is still too limited in its applicability for the complex expressions arising in these mechanical system models (Padhi and Ali, 2009). From a control point of view, however, it is well known that continuum models, unlike discretized approximations, do not face so called spillover phenomena that can lead to instabilities due to unmodeled high frequency dynamics (Meirovitch and Baruh, 1983). For these reasons, as a first contribution toward filling the discussed gap in literature, this work proposes a geometrically exact model of planar Euler–Bernoulli beam dynamics for arbitrarily large deflections that admit a surprisingly concise PDE formulation.

Although there are special purpose models for large deflection models, such as bullwhip dynamics by Bernstein et al. (1958) and McMillen and Goriely (2003) or dynamics of a fly fishing line by Spolek (1986), we are interested in more general beam dynamics that admit rope and cable dynamics as a special case of very low elastic stiffness. A vast literature body exists on 1D analytic beam theories alone. Therefore, because of the considerable computational advantage of 1D beam theories over 3D continuum mechanics theories, the latter will not be considered in this concise literature review. Starting with the first mathematical treatment of static elasticity by Galileo already in 1638, Hooke’s treatise of linear elasticity is followed in 1678. The first precise definition of the elastica (a thin strip of elastic material) problem was carried out by Jakob Bernoulli (1691), and he published its first solution in Bernoulli (1694). His nephew Daniel Bernoulli (1742) did not himself solve the problem, but he suggested Euler the use of variational analysis, who delivers a closed-form solution of the elastica in the study by Euler (1744). A more detailed and insightful mathematical historical overview of the elastica can be found in the study by Levien (2008). While Euler’s early work already predicts slender beam deformations with astounding precision, many authors built on this work to include further effects to account for more general conditions as well as geometries. It is said that Rayleigh (1877) added rotary inertia effects, and Timoshenko (1921) further enhanced the theory to account for shear effects. However, Elishakoff (2020) discusses original authors and naming of linear beam theories. He, e.g., mentioned that Bresse (1859) already included rotary inertia effects before Rayleigh, though the works were developed independently. Furthermore, a beam theory including shear effects was originally published in Timoshenko (1916), an earlier book in Russian language, where Timoshenko mentions to have developed the theory together with P. Ehrenfest. Elishakoff, therefore, suggests the historically justifiable name *Bresse–Rayleigh–Timoshenko–Ehrenfest beam theory*.

For large deflections, also referred to as *finite strain*, geometrically nonlinear models are necessary. Kirchhoff (1859) is a spatial generalization of the Euler–Bernoulli beam and allows modeling of 3D deformations through bending and torsion. The theory was later extended by Love (1892) to further account for axial tensions and is referred to as *Kirchhoff–Love beam theory*. Reissner (1972) added further measures to Kirchhoff’s theory, accounting for shear deformations in planar curves and later for space curves in Reissner (1981). Simo (1985) enhanced Reissner’s work in terms of approximations to what is nowadays known as *Simo–Reissner beam theory* or *geometrically exact beam theory*. To complete the overview, it is also worth mentioning that reduced versions of these two well-known theories have been proposed that neglect torsion modes, e.g., the study by Meier et al. (2015) for the Kirchhoff–Love and Romero et al. (2014) for the Simo–Reissner case. Meier et al. (2019) gave a more in-depth overview and analysis of these nonlinear beam theories. Table 1 shows how our proposed model fills the current literature gap of a concise model for geometrically exact descriptions.

**TABLE 1**. Comparison of this work to the most commonly used beam theories. While the proposed model is limited to planar bending, it offers a geometrically exact description in a single-dimensional equation.

### 1.2 Contribution

The main contribution of this work is twofold. We provide

• the first single-dimensional, geometrically exact, dynamic beam model, and

• a method for incorporating boundary conditions in a FEM formulation, for cases where the FEM model has only descriptive variables of higher-order derivatives than the boundary condition itself

### 1.3 Outline

In Section 2, we derive the model via step-by-step reduction, starting from the general Simo–Reissner beam theory, extracting first the Kirchhoff–Love beam theory, followed by special case assumptions — isotropic, torsion-free, inextensible, and planar bending. The translation into a FEM description is explained in Section 3, including our proposal for incorporating position-level boundary conditions into the tangent angle PDE system. This FEM description is applied for the simulation verification in Section 4. The work is concluded in the final Section 5.

## 2 Modeling

This section outlines the model derivation, starting from a general 3D theory. After reducing the model by gradually introducing assumptions, we couple the Cartesian coordinates to achieve an expression in a single PDE.

### 2.1 Nomenclature

This work follows the convention of using lowercase bold variables for vectors and uppercase bold variables for matrices. All nonbold variables are scalars. Furthermore, subscript annotations are reserved for index notation of multidimensional variables as well as expressing partial derivatives, whereas superscript annotations are part of the variable specification. Also, note that we omit explicit listing of function parameters whenever it is clear from the context to not unnecessarily clutter the notation. A list of the most frequently used variables at the end of this work is given.

### 2.2 Model Reduction

The presented geometrically exact model for a planar Euler–Bernoulli beam is found via reduction of the general 3D Simo–Reissner beam theory. We gradually introduce further assumptions to simplify the dynamic governing equations of the beam. Please note that we only define the individual components of the equations once necessary to keep the derivation clear and easy to follow. The resulting model forms a special case of a planar Kirchhoff–Love beam theory parametrized solely in the curve tangent angle, in an analog to Euler’s elastica (Euler, 1744). A more in-detail analysis and discussion of the general Simo–Reissner and Kirchhoff–Love beam theories, including the special cases of isotropic cross-sections and torsion-free formulations, can be found in Meier et al. (2019).

##### General Simo–Reissner Beam Theory

The first theory that accounts for very general 3D beam deformations including spatial bending, torsion, axial tension, and shear deformation was published by Simo (1985). The consideration of shear effects makes it an adequate theory for thick rod dynamics. It is also been denoted by *geometrically exact beam theory* because the description is consistent at the deformed state regardless of the magnitude of displacements, rotations, and strains, cf. Crisfield and Jelenić (1999). Simo himself also used the term *finite strain beam formulation*.

The strong form of the Simo–Reissner beam theory, cf. Simo (1985), is a system of six coupled PDEs and can be stated with the equilibrium equations

where the internal force * r* of the beam center line. The quantities

In the following, all objective deformation measures are chosen to be work conjugated to the material stress resultants in Eq. 1. We further assume a hyperelastic constitutive relation between these kinetic and kinematic quantities.

##### Assumption: Vanishing Shear Strains (Kirchhoff–Love Beam Theory)

Neglecting shear deformations and assuming that the cross-section is always perpendicular to the center line of the beam, the change in internal forces

where Eq. 2b is now a scalar expression, and the beam model is thus reduced to 4 PDEs.

##### Assumption: Initially Straight and Isotropic

If now an initially straight beam with an isotropic cross-section is assumed for the hyperelastic beam, the components of Eq. 2 are given with

with the Young modulus *E*, inertia *I*, density *ρ*, cross-section area *A,* and the curvature vector

The components

##### Assumption: Torsion-Free

Assuming pure bending and no torsional effects, the inertia moment and the moment balance Eq. 3b vanish completely. Moreover, only the perpendicular component of the external moment affects the force balance equation

with the axial tension parameter

of the undeformed beam, does not in general hold for the deformed case

due to possible elongations in the beam structure.

##### Assumption: Inextensible Beam

If it is assumed that the beam does not undergo axial elongations, the gradient Eq. 7 does always equals 1. Hence, the axial tension parameter *ϵ* evaluates to

and thus, ^{∥}

However, unlike the previous assumptions that can be incorporated implicitly with an adequate choice of parametrization variables, it is in general difficult to find such a set of variables that fulfill the inextensibility constraint Eq. 8a by construction. A common practice to enforce the equality constraint Eq. 8a on the simulation result in a weak sense, i.e., in an integral form instead of point-wise, is by means of extending the weak form of the model Eq. 8b with a Lagrange multiplier potential, cf. Meier et al. (2019).

The following assumption of pure planar bending, however, does again permit a parametrization that fulfills this constraint directly in the strong sense, i.e., for every point along the beam.

##### Assumption: Pure Planar Bending

The last step in the model reduction is the restriction to pure planar bending. For the remainder of this section, we will switch to a component-wise notation in Cartesian coordinates. The beam model from Eq. 8b reduces to two coupled PDEs and is fully described by

and the additional inextensibility constraint

because the third row of the vector equation Eq. 8b vanishes for the purely planar problem. Thus, the remaining external inputs are forces *xy*-plane, as well as the external moment

With the curve tangent angle

for the shear-free, planar, inextensible case. Furthermore, it allows stating the geometric identities

The planar beam model Eq. 9 can thus be stated as

a beam model in a mixed form containing Cartesian coordinates as well as the curve tangent angle as describing variables. In the remainder of this Section, Eq. 12 is the starting point to first develop the static case followed by the general dynamic case, both entirely expressed in the curve tangent angle.

##### 2.2.1 Static Beam Model Expressed in the Curve Tangent Angle

Only considering solutions in a static equilibrium, the Cartesian acceleration terms in Eq. 12 vanish and

remains. By computing all derivatives,

and rotating the equations from their Cartesian *φ* around the *z*-axis, by premultiplying Eq. 14 with the rotation matrix

which allows extracting the components perpendicular

The fact that no Cartesian

which can be integrated twice and yields a unique solution if boundary conditions are applied.

##### 2.2.2 Dynamic Beam Model Expressed in the Curve Tangent Angle

To also fully state the dynamic planar beam model in terms of the curve tangent angle, the Cartesian *φ*. This is achieved by differentiating the system of equations Eq. 12 w.r.t. the beam parameter *l*. Assuming no buckling of the object, *x* and *y* have continuous derivatives, and thus, Schwarz’s theorem allows changing the order of the derivations. Applying the geometric identities (Eq. 11) now also to the acceleration terms, leads to

with the auxiliary variable

a PDE system entirely expressed in terms of the curve tangent angle *φ*. As for the static case (Eq. 16), the geometric identities fulfill the inextensibility constraint Eq. 9b by construction; thus, no special consideration is necessary. Expanding all partial derivatives and grouping the trigonometric terms, Eq. 18a yields

Similar to the static case, premultiplying the entire system with the rotation matrix

which in the case of no external inputs admits the very concise strong form

as a single PDE governing the beam dynamics in a single parameter *φ*.

While this reduced model is relevant for PDE controller development, it is not directly applicable for use in simulations. We therefore present in the following section a respective approximation with a system of ordinary differential equations (ODEs), in terms of a FEM formulation.

## 3 FEM Formulation

In this section, we outline the development of a FEM simulation procedure, starting from the development of the weak form of the beam model (Eq. 21), without considering external forces. After transforming the integrodifferential weak form into a system of nonlinear ODEs of thee second order in time via a Bubnov–Galerkin approximation, a finite element discretization leads to a simulation procedure.

### 3.1 Weak Form of Large Deformation in the Curve Tangent Angle

The weak form of Eq. 21 is found by multiplying the equation with the test function

A sequence of integrations by parts will lead to the final weak form. In the first intermediate step,

allows to solve the integral

Using this result, the terms on the right-hand side of Eq. 22 result in

and

which lead to the final dynamic equations in a weak form

that builds the basis for the following FEM formulation. The function *φ* as well as the variation

with the function space of square integrable functions

such that they are twice continuously differentiable in *l*, cf. Reddy (2013). A common method to choose candidates for *φ* and its variation

### 3.2 Bubnov–Galerkin Approximation

In the sense of the Bubnov–Galerkin method, the function *φ* as well as the test function

using the same set of *n* weighted orthogonal spatial basis functions *n* time-dependent scaling coefficients

As the coefficients *n* equations

for

with the component-wise definitions

Note the index order of the matrix component definitions

### 3.3 Finite Element Discretization

The core idea of the FEM is to discretize the continuous body into a finite number of *N* elements

For the weak formulation (Eq. 24), the orthogonal basis functions need to be at least of class

They have the special property that those coefficients that build the spline element,

directly correspond to the function value *φ* and its derivative

**FIGURE 2**. FEM nodes and elements. The local basis functions (Eq. 29) of a unit element are shown on the left. The right plot shows the assembly of *φ* and

### 3.4 Boundary Conditions Expressed in the Curve Tangent Angle

Incorporate the boundary conditions on the curve tangent angle *φ*, and its derivatives follow standard procedures, e.g., Hughes (2012). Because treating position-based boundary conditions is not directly possible in the curve tangent angle beam model (Eq. 21), this chapter focuses on the development of a strategy to express boundary conditions in terms of higher order derivatives only. Without loss of generality, we consider for our beam model a fixed end,

and, respectively, its dynamic counterpart

Fixing an elastic beam in its position at one end introduces point-wise reaction forces from the mounting onto the beam in *x* and/or *y* directions. While these boundary conditions are straight forward to be incorporated in a FEM formulation for a model in the parameters *x* and *y*, e.g., Eq. 9a, the FEM description of the reduced model (Eq. 28a) directly acts on the tangent angle function *φ* and its derivatives, thus not offering any parameter to incorporate position boundary conditions. However, except for the a-priori known position at the boundary condition itself, the curve tangent angle dynamics Eq. 21 do not require any position parameters to govern the beam profile.

We thus propose to transform the position boundary condition at node

Formulating the weak forms and applying another integration by parts reads

and the Bubnov–Galerkin approximation

with the same set of orthogonal functions ψ that leads to the systems of equations

with the components

Note that for both systems of equations in Eq. 36, the matrices *φ* from Eq. 28a.

Again using Hermite cubic spline basis functions (Eq. 30), the equations of interest in the two systems of *N* equations (Eq. 36) are the ones relating to the position basis function of

As illustrated in Figure 2, this requires an integration from

relating to the function values and first derivatives at the first three nodes. While the right hand side of the FEM formulations (Eq. 36) is already defined in the curve tangent angle *φ*, what remains is to also rewrite coefficient vectors *φ* instead of *x* and *y*. Starting again from the geometric identities

containing Volterra integrals with an upper limit *l*. Recalling the spline approximation

allows to express the entire FEM balance equations (Eq. 38) in terms of *φ* with the left hand sides

and

where ^{○}, denotes an element-wise product, considering the piecewise function element definitions *x* direction of *φ* reads

and the load value in *y* direction reads

These expressions are used to incorporate the position boundary conditions into the FEM formulation of the curve tangent angle beam equation.

The resulting system of nonlinear ODEs is used to verify the proposed beam model as well as the presented strategy for incorporating position boundary conditions.

## 4 Simulation Verification

The proposed beam model of Section 2 together with the method for incorporating boundary conditions on lower level derivatives than the descriptive variables of the FEM formulation from Section 3.2 is verified in simulation. First, the developed FEM strategy is tested in two different initial value problems and demonstrated for different magnitudes of external forces and momenta. The second part verifies a dynamic simulation in terms of energy consistency of the beam profile. For both cases, we consider the case of a clamped end at

While the essential boundary condition Eq. 43c and natural boundary conditions Eq. 43d and Eq. 43e are directly considered in the FEM formulation with conventional techniques, conditions Eq. 43a and Eq. 43b are incorporated with the expressions developed in Section 3.4. All simulations have been conducted in MATLAB R2020a.

### 4.1 Initial Value Problem

The static solutions of the FEM formulation (Eq. 28a) with external nodal forces

are verified in two scenarios. First, the analytic solution of Eq. 17 is replicated with an external moment at the free end, and in the second example, an external nodal force

where the last two rows impose the boundary conditions Eq. 43a and Eq. 43b using the expressions from Eq. 42a and Eq. 42b. The results show that stem from MATLAB’s nonlinear least-squares solver lsqnonlin() for a nominal beam of

##### 4.1.1 External Moment

Applying an external nodal moment at the free end of the beam results in

and thus a constant curvature *φ* is a purely linear function in this special case and thus can be approximated with arbitrary accuracy, even for a low number of elements. Note that the analytical trigonometric solution is not visually distinguishable from the simulation result and thus is not shown in Figure 3.

**FIGURE 3**. Results of the initial value problem, simulated in the beam tangent angle *φ* for different external momenta at the free end *φ* profiles are shown on the left. Additionally, the middle and right columns show the resulting *x* and *y* profiles (dashed lines). The large red dots mark the imposed boundary conditions on the curve tangent angle profile and their respective impact in *x* and *y* direction. Note that the *x* and *y* profiles, as well as the 2D visualization on the top left, are evaluated during postprocessing of the actual simulation results in the curve tangent angle *φ*, applying the geometric identities Eq. 11 and respective derivatives.

##### 4.1.2 External Force

This second example shows the beam deflection under different external nodal forces *y* direction correspond to results of a common linearized Euler–Bernoulli Beam model, the geometrical nonlinearities show the full effect for the two larger deflection cases. Note that the initial curvature at the clamped end

**FIGURE 4**. Results of the initial value problem simulated in the beam tangent angle *φ* for different external nodal forces at *φ* profiles are shown on the left. A detailed explanation of the illustration is given in the caption of Figure 3.

### 4.2 Dynamic Simulation

For the dynamic case, the FEM formulation (Eq. 28a) is again extended by external nodal forces

which allow external perturbation of the system. Assuming that the initial beam configuration

where the position boundary conditions from Eq. 42a are encoded in

with the first six columns of

The resulting acceleration vector of curve tangent angle coefficients *t*.

##### 4.2.1 Energy Analysis

The dynamic simulation result is verified by confirming energy consistency of the simulation result. This is achieved by calculating the total energy,

consisting of the potential energy

It can be directly evaluated from the simulation results. The kinetic energy

however, needs to be rewritten in terms of *φ* first. The geometric identities *l* of the beam center line, with the Volterra integrals

and thus, the kinetic energy of the entire beam can be evaluated with

Note that the rotational component of the kinetic energy is not considered, due to the Euler–Bernoulli assumption of slender beams.

##### 4.2.2 Simulation Results

For the dynamic case, an initially unloaded beam with a normalized beam material parameter

and again at the beam center, with a maximal unit force

**FIGURE 5**. Energy conservation during an FEM simulation of the dynamic curve tangent angle model in a clamped-free scenario, where an external torque is applied to the center of an initially straight beam for 0.5 s. The left plot shows the energy distribution. The right plot shows the energy along the beam at a snapshot taken at *t* = 1 s.

## 5 Conclusion

A model for large planar deformation dynamics of Euler–Bernoulli beams was presented and put into context with well-known more general beam models. Literature does already offer various models that account for arbitrarily large deformations; however, they typically result in a system of coupled nonlinear PDEs expressions. Whereas, the presented approach admits a single-dimensional PDE in one variable, i.e., the curve tangent angle of the beam center line to describe planar beam dynamics under the common Euler–Bernoulli assumptions of shear-free constant cross-sections.

While boundary conditions on the beam profile derivatives—which is sufficient for sliding and/or free ends—can be directly encoded in a simulation algorithm of the curve tangent angle beam model, there is no descriptive variable available to directly incorporate boundary conditions on the beam position—needed for clamped and/or hinged ends. These cases are, however, of course of highly practical relevance. To also address these cases, we additionally outlined a novel method that allows incorporating boundary conditions in FEM formulations, using solely descriptive variables of higher order derivatives. We apply this method to impose position-based boundary conditions in the FEM formulation of the presented beam model, but it is not limited to solely use this case. The strategy is verified in initial value problems, where it replicates analytical solutions and a time-variant FEM simulation by evaluating energy conservation. The presented beam model does not only reduce the computational effort due to the dimensional reduction to a single parameter, but the beam profiles are at the same time less complex in this curve tangent angle parametrization and thus require fewer elements in the FEM description.

Although being nonlinear, the derived beam model provides a concise continuum model for future control theory applications of large deformations, where other more complex model descriptions are not appropriate for current model-based PDE controller development. The model reduction process in this work, however, also outlines more general beam models considering, e.g., shearing, axial torsion, elongation, and/or 3D spatial deformations. This work is thus also a good source for extended beam models including these additional effects, which might be relevant with the advancement of PDE controller development.

## Data Availability Statement

The original contributions presented in the study are included in the article/Supplementary Material, further inquiries can be directed to the corresponding authors.

## Author Contributions

GH, DW and MB initiated the study. GH formulated, implemented and conducted the research. All authors contributed to manuscript revision, read, and approved the submitted version.

## Funding

The research leading to these results has received funding from the Horizon 2020 research and innovation program under grant agreement

## Conflict of Interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

## Acknowledgments

The authors would like to thank Christoph Meier for valuable discussions on the methodology and FEM formulation, as well as Armin Lederer for constructive discussions and implementations at an early stage of the project.

## Nomenclature

#### Beam parameters

*A* Area of cross-section surface

*E* Young’s modulus

*G* Shear modulus

*L* Total length

ρ Density

*c* Constant scaling factor of the PDE

#### Vectors

** f** Force vector

** m** Moment vector

** r** Vector from the origin to a point on the deformed beam

#### Superscripts

^{║}( ) Parallel component

^{┴}( ) Perpendicular component

#### Modeling variables

φ Curve tangent angle between deformed

*l* Beam parameter along beam axis with range

*t* Time

*x* First component of the location vector along

*y* Second component of the location vector along

#### Finite element description

*λ* Local length parameter with range

*N* Number of elements

*e* Element number with range

#### Operators

*x*, i.e.,

## References

Antman, S. S. (2005). Planar steady-state problems for elastic rods. *Nonlinear problems of elasticity*. 2 edn. New York, United States: Springer Science & Business Media Chap. 4, 95–134.

Bernoulli, J. (1691). Quadratura curvae, e cujus evolutione describitur inflexae laminae curvatura. *Die Werke von Jakob Bernoulli* 5, 223–227.

Bernstein, B., Hall, D. A., and Trent, H. M. (1958). On the dynamics of a bull whip. *The J. Acoust. Soc. America* 30, 1112–1115. doi:10.1121/1.1909473

Bicchi, A., and Tonietti, G. (2004). Fast and “Soft-Arm” tactics. *IEEE Robot. Automat. Mag.* 11, 22–33. doi:10.1109/mra.2004.1310939

Bresse, J. A. C. (1859). *Cours de mecanique appliquee: Re’sistance des mate’riaux et stabilite’des constructions*, (Paris, France: Mallet-Bachelier). 1.

Bretl, T., and McCarthy, Z. (2014). Quasi-static manipulation of a Kirchhoff elastic rod based on a geometric analysis of equilibrium configurations. *Int. J. Robotics Res.* 33, 48–68. doi:10.1177/0278364912473169

Crisfield, M. A., and Jelenić, G. (1999). Objectivity of strain measures in the geometrically exact three-dimensional beam theory and its finite-element implementation. *Proc. R. Soc. Lond. A.* 455, 1125–1147. doi:10.1098/rspa.1999.0352

Duriez, C. (2013). Control of elastic soft robots based on real-time finite element method. In International Conference on Robotics and Automation *(ICRA) (IEEE)*, 3982–3987. doi:10.1109/ICRA.2013.6631138

Elishakoff, I. (2020). Who developed the so-called timoshenko beam theory? *Math. Mech. Sol.* 25, 97–116. doi:10.1177/1081286519856931

Euler, L. (1744). *Methodus inveniendi lineas curvas maximi minimive proprietate gaudentes sive solutio problematis isoperimetrici latissimo sensu accepti*. Genevae, Switzerland: Marcus-Michael Bousquet), 1.

Haddadin, S., Krieger, K., and Albu-Schäffer, A. (2011). Exploiting elastic energy storage for cyclic manipulation: modeling, stability, and observations for dribbling. *(CDC-ECC) (IEEE)*, 690–697. doi:10.1109/CDC.2011.6161022

Haddadin, S., Laue, T., Frese, U., Wolf, S., Albu-Schäffer, A., and Hirzinger, G. (2009). Kick it with elasticity: safety and performance in human-robot soccer. *Robotics Autonomous Syst.* 57, 761–775. doi:10.1016/j.robot.2009.03.004

Hegarty, G., and Taylor, S. (2012). Classical solutions of nonlinear beam equations: existence and stabilization. *SIAM J. Control. Optim.* 50, 703–719. doi:10.1137/100793864

Hughes, T. J. (2012). *The finite element method: linear static and dynamic finite element analysis*. (Chelmsford, MA: Courier Corporation).

Ito, K. (2001). Numerical study for the stability of a geometrically nonlinear elastic beam with velocity feedback. *Nonlinear Anal. Theor. Methods Appl.* 47, 3813–3821. doi:10.1016/s0362-546x(01)00500-4

Kirchhoff, G. (1859). Ueber das gleichgewicht und die bewegung eines unendlich dünnen elastischen stabes. *J. für die reine Angew. Mathematik* 1859, 285–313.

Krstic, M., Siranosian, A. A., and Smyshlyaev, A. (2006a). Backstepping boundary controllers and observers for the slender timoshenko beam: Part Ⅰ-design. American control conference. (ACC) (IEEE), Minneapolis, MN, 24 July, 2006, 2412–2417. doi:10.1109/ACC.2006.1656581

Krstic, M., Siranosian, A. A., Smyshlyaev, A., and Bement, M. (2006b). Backstepping boundary controllers and observers for the slender timoshenko beam: Part Ⅱ—stability and simulations. Conference on Decision and Control (CDC) (IEEE), San Diego, CA, 13-15 December, 2006, 3938–3943. doi:10.1109/CDC.2006.377717

Levien, R. (2008). *The elastica: a mathematical history*. Berkeley, CA: University of California. Technical Report No. UCB/EECS-2008-103.

Love, A. E. H. (1892). *A treatise on the mathematical theory of elasticity*. Cambridge, United Kingdom: The University Press.

Malzahn, J., Schloss, R., and Bertram, T. (2015). Link elasticity exploited for payload estimation and force control. International Conference on Intelligent Robots and Systems (IROS) (IEEE), Hamburg, Germany, 28 September–2 October, 2015, 1508–1513. doi:10.1109/IROS.2015.7353567

McCarragher, B. J. (2000). “A hybrid position/force approach to the exploitation of elasticity in manipulation,” in *Robot manipulation of deformable objects*. Editors D. Henrich, and H. Wörn (London, United Kingdom: Springer London), 91–109. doi:10.1007/978-1-4471-0749-1_7

McMillen, T., and Goriely, A. (2003). Whip waves. *Physica D: Nonlinear Phenomena* 184, 192–225. doi:10.1016/s0167-2789(03)00221-5

Meier, C., Popp, A., and Wall, W. A. (2015). A locking-free finite element formulation and reduced models for geometrically exact Kirchhoff rods. *Comp. Methods Appl. Mech. Eng.* 290, 314–341. doi:10.1016/j.cma.2015.02.029

Meier, C., Popp, A., and Wall, W. A. (2019). Geometrically exact finite element formulations for slender beams: Kirchhoff-Love theory versus simo-reissner theory. *Arch. Computat Methods Eng.* 26, 163–243. doi:10.1007/s11831-017-9232-5

Meirovitch, L., and Baruh, H. (1983). On the problem of observation spillover in self-adjoint distributed-parameter systems. *J. Optim Theor. Appl* 39, 269–291. doi:10.1007/bf00934533

Padhi, R., and Ali, S. F. (2009). An account of chronological developments in control of distributed parameter systems. *Annu. Rev. Control.* 33, 59–68. doi:10.1016/j.arcontrol.2009.01.003

Pekarovskiy, A., Saluja, K., Sarkar, R., and Buss, M. (2014). Resonance-driven dynamic manipulation: dribbling and juggling with elastic beam. International conference on robotics and automation (ICRA) (IEEE), Hong Kong, China, 31 May–7 June 2014, 943–948. doi:10.1109/ICRA.2014.6906967

Pham, H., and Pham, Q.-C. (2017). Robotic manipulation of a rotating chain. *IEEE Trans. Robotics* 34, 139–150. doi:10.1109/TRD.2017.2775650

Reddy, B. D. (2013). *Introductory functional analysis: with applications to boundary value problems and finite elements*, 27. New York, NY, United States: Springer Science & Business Media.

Reissner, E. (1981). On finite deformations of space-curved beams. *Z. Angew. Math. Phys.* 32, 734–744. doi:10.1007/bf00946983

Reissner, E. (1972). On one-dimensional finite-strain beam theory: the plane problem. *J. Appl. Math. Phys. (Zamp)* 23, 795–804. doi:10.1007/bf01602645

Romero, I., Urrecha, M., and Cyron, C. J. (2014). A torsion-free non-linear beam model. *Int. J. Non-Linear Mech.* 58, 1–10. doi:10.1016/j.ijnonlinmec.2013.08.008

Simo, J. C. (1985). A finite strain beam formulation. The three-dimensional dynamic problem. Part I. *Comp. Methods Appl. Mech. Eng.* 49, 55–70. doi:10.1016/0045-7825(85)90050-7

Spolek, G. A. (1986). The mechanics of flycasting: the flyline. *Am. J. Phys.* 54, 832–836. doi:10.1119/1.14425

Tantanawat, T., and Kota, S. (2007). Design of compliant mechanisms for minimizing input power in dynamic applications. *J. Mech. Des.* 129, 1064–1075. doi:10.1115/1.2756086

Tavasoli, A. (2015). Dynamic modeling and nonlinear boundary control of hybrid Euler-Bernoulli beam system with a tip mass. *Proc. IMechE* 229, 3–15. doi:10.1177/1464419314541790

Timoshenko, S. P. (1916). *A course of elasticity theory. Part 2: rods and plates*. St. Petersburg, Russia: A.E. Collins Publishers.

Keywords: Kirchhoff–Love beam, elastic deformation, large deformation, planar deflection, geometrically exact, PDE, FEM

Citation: Huber G, Wollherr D and Buss M (2021) A Concise and Geometrically Exact Planar Beam Model for Arbitrarily Large Elastic Deformation Dynamics. *Front. Robot. AI* 7:609478. doi: 10.3389/frobt.2020.609478

Received: 23 September 2020; Accepted: 29 December 2020;

Published: 25 May 2021.

Edited by:

Ziwei Liu, The Chinese University of Hong Kong, ChinaReviewed by:

Renate Sachse, University of Stuttgart, GermanyLeopoldo Greco, University of Catania, Italy

Copyright © 2021 Huber, Wollherr and Buss. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.

*Correspondence: Gerold Huber, gerold.huber@tum.de