# Closed-forms of Kirchhoff elastic rods shape and sensitivity in the planar case

Olivier Roussel, Marc Renaud, Michel Taïx

► To cite this version:

Olivier Roussel, Marc Renaud, Michel Taïx. Closed-forms of Kirchhoff elastic rods shape and sensitivity in the planar case. [Research Report] LAAS-CNRS. 2015. hal-01133395v3

HAL Id: hal-01133395

https://hal.science/hal-01133395v3

Submitted on 7 Apr 2015

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution - ShareAlike 4.0 International License
# Closed-forms of Kirchhoff elastic rods shape and sensitivity in the planar case

Olivier Roussel¹, Marc Renaud¹ and Michel Taïx¹*

## Abstract

In this report we give closed-forms of Kirchhoff 3-D elastic rods curvature in terms of elliptic functions and, by treating planar rods as a special case, we show we can also obtain closed-forms of planar rods shape, sensitivity and total elastic energy.

## 1 General case of 3-D rods

Consider an inextensible, non-shearable and unit length linearly elastic rod. The shape of the rod traces a curve that we will describe by the mapping $q: [0, 1] \to SE(3)$. The position along the rod is parametrized by $t \in [0, 1]$ and we will name "base" and "tip" of the rod its extremity at $t = 0$ and $t = 1$ respectively. Let the mappings $u_1(t), u_2(t), u_3(t)$ such that $u_i: [0, 1] \to \mathbb{R}$ be axial and bending rod strains respectively, and $c_1, c_2, c_3$ be the constants that reflect its elasticity properties. As in [3], we say the elastic rod is in static equilibrium in the sense of Kirchhoff if it locally minimizes the elastic energy defined by

$$E_{el} = \frac{1}{2} \int_{0}^{1} \sum_{i=1}^{3} c_i u_i^2 dt. \qquad (1.1)$$

Without loss of generality, we will also assume that the base of the rod is held fixed at the origin, i.e. $q(0) = e$ where $e$ is the identity element of $SE(3)$. Under these assumptions, we will denote by $\mathcal{B}$ the set of positions that the other extremity of the rod $q(1)$ can reach. As shown in [2], the problem of static equilibrium of such rods can be formulated as an optimal control problem by

$$\begin{array}{ll} \underset{q,u}{\text{minimize}} & \frac{1}{2} \int_{0}^{1} \sum_{i=1}^{3} c_i u_i^2 dt \\ \text{subject to} & \dot{q} = q \left( \sum_{i=1}^{3} u_i X_i + X_4 \right) \\ & q(0) = e, \quad q(1) = b \end{array} \qquad (1.2)$$

for some $b \in \mathcal{B}$ and where

$$\begin{aligned} X_1 &= \begin{bmatrix} 0 & 0 & 0 \\ 0 & 0 & -1 \\ 0 & 1 & 0 \\ 0 & 0 & 0 \\ 0 & 0 & 0 \end{bmatrix} & X_2 &= \begin{bmatrix} 0 & 0 & 1 \\ 0 & 0 & 0 \\ -1 & 0 & 0 \\ 0 & 0 & 0 \\ 0 & 0 & 0 \end{bmatrix} & X_3 &= \begin{bmatrix} 0 & -1 & 0 \\ 1 & 0 & 0 \\ 0 & 0 & 0 \\ 0 & 0 & 0 \\ 0 & 0 & 0 \end{bmatrix} \\ X_4 &= \begin{bmatrix} 0 & 0 & 0 \\ 0 & 0 & 0 \\ 0 & 0 & 0 \\ 0 & 0 & 0 \end{bmatrix} & X_5 &= \begin{bmatrix} 0 & 0 & 0 \\ 0 & 0 & 1 \\ 0 & 0 & 0 \\ 0 & 0 & 0 \end{bmatrix} & X_6 &= \begin{bmatrix} 0 & 0 & 0 \\ 0 & 0 & 0 \\ 0 & 0 & 1 \\ 0 & 0 & 0 \end{bmatrix} \end{aligned}$$

is a basis for $\mathfrak{se}(3)$, the Lie algebra of $SE(3)$. Note that when solving this optimal control problem, the rod tip position $b$ is not an input.

<footnote>*¹CNRS, LAAS, Univ. de Toulouse, UPS, 7 avenue du Colonel Roche, F-31400 Toulouse, France</footnote>
In these conditions, the Maximum Principle states that solutions to this optimal control problem are the projections of extremal curves defined on the cotangent bundle $T^*SE(3)$ onto $SE(3)$. Thanks to the Lie Group structure of $SE(3)$, the Hamiltonian can be reduced on the dual of the Lie algebra $\mathfrak{se}(3)^*$ and the corresponding (time-varying) Hamiltonian vector fields $\mu: [0, 1] \to \mathfrak{se}(3)^*$ can be expressed by

$$ \left\{ \begin{aligned} \dot{\mu_1} &= \frac{\mu_3\mu_2}{c_3} - \frac{\mu_2\mu_3}{c_2} \\ \dot{\mu_2} &= \mu_6 + \frac{\mu_1\mu_3}{c_1} - \frac{\mu_1\mu_3}{c_3} \\ \dot{\mu_3} &= -\mu_5 + \frac{\mu_1\mu_2}{c_2} - \frac{\mu_1\mu_2}{c_1} \\ \dot{\mu_4} &= \frac{\mu_3\mu_5}{c_3} - \frac{\mu_2\mu_6}{c_2} \\ \dot{\mu_5} &= \frac{\mu_1\mu_6}{c_1} - \frac{\mu_3\mu_4}{c_3} \\ \dot{\mu_6} &= \frac{\mu_2\mu_4}{c_2} - \frac{\mu_1\mu_5}{c_1} \end{aligned} \right. \qquad (1.3) $$

where vector fields $\mu$ are related to controls $u_i$ by $u_i = c_i^{-1}\mu_i$ for $i \in \{1, 2, 3\}$.

Let $\mathcal{A}$ be the set homeomorphic to $\mathbb{R}^6$ and $a \in \mathcal{A}$ such that $a_i \triangleq \mu_i(0)$, $i \in \{1, \dots, 6\}$. It has been shown in [2] that coordinates in $\mathcal{A}$ offer a global parameterization to the set of static equilibrium configuration for the rod. In other words, we can describe configurations of quasi-static 3-D elastic rods using the 6-dimensional configuration space $\mathcal{A}$.

Assuming isotropy and normalized elasticity constants such that $c_i = 1$ for $i \in \{1, 2, 3\}$, we have from (1.3) $\dot{\mu}_1 = 0$. Then $\mu_1$ is a constant of motion with $\mu_1 = a_1$ and

$$ \left\{ \begin{aligned} \dot{\mu}_2 &= \mu_6 \\ \dot{\mu}_3 &= -\mu_5 \\ \dot{\mu}_4 &= \mu_3\mu_5 - \mu_2\mu_6 \\ \dot{\mu}_5 &= a_1\mu_6 - \mu_3\mu_4 \\ \dot{\mu}_6 &= \mu_2\mu_4 - a_1\mu_5 \end{aligned} \right. \qquad (1.4) $$

The signed curvature $\kappa$ and the torsion $\tau$ of the curve can be expressed in terms of $\mu$ by

$$ \kappa^2 = \mu_2^2 + \mu_3^2 \qquad \tau = \mu_1 - \frac{\mu_2\mu_5 + \mu_3\mu_6}{\mu_2^2 + \mu_3^2} $$

and, as mentioned in [2], the differential system (1.4) is equivalent to

$$ 2\ddot{\kappa} + \kappa^3 - 2\kappa(\tau - \lambda_1)^2 = \lambda_2\kappa \qquad (1.5a) $$

$$ \kappa^2(\tau - \lambda_1) = \lambda_3 \qquad (1.5b) $$

where the constants of integration are given by

$$ \left\{ \begin{aligned} \lambda_1 &\triangleq \frac{a_1}{2} \\ \lambda_2 &\triangleq a_2^2 + a_3^2 + 2a_4 - \frac{a_1^2}{2} \\ \lambda_3 &\triangleq \frac{a_1}{2}(a_2^2 + a_3^2) - (a_2a_5 + a_3a_6). \end{aligned} \right. $$

Substituting (1.5b) into (1.5a) and integrating, we obtain

$$ \dot{\kappa}^2 + \frac{1}{4}\kappa^4 + \lambda_3^2\kappa^{-2} - \frac{\lambda_2}{2}\kappa^2 = \lambda_4 \qquad (1.6) $$

where the constant of integration $\lambda_4$ is given by

$$ \lambda_4 \triangleq a_5^2 + a_6^2 - \frac{1}{4}(a_2^2 + a_3^2)^2 + \frac{1}{2}(a_2^2 + a_3^2)(a_1^2 - 2a_4) - a_1(a_2a_5 + a_3a_6) $$

By making the change of variable $v = \kappa^2$, (1.6) transforms to

$$ \dot{v}^2 + v^3 - 2\lambda_2v^2 - 4\lambda_4v + 4\lambda_3^2 = 0 \qquad (1.7) $$
<figure><img src="image_2.png" /><figcaption>Figure 1: Plot of the cubic polynomial $P(v)$ with respect to the squared curvature $v$. Hatched regions correspond to impossible values illustrating that the only valid range for $v$ is given by $\alpha_2$ and $\alpha_3$, the zeros of $P(v)$.</figcaption></figure>

As already stated in [4], this equation is in the form $\dot{v}^2 = P(v)$ with $P$ is the cubic polynomial

$$P(v) \triangleq -v^3 + 2\lambda_2 v^2 + 4\lambda_4 v - 4\lambda_3^2. \qquad (1.8)$$

Let $-\alpha_1, \alpha_2, \alpha_3$ be the zeros of the polynomial $P(v)$ such that

$$-\alpha_1 \le 0 \le \alpha_2 \le \alpha_3. \qquad (1.9)$$

As $P(\pm\infty) = \mp\infty$ and $P(0) = -4\lambda_3^2 \le 0$, $P(v)$ is in the form illustrated in figure 1.

Also, we have $v \ge 0$ and $P(v) \ge 0$ as they are both squares, so $v \in [\alpha_2, \alpha_3]$.

The polynomial $P(v)$ can be rewritten for its zeros by

$$P(v) = -(v + \alpha_1)(v - \alpha_2)(v - \alpha_3).$$

We can express the polynomial zeros $-\alpha_1, \alpha_2, \alpha_3$ from the constants of integrations $\lambda_i$ by

$$\begin{aligned} \alpha_1 - \alpha_2 - \alpha_3 &= -2\lambda_2 \\ \alpha_1\alpha_2 + \alpha_1\alpha_3 - \alpha_2\alpha_3 &= 4\lambda_4 \\ \alpha_1\alpha_2\alpha_3 &= 4\lambda_3^2. \end{aligned} \qquad (1.10)$$

The squared curvature $v$ can be expressed in terms of elliptic functions by

$$v(t) = \alpha_3 \left(1 - n \operatorname{sn}^2 \left(r t + \varphi | m|\right)\right) \qquad (1.11)$$

the parameter $m$, the characteristic $n$ and $r$ can be expressed from the polynomial zeros by

$$m = \frac{\alpha_3 - \alpha_2}{\alpha_3 + \alpha_1} \qquad n = \frac{\alpha_3 - \alpha_2}{\alpha_3} \qquad r = \frac{1}{2} \sqrt{\alpha_3 + \alpha_1} \qquad (1.12)$$
Given

$$ \varpi \triangleq \sqrt{\frac{1}{n}\left(1 - \frac{a_2^2 + a_3^2}{\alpha_3}\right)} \qquad (1.13) $$

the phase $\varphi$ can be retrieved from $a_2^2 + a_3^2 = \alpha_3(1 - n \operatorname{sn}^2(\varphi|m))$, and is given by

$$ \varphi = \operatorname{sgn}(a_3 a_5 - a_2 a_6) \operatorname{arcsn}(\varpi|m) \qquad (1.14) $$

where $\operatorname{arcsn}$ is the inverse of the Jacobi elliptic function sn.

Note that from (1.9), we have $0 \le m \le n \le 1$.

As outlined in [3], it has been shown the Hamiltonian vector fields in (1.4) is integrable and we have proved it can be expressed in the following form

$$ \left\{ \begin{aligned} \mu_2 &= \kappa \sin \psi \\ \mu_3 &= \kappa \cos \psi \\ \mu_4 &= \frac{1}{2}(\lambda_2 + \frac{a_1^2}{2} - \upsilon) \\ \mu_5 &= -\dot{\kappa} \cos \psi + \kappa \dot{\psi} \sin \psi \\ \mu_6 &= \dot{\kappa} \sin \psi + \kappa \dot{\psi} \cos \psi. \end{aligned} \right. \qquad (1.15) $$

where

$$ \psi(t) = \lambda_1 t - \frac{\lambda_3}{\alpha_3 r} \left( \Pi(n, \operatorname{am}(rt + \varphi|m)|m) - \Pi(n, \operatorname{am}(\varphi|m)|m) \right) + \psi(0) $$

with $\Pi(n, u|m)$ the elliptic integral of the third kind and $\operatorname{am}(u|m)$ is the Jacobi amplitude.

## 2 Planar case

Although neither the curve $q(t)$ nor the rod sensitivity $\frac{\partial q(t)}{\partial a}$ can be explicitly expressed in the general 3-D case, we will show in this section that closed forms can be obtained in the planar case which can be treated as a particular case of the previously presented model.

### 2.1 Curvature and internal wrenches

Considering only planar curves $q(t)$ in the $xy$-plane with $q = (0, 0, \theta, x, y, 0)^T$, Hamiltonian vector fields defined in (1.3) simplify to

$$ \left\{ \begin{aligned} \dot{\mu}_1 &= 0 \\ \dot{\mu}_2 &= 0 \\ \dot{\mu}_3 &= -\mu_5 \\ \dot{\mu}_4 &= \mu_3 \mu_5 \\ \dot{\mu}_5 &= -\mu_3 \mu_4 \\ \dot{\mu}_6 &= 0 \end{aligned} \right. \qquad (2.1) $$

Closed-forms of rod internal wrenches $\mu(t)$ defined in (1.15) reduce to

$$ \left\{ \begin{aligned} \mu_1 &= 0 \\ \mu_2 &= 0 \\ \mu_3 &= \kappa \\ \mu_4 &= -\frac{1}{2}(\kappa^2 + \lambda_2) \\ \mu_5 &= -\dot{\kappa} \\ \mu_6 &= 0 \end{aligned} \right. \qquad (2.2) $$
And constants of integration defined in (1) simplify to

$$
\left\{
\begin{array}{l}
\lambda_1 = 0 \\
\lambda_2 = a_3^2 + 2a_4 \\
\lambda_3 = 0 \\
\lambda_4 = a_5^2 - a_3^2 (\frac{1}{4}a_3^2 + a_4)
\end{array}
\right.
$$

We retrieve the same results as we would have obtain by applying the same problem formu-
lation on the Lie Group SE(2) rather than SE(3). Therefore, in the rest of this section we will
restrict to solutions of (1.2) that are similar to trajectories on SE(2), which are generated by the
subset of initial conditions {a ∈ A : (a₁, a₂, a₆) = (0, 0, 0)}.

In the following equations, when referring to an elliptic function pq, we will simplify the
notation pq(u|m) to pq u. Also, let us define

$$
\Gamma(t) \triangleq rt + \varphi
$$

and the following constants of motion that be needed in the following developments

$\epsilon \triangleq \operatorname{sgn}(a_3)$

$$
\Gamma_0 \triangleq \Gamma(0)
$$

$$
\delta \triangleq \lambda_2^2 + 4\lambda_4
$$

The expression of the phase $\varphi$ given in (1.14) simplifies to

$$
\varphi = \operatorname{sgn}(a_3 a_5) \operatorname{arcsn}(\varpi | m)
$$

where $\varpi$ given in (1.13) reduces to

$$
\varpi = \sqrt{\frac{1}{n} \left( 1 - \frac{a_3^2}{\alpha_3} \right)}.
$$

**2.1.1 Expression of the curvature**

In the planar case, as $\lambda_3 = 0$, the polynomial $P(v)$ simplifies to

$$
P(v) = -v^3 + 2\lambda_2 v^2 + 4\lambda_4 v
$$

so $P(v)$ has one trivial zero at $v = 0$.

From (1.9), we can distinguish three cases as outlined in [4] and [5]:

• Case I: $\lambda_4 > 0$

Using (1.10), we have that

$$
\lambda_1(\lambda_2 + \lambda_3) > \lambda_2\lambda_3.
$$

This imposes the choice for the zeros to

$$
\begin{cases}
\alpha_1 = -\lambda_2 + \sqrt{\delta} \\
\alpha_2 = 0 \\
\alpha_3 = \lambda_2 + \sqrt{\delta}.
\end{cases}
$$

From (1.12) we get $n = 1$, so the squared curvature formula in (1.11) simplifies to

$$
\begin{align*}
\upsilon(t) &= \alpha_3 (1 - \operatorname{sn}^2 \Gamma(t)) \\
            &= \alpha_3 \operatorname{cn}^2 \Gamma(t).
\end{align*}
$$
Then the signed curvature is given by

$$ \kappa(t) = \varepsilon\sqrt{\alpha_3} \operatorname{cn} \Gamma(t). \qquad (2.3) $$

The curvature $\kappa(t)$ oscillates between $\sqrt{\alpha_3}$ and $-\sqrt{\alpha_3}$ and the resulting curve $q(t)$ is called a "wavelike" elastica.

* Case II: $\lambda_4 < 0$ Using (1.10), we have

$$ \lambda_1(\lambda_2 + \lambda_3) < \lambda_2\lambda_3. $$

This imposes the choice for the zeros to

$$ \begin{cases} \alpha_1 = 0 \\ \alpha_2 = \lambda_2 - \sqrt{\delta} \\ \alpha_3 = \lambda_2 + \sqrt{\delta}. \end{cases} $$

From (1.12) we get $n = m$, so the squared curvature formula in (1.11) simplifies to

$$ \begin{aligned} v(t) &= \alpha_3 (1 - m \operatorname{sn}^2 \Gamma(t)) \\ &= \alpha_3 \operatorname{dn}^2 \Gamma(t). \end{aligned} $$

Then the signed curvature is given by

$$ \kappa(t) = \varepsilon\sqrt{\alpha_3} \operatorname{dn} \Gamma(t) \qquad (2.4) $$

The curvature $\kappa(t)$ is non-vanishing and the resulting curve $q(t)$ is called a "orbit-like" elastica.

* Case III: $\lambda_4 = 0$ This borderline case implies the polynomial $P(v)$ reduces to

$$ P(v) = -v^3 + 2\lambda_2 v^2 $$

which has a double zero.

Using first equation of (1.10), only one choice is possible for the zeros $\alpha_i$:

$$ \begin{cases} \alpha_1 = \alpha_2 = 0 \\ \alpha_3 = |2\lambda_2|. \end{cases} $$

which leads to the signed curvature

$$ \kappa(t) = \varepsilon\sqrt{\alpha_3} \operatorname{sech} \Gamma(t) \qquad (2.5) $$

This corresponds to the borderline case where the curvature is non-periodic.

### 2.1.2 Reduction to a unique formulation of the curvature

These cases can be reduced to a single formulation of the curvature by allowing the parameter $m$ to be any positive or null real and applying the Jacobi's real transformation (see [1] §16.11). By relaxing the constraint on the zeros $\alpha_i$ given in (1.9), and keeping only one fixed choice on the zeros that we will denote by $\alpha'$ such that
$$ \left\{ \begin{array}{l} \alpha'_{1} = 0 \\ \alpha'_{2} = \lambda_{2} - \sqrt{\delta} \\ \alpha'_{3} = \lambda_{2} + \sqrt{\delta} \end{array} \right. $$

In this form, $\alpha'_3$ is positive as $\sqrt{\delta} > |\lambda_2|$ but $\alpha'_2$ can now be negative. Note that we still have $\alpha'_3 \ge \alpha'_2$.

Using same forms as in (1.12), the elliptic parameter $m'$ and $r'$ by

$$ m' = \frac{\alpha'_{3} - \alpha'_{2}}{\alpha'_{3}} \qquad r' = \frac{1}{2} \sqrt{\alpha'_{3}} $$

but as mentioned before, the new elliptic parameter $m'$ is only constrained to in $[0, \infty)$. Then, the signed curvature can be expressed by a unique expression by

$$ \kappa(t) = \varepsilon \sqrt{\alpha'_{3}} \operatorname{dn} (r' (t + \phi) |m') \qquad (2.6) $$

When $m' > 1$, the Jacobi's real transformation can be applied to reduce to a parameter $m$ such that $0 \le m \le 1$ and we retrieve the previously described cases.

### 2.1.3 Explicit formulation of rod total elastic energy

Recall from (1.1) the total elastic energy of the rod is given by

$$ \begin{aligned} E_{el} &= \frac{1}{2} \int_{0}^{1} u_{3}(t)^{2} dt \\ &= \frac{1}{2} \int_{0}^{1} \kappa(t)^{2} dt \end{aligned} $$

Using the unique formulation of the curvature $\kappa(t)$ given in (2.6) can be integrated to give an explicit formulation in terms of the elliptic integral of the second kind by

$$ E_{el} = \frac{\alpha'_{3}}{2} E (r' (t + \phi) |m') $$

## 2.2 Integration of the curve $q(t)$

From the differential system defined in (1.2), it follows that

$$ \dot{\theta} = u_3 = \kappa \qquad \dot{x} = \cos \theta \qquad \dot{y} = \cos \theta. $$

Using (2.2), the integration of the curvature is given by

$$ \cos \theta(t) = \beta_{1}(0)\beta_{1}(t) + 4\beta_{2}(0)\beta_{2}(t) \qquad (2.7a) $$

$$ \sin \theta(t) = 2 \varepsilon (\beta_{1}(0)\beta_{2}(t) - \beta_{2}(0)\beta_{1}(t)) \qquad (2.7b) $$

$$ x(t) = \beta_{1}(0) \int \beta_{1}(t) + 4\beta_{2}(0) \int \beta_{2}(t) \qquad (2.7c) $$

$$ y(t) = 2\varepsilon \left( \beta_{1}(0) \int \beta_{2}(t) - \beta_{2}(0) \int \beta_{1}(t) \right) \qquad (2.7d) $$

The functions $\beta_1(t)$ and $\beta_2(t)$ can be explicitly given using Jacobi elliptic functions and the elliptic integral of second kind $E(u|m)$ in the three cases previously described as follows
• Case I: $\lambda_4 > 0$

Integrating the curvature in (2.3) (see [1] §16.24) leads to

$$ \theta(t) = 2 \varepsilon (\arccos (\mathrm{dn} \, \Gamma(t)) - \arccos (\mathrm{dn} \, \Gamma_0)) $$

Let $A(t) \triangleq \arccos (\mathrm{dn} \, \Gamma(t))$ and $A(0) \triangleq \arccos (\mathrm{dn} \, \Gamma_0)$, then

$$ \cos A(t) = \mathrm{dn} \, \Gamma(t) $$

$$ \sin A(t) = \pm \sqrt{1 - \mathrm{dn}^2 \Gamma(t)} $$

$$ = \sqrt{m} \sin \Gamma(t) $$

Given that $\frac{\theta(t)}{2} = \varepsilon (A(t) - A(0))$, we have

$$ \cos \frac{\theta(t)}{2} = \cos A(t) \cos A(0) + \sin A(t) \sin A(0) $$

$$ = \mathrm{dn} \, \Gamma(t) \mathrm{dn} \, \Gamma_0 + m \sin \Gamma(t) \sin \Gamma_0 $$

$$ \sin \frac{\theta(t)}{2} = \sin A(t) \cos A(0) - \cos A(t) \sin A(0) $$

$$ = \varepsilon \sqrt{m} (\sin \Gamma(t) \mathrm{dn} \, \Gamma_0 - \sin \Gamma_0 \mathrm{dn} \, \Gamma(t)) $$

Using half-angle formulas, we get

$$ \cos \theta(t) = \cos^2 \frac{\theta(t)}{2} - \sin^2 \frac{\theta(t)}{2} $$

$$ = (2 \mathrm{dn}^2 \Gamma_0 - 1) (2 \mathrm{dn}^2 \Gamma(t) - 1) + 4m \mathrm{dn} \, \Gamma(t) \sin \Gamma(t) \mathrm{dn} \, \Gamma_0 \sin \Gamma_0 $$

$$ \sin \theta(t) = 2 \cos \frac{\theta(t)}{2} \sin \frac{\theta(t)}{2} $$

$$ = 2 \varepsilon \sqrt{m} \left( (2 \mathrm{dn}^2 \Gamma_0 - 1) \mathrm{dn} \, \Gamma(t) \sin \Gamma(t) - (2 \mathrm{dn}^2 \Gamma(t) - 1) \mathrm{dn} \, \Gamma_0 \sin \Gamma_0 \right) $$

which is in the form (2.7) with $\beta_1$ and $\beta_2$ given by

$$ \beta_1(t) \triangleq 2 \mathrm{dn}^2 \Gamma(t) - 1 \qquad (2.8a) $$

$$ \beta_2(t) \triangleq \sqrt{m} \sin \Gamma(t) \mathrm{dn} \, \Gamma(t) \qquad (2.8b) $$

and can be integrated to

$$ \int \beta_1(t) = 2r^{-1} (E(am \Gamma(t)) - E(am \Gamma_0)) - t \qquad (2.9a) $$

$$ \int \beta_2(t) = -r^{-1} (\mathrm{cn} \, \Gamma(t) - \mathrm{cn} \, \Gamma_0) \qquad (2.9b) $$

(2.9c)

• Case II: $\lambda_4 < 0$

Integrating the curvature in (2.4) (see [1] §16.24) leads to

$$ \theta(t) = 2 \varepsilon (\arcsin (\mathrm{sn} \, \Gamma(t)) - \arcsin (\mathrm{sn} \, \Gamma_0)) $$

Let $A(t) \triangleq \arcsin (\mathrm{sn} \, \Gamma(t))$ and $A(0) \triangleq \arcsin (\mathrm{sn} \, \Gamma_0)$, then

$$ \cos A(t) = \pm \sqrt{1 - \mathrm{sn}^2 \Gamma(t)} $$

$$ = \mathrm{cn} \, \Gamma(t) $$

$$ \sin A(t) = \mathrm{sn} \, \Gamma(t) $$
Given that $\frac{\theta(t)}{2} = \varepsilon (A(t) - A(0))$, we have

$$ 
\begin{aligned}
\cos \frac{\theta(t)}{2} &= \cos A(t) \cos A(0) + \sin A(t) \sin A(0) \\
&= \operatorname{cn} \Gamma(t) \operatorname{cn} \Gamma_0 + \operatorname{sn} \Gamma(t) \operatorname{sn} \Gamma_0
\end{aligned}
 $$

$$ 
\begin{aligned}
\sin \frac{\theta(t)}{2} &= \sin A(t) \cos A(0) - \cos A(t) \sin A(0) \\
&= \varepsilon (\operatorname{sn} \Gamma(t) \operatorname{cn} \Gamma_0 - \operatorname{cn} \Gamma_0 \operatorname{sn} \Gamma(t))
\end{aligned}
 $$

Using half-angle formulas, we get

$$ 
\begin{aligned}
\cos \theta(t) &= \cos^2 \frac{\theta(t)}{2} - \sin^2 \frac{\theta(t)}{2} \\
&= (1 - 2 \operatorname{dn}^2 \Gamma(t)) (1 - 2 \operatorname{dn}^2 \Gamma_0) + 4 \operatorname{cn} \Gamma(t) \operatorname{sn} \Gamma(t) \operatorname{cn} \Gamma_0 \operatorname{sn} \Gamma_0
\end{aligned}
 $$

$$ 
\begin{aligned}
\sin \theta(t) &= 2 \cos \frac{\theta(t)}{2} \sin \frac{\theta(t)}{2} \\
&= 2\varepsilon \left( (1 - 2 \sin^2 \Gamma_0) \operatorname{cn} \Gamma(t) \sin \Gamma(t) - (1 - 2 \operatorname{dn}^2 \Gamma(t)) \operatorname{cn} \Gamma_0 \sin \Gamma_0 \right)
\end{aligned}
 $$

which is in the form (2.7) with $\beta_1$ and $\beta_2$ given by

$$ \beta_1(t) \triangleq 1 - 2 \sin^2 \Gamma(t) \qquad (2.10a) $$

$$ \beta_2(t) \triangleq \operatorname{sn} \Gamma(t) \operatorname{cn} \Gamma(t) \qquad (2.10b) $$

and can be integrated to

$$ \int \beta_1(t) = m^{-1} (t(m-2) + 2r^{-1} (E(am \Gamma(t)) - E(am \Gamma_0))) \qquad (2.11a) $$

$$ \int \beta_2(t) = -r^{-1} (\operatorname{cn} \Gamma(t) - \operatorname{cn} \Gamma_0) \qquad (2.11b) $$

(2.11c)

### Case III: $\lambda_4 = 0$

Integrating the curvature in (2.5) (see [1] §16.24) leads to

$$ \theta(t) = 2 \varepsilon (\operatorname{arctan}(\sinh \Gamma(t)) - \operatorname{arctan}(\sinh \Gamma_0)) $$

Let $A(t) \triangleq \operatorname{arctan}(\sinh \Gamma(t))$ and $A(0) \triangleq \operatorname{arctan}(\sinh \Gamma_0)$, then

$$ 
\begin{aligned}
\cos A(t) &= \pm (1 + \sinh^2 \Gamma(t))^{-\frac{1}{2}} \\
&= \operatorname{sech} \Gamma(t) \\
\sin A(t) &= \operatorname{tanh} \Gamma(t)
\end{aligned}
 $$

Given that $\frac{\theta(t)}{2} = \varepsilon (A(t) - A(0))$, we have

$$ 
\begin{aligned}
\cos \frac{\theta(t)}{2} &= \cos A(t) \cos A(0) + \sin A(t) \sin A(0) \\
&= \operatorname{sech} \Gamma(t) \operatorname{sech} \Gamma_0 + \operatorname{tanh} \Gamma(t) \operatorname{tanh} \Gamma_0
\end{aligned}
 $$

$$ 
\begin{aligned}
\sin \frac{\theta(t)}{2} &= \sin A(t) \cos A(0) - \cos A(t) \sin A(0) \\
&= \varepsilon (\operatorname{tanh} \Gamma(t) \operatorname{sech} \Gamma_0 - \operatorname{tanh} \Gamma_0 \operatorname{sech} \Gamma(t))
\end{aligned}
 $$
• The elliptic parameters $m, n$ and $r$

$$ \frac{\partial m}{\partial a} = \frac{1}{(\alpha_3 + \alpha_1)^2} \left( \left( \frac{\partial \alpha_3}{\partial a} - \frac{\partial \alpha_2}{\partial a} \right) (\alpha_3 + \alpha_1) - \left( \frac{\partial \alpha_3}{\partial a} + \frac{\partial \alpha_1}{\partial a} \right) (\alpha_3 - \alpha_2) \right) \quad (2.17a) $$

$$ \frac{\partial r}{\partial a} = \frac{1}{4\sqrt{\alpha_3 + \alpha_1}} \left( \frac{\partial \alpha_3}{\partial a} + \frac{\partial \alpha_1}{\partial a} \right) \quad (2.17b) $$

$$ \frac{\partial n}{\partial a} = \frac{1}{\alpha_3^2} \left( \frac{\partial \alpha_3}{\partial a} \alpha_2 - \frac{\partial \alpha_2}{\partial a} \alpha_3 \right) \quad (2.17c) $$

• The phase $\varphi$

Given

$$ \frac{\partial \varphi}{\partial a} = \frac{1}{2n\varphi} \left( \frac{a_3^2}{\alpha_3^2} \left( \begin{matrix} \frac{\partial \alpha_3}{\partial a_3} - 2 \frac{\alpha_3}{a_3} \\ \frac{\partial \alpha_3}{\partial a_4} \\ \frac{\partial \alpha_3}{\partial a_5} \end{matrix} \right) - \frac{1}{n} \left( 1 - \frac{a_3^2}{\alpha_3^2} \right) \frac{\partial n}{\partial a} \right) \quad (2.18) $$

and the first order derivatives of the function $\arcsin(z|m)$

$$ \frac{\partial \arcsin(z|m)}{\partial z} = \frac{1}{\sqrt{1 - z^2 \sqrt{1 - m z^2}}} $$

$$ \frac{\partial \arcsin(z|m)}{\partial m} = \frac{1}{2(m-1)m} \left( \frac{m\sqrt{1-z^2}z}{\sqrt{1-mz^2}} - E(\arcsin z|m) - (m-1)F(\arcsin z|m) \right) $$

with $cd(z|m)$ is the Jacobi elliptic function defined by

$$ cd z = \frac{cn z}{dn z}, $$

we can express the derivative of the function $\arcsin(\varphi|m)$ with respect to $a$ using the chain rule

$$ \begin{aligned} \frac{\partial \arcsin \varphi}{\partial a} &= \frac{\partial \arcsin \varphi}{\partial \varphi} \frac{\partial \varphi}{\partial a} + \frac{\partial \arcsin \varphi}{\partial m} \frac{\partial m}{\partial a} \\ &= \frac{1}{\sqrt{1 - m z^2}} \left( \frac{1}{\sqrt{1 - z^2}} \frac{\partial \varphi}{\partial a} + \frac{m z \sqrt{1 - z^2} - E(\arcsin \varphi) - (m-1)F(\arcsin \varphi)}{2(m-1)m} \frac{\partial m}{\partial a} \right) \end{aligned} $$

Then, the general expression of the derivative of the phase $\varphi$ with respect to $a$ is

$$ \frac{\partial \varphi}{\partial a} = \operatorname{sgn}(a_3 a_5) \frac{\partial \arcsin(\varphi|m)}{\partial a} $$

• The function $\Gamma(t)$

$$ \frac{\partial \Gamma(t)}{\partial a} = t \frac{\partial r}{\partial a} + \frac{\partial \varphi}{\partial a} $$

• The Jacobi elliptic function $sn(\Gamma(t)|m)$

Given the first order derivatives of the function $sn(z|m)$

$$ \frac{\partial \operatorname{sn}(z|m)}{\partial z} = \operatorname{cn}(u|m) \operatorname{dn}(u|m) $$

$$ \frac{\partial \operatorname{sn}(z|m)}{\partial m} = \frac{\operatorname{dn}(z|m) \operatorname{cn}(z|m) ((1-m)z - E(\operatorname{am}(z|m)|m) + m \operatorname{cd}(z|m) \operatorname{sn}(z|m))}{2m(1-m)} $$
we can compute directly

$$ \begin{aligned} \frac{\partial \operatorname{sn} \Gamma(t)}{\partial a} &= \frac{\partial \operatorname{sn} \Gamma(t)}{\partial \Gamma(t)} \frac{\partial \Gamma(t)}{\partial a} + \frac{\partial \operatorname{sn} \Gamma(t)}{\partial m} \frac{\partial m}{\partial a} \\ &= \operatorname{cn} \Gamma(t) \operatorname{dn} \Gamma(t) \left( \frac{\partial \Gamma(t)}{\partial a} + \frac{(m-1)\Gamma(t) + \operatorname{E}(\operatorname{am} \Gamma(t)) - m \operatorname{cd} \Gamma(t) \operatorname{sn} \Gamma(t)}{2m(m-1)} \frac{\partial m}{\partial a} \right). \end{aligned} $$

* The Jacobi elliptic function $\operatorname{cn}(\Gamma(t)|m)$

Given the first order derivatives of the function $\operatorname{cn}(z|m)$

$$ \begin{aligned} \frac{\partial \operatorname{cn}(z|m)}{\partial z} &= -\operatorname{sn}(u|m) \operatorname{dn}(u|m) \\ \frac{\partial \operatorname{cn}(z|m)}{\partial m} &= \frac{\operatorname{dn}(z|m) \operatorname{sn}(z|m) ((m-1)z + \operatorname{E}(\operatorname{am}(z|m)|m) - m \operatorname{cd}(z|m) \operatorname{sn}(z|m))}{2m(1-m)} \end{aligned} $$

we can compute directly

$$ \begin{aligned} \frac{\partial \operatorname{cn} \Gamma(t)}{\partial a} &= \frac{\partial \operatorname{cn} \Gamma(t)}{\partial \Gamma(t)} \frac{\partial \Gamma(t)}{\partial a} + \frac{\partial \operatorname{cn} \Gamma(t)}{\partial m} \frac{\partial m}{\partial a} \\ &= -\operatorname{sn} \Gamma(t) \operatorname{dn} \Gamma(t) \left( \frac{\partial \Gamma(t)}{\partial a} + \frac{(m-1)\Gamma(t) + \operatorname{E}(\operatorname{am} \Gamma(t)) - m \operatorname{cd} \Gamma(t) \operatorname{sn} \Gamma(t)}{2m(m-1)} \frac{\partial m}{\partial a} \right). \end{aligned} $$

* The Jacobi elliptic function $\operatorname{dn}(\Gamma(t)|m)$

Given the first order derivatives of the function $\operatorname{dn}(z|m)$

$$ \begin{aligned} \frac{\partial \operatorname{dn}(z|m)}{\partial z} &= -m \operatorname{cn}(u|m) \operatorname{sn}(u|m) \\ \frac{\partial \operatorname{dn}(z|m)}{\partial m} &= \frac{\operatorname{sn}(z|m) \operatorname{cn}(z|m) ((m-1)z + m \operatorname{E}(\operatorname{am}(z|m)|m) - m \operatorname{dn}(z|m) \operatorname{sc}(z|m))}{2(1-m)} \end{aligned} $$

with $\operatorname{sc}(z|m)$ is the Jacobi elliptic function defined by

$$ \operatorname{sc} z = \frac{\operatorname{sn} z}{\operatorname{cn} z}, $$

we can compute directly

$$ \begin{aligned} \frac{\partial \operatorname{dn} \Gamma(t)}{\partial a} &= \frac{\partial \operatorname{dn} \Gamma(t)}{\partial \Gamma(t)} \frac{\partial \Gamma(t)}{\partial a} + \frac{\partial \operatorname{dn} \Gamma(t)}{\partial m} \frac{\partial m}{\partial a} \\ &= -m \operatorname{cn} \Gamma(t) \operatorname{sn} \Gamma(t) \left( \frac{\partial \Gamma(t)}{\partial a} + \frac{(m-1)\Gamma(t) + m \operatorname{E}(\operatorname{am} \Gamma(t)) - m \operatorname{dn} \Gamma(t) \operatorname{sc} \Gamma(t)}{2m(m-1)} \frac{\partial m}{\partial a} \right). \end{aligned} $$

* The elliptic integral of the second kind $\operatorname{E}(\operatorname{am}(\Gamma(t)|m)|m)$

We first need to express the derivative of the amplitude $\operatorname{am}(\Gamma(t)|m)$ with respect to $a$. Given first derivatives of the Jacobi amplitude $\operatorname{am}(z|m)$

$$ \begin{aligned} \frac{\partial \operatorname{am}(z|m)}{\partial z} &= \operatorname{dn}(u|m) \\ \frac{\partial \operatorname{am}(z|m)}{\partial m} &= \frac{\operatorname{dn}(z|m) ((m-1)z + \operatorname{E}(\operatorname{am}(z|m)|m)) - m \operatorname{cn}(z|m) \operatorname{sn}(z|m)}{2m(1-m)} \end{aligned} $$
we get the derivative of am(Γ(t)|m) with respect to a by applying the chain rule as usual

$$ \begin{aligned} \frac{\partial \operatorname{am} \Gamma(t)}{\partial a} &= \frac{\partial \operatorname{am} \Gamma(t)}{\partial \Gamma(t)} \frac{\partial \Gamma(t)}{\partial a} + \frac{\partial \operatorname{am} \Gamma(t)}{\partial m} \frac{\partial m}{\partial a} \\ &= \operatorname{dn} \Gamma(t) \left( \frac{\partial \Gamma(t)}{\partial a} + \frac{(m-1)\Gamma(t) + \operatorname{E}(\operatorname{am} \Gamma(t)) - m \operatorname{cd} \Gamma(t) \operatorname{sn} \Gamma(t)}{2m(m-1)} \frac{\partial m}{\partial a} \right). \end{aligned} $$

Then, given first derivatives of the elliptic integral of the second kind E(z|m)

$$ \begin{aligned} \frac{\partial \operatorname{E}(z|m)}{\partial z} &= \sqrt{1 - m \sin^2 z} \\ \frac{\partial \operatorname{E}(z|m)}{\partial m} &= \frac{\operatorname{E}(z|m) - \operatorname{F}(z|m)}{2m} \end{aligned} $$

Noting the simplification

$$ \begin{aligned} \frac{\partial \operatorname{E}(\operatorname{am} \Gamma(t))}{\partial \operatorname{am} \Gamma(t)} &= \sqrt{1 - m \sin^2 (\operatorname{am} \Gamma(t))} \\ &= \sqrt{1 - \operatorname{sn}^2 \Gamma(t)} \\ &= \operatorname{dn} \Gamma(t), \end{aligned} $$

and that

$$ F(\operatorname{am} \Gamma(t)) = \Gamma(t), $$

we finally have all the expressions to compute the derivative of E(am Γ(t)) with respect to a by applying the chain rule

$$ \begin{aligned} \frac{\partial \operatorname{E}(\operatorname{am} \Gamma(t))}{\partial a} &= \frac{\partial \operatorname{E}(\operatorname{am} \Gamma(t))}{\partial \operatorname{am} \Gamma(t)} \frac{\partial \operatorname{am} \Gamma(t)}{\partial a} + \frac{\partial \operatorname{E}(\operatorname{am} \Gamma(t))}{\partial m} \frac{\partial m}{\partial a} \\ &= \operatorname{dn}^2 \Gamma(t) \left( \frac{\partial \Gamma(t)}{\partial a} + \frac{\operatorname{E}(\operatorname{am} \Gamma(t)) - \operatorname{cd} \Gamma(t) \operatorname{sn} \Gamma(t)}{2(m-1)} \frac{\partial m}{\partial a} \right) \end{aligned} $$

Then, most of these forms simplify in the three previously introduced cases and we can give the explicit forms of the derivatives of functions β₁(t) and β₂(t) (and their respective integrals) with respect to a:

* Case I: λ₄ > 0

From (2.8) and (2.9), we get

$$ \begin{aligned} \frac{\partial \beta_1(t)}{\partial a} &= 4 \operatorname{dn} \Gamma(t) \frac{\partial \operatorname{dn} \Gamma(t)}{\partial a} \\ \frac{\partial \beta_2(t)}{\partial a} &= \operatorname{dn} \Gamma(t) \frac{\partial \operatorname{sn} \Gamma(t)}{\partial a} + \operatorname{sn} \Gamma(t) \frac{\partial \operatorname{dn} \Gamma(t)}{\partial a} \\ \frac{\partial \int \beta_1(t)}{\partial a} &= \frac{2}{r} \left( \left( \frac{\partial \operatorname{E}(\Gamma(t))}{\partial a} - \frac{\partial \operatorname{E}(\Gamma_0)}{\partial a} \right) - \frac{1}{r} \left( \operatorname{E}(\Gamma(t)) - \operatorname{E}(\Gamma_0) \right) \frac{\partial r}{\partial a} \right) \\ \frac{\partial \int \beta_2(t)}{\partial a} &= \frac{1}{r} \left( \frac{1}{r} (\operatorname{cn} \Gamma(t) - \operatorname{cn} \Gamma_0) \frac{\partial r}{\partial a} - \left( \frac{\partial \operatorname{cn} \Gamma(t)}{\partial a} - \frac{\partial \operatorname{cn} \Gamma_0}{\partial a} \right) \right) \end{aligned} $$

* Case II: λ₄ < 0
From (2.10) and (2.11), we get

$$
\begin{align*}
\frac{\partial \beta_1(t)}{\partial a} &= -4 \operatorname{sn} \Gamma(t) \frac{\partial \operatorname{sn} \Gamma(t)}{\partial a} \\
\frac{\partial \beta_2(t)}{\partial a} &= \operatorname{cn} \Gamma(t) \frac{\partial \operatorname{sn} \Gamma(t)}{\partial a} + \operatorname{sn} \Gamma(t) \frac{\partial \operatorname{cn} \Gamma(t)}{\partial a} \\
\frac{\partial \int \beta_1(t)}{\partial a} &= \frac{1}{m} \left( \frac{t(1-m)}{m} \frac{\partial m}{\partial a} + \frac{2}{r} (\operatorname{E}(\operatorname{am} \Gamma(t)) - \operatorname{E}(\operatorname{am} \Gamma_0)) \left( \frac{1}{r} \frac{\partial r}{\partial a} - \frac{1}{m} \frac{\partial m}{\partial a} \right) \right. \\
&\qquad \left. + \frac{2}{r} \left( \frac{\partial \operatorname{E}(\Gamma(t))}{\partial a} - \frac{\partial \operatorname{E}(\Gamma_0)}{\partial a} \right) \right) \\
\frac{\partial \int \beta_2(t)}{\partial a} &= \frac{1}{rm} \left( (\operatorname{dn} \Gamma(t) - \operatorname{dn} \Gamma_0) \left( \frac{1}{m} \frac{\partial m}{\partial a} + \frac{1}{r} \frac{\partial r}{\partial a} \right) - \left( \frac{\partial \operatorname{dn} \Gamma(t)}{\partial a} - \frac{\partial \operatorname{dn} \Gamma_0}{\partial a} \right) \right)
\end{align*}
$$

• Case III: $\lambda_4 = 0$

From (2.12) and (2.13), we get

$$
\begin{align*}
\frac{\partial \beta_1(t)}{\partial a} &= -2 \operatorname{sech}^2 \Gamma(t) \tanh \Gamma(t) \\
\frac{\partial \beta_2(t)}{\partial a} &= \operatorname{sech} \Gamma(t) (1 - 2 \tanh^2 \Gamma(t)) \\
\frac{\partial \int \beta_1(t)}{\partial a} &= \frac{2}{r} \left( \frac{1}{r} \frac{\partial r}{\partial a} (\tanh \Gamma(t) - \tanh \Gamma_0) - (\tanh^2 \Gamma(t) - \tanh^2 \Gamma_0) \right) \\
\frac{\partial \int \beta_2(t)}{\partial a} &= \frac{1}{r} \left( \tanh \Gamma(t) \operatorname{sech} \Gamma(t) - \tanh \Gamma_0 \operatorname{sech} \Gamma_0 - \frac{1}{r} (\operatorname{sech} \Gamma(t) - \operatorname{sech} \Gamma_0) \right)
\end{align*}
$$

Acknowledgments

This work was supported by the French National Research Agency under the project Flecto
(ANR- Digital Models).

References

[1] M. Abramowitz and I. Stegun. *Handbook of Mathematical Functions*. Dover, 1964.

[2] T. Bretl and Z. McCarthy. Quasi-static manipulation of a kirchhoff elastic rod based on a geometric analysis of equilibrium configurations. *I. J. Robotic Res.*, 33(1):48–68, 2014.

[3] V. Jurdjevic. *Integrable Hamiltonian systems on complex Lie groups*. Memoirs of the American Mathematical Society. 2005.

[4] J. Langer and D. A. Singer. The total squared curvature of closed curves. *J. Differential Geom.*, 20(1):1–22, 1984.

[5] D. A. Singer. Lectures on elastic curves and rods, 2007.