# Large deformation analysis of spontaneous twist and contraction in nematic
elastomer fibres with helical director.

Andrea Giudici<EMAIL_ADDRESS>John S. Biggins  Department of Engineering,
University of Cambridge, Trumpington St., Cambridge CB21PZ, U.K.

###### Abstract

A cylindrical rubber fibre subject to twist will also elongate: a
manifestation of Poynting’s effect in large strain elasticity. Here, we
construct an analogous treatment for an active rubber fibre actuated via an
axisymmetric pattern of spontaneous distortion. We start by constructing an
exact large-deformation solution to the equations of elasticity for such fibre
subject to imposed twist and stretch, which reveals spontaneous warping and
twisting of the fibre cross-section absent in passive rubbers. We then compute
the corresponding non-linear elastic energy, which encompasses the Poynting
effect, but is minimized by a finite spontaneous twist and stretch. In the
second half of the paper, we apply these results to understand the twist-
contraction actuation of nematic elastomer fibres fabricated with director-
fields that encode helical patterns of contraction on heating. We first
consider patterns making a constant angle with respect to the local
cylindrical coordinate system (conical spiral director curves) and verify the
predicted spontaneous twist, contraction and cross-section deformation via
finite elements. Secondly, we consider realistic director distributions for
the experimentally reported fibres fabricated by cross-linking while
simultaneously applying stretch and twist. Counter-intuitively, we find that
maximum actuation twist is produced by applying a finite optimal twist during
fabrication. Finally, we illustrate that spontaneously twisting fibres will
coil into spring-like shapes on actuation if the ends are prevented from
twisting relative to each other. Such twist-torsion coupling would allow to
make a tendril-like “soft-spring” actuator with low force and high linear
stroke compared to the intrinsic contraction of the elastomer itself.

††preprint: AIP/123-QED

##  I Introduction

John Harrison invented the bi-metalic thermostat in 1759. Ever since,
scientists and engineers have been deploying spatial patterns of spontaneous
deformation to induce complex and dramatic actuation in solid materials. The
thermal strains in metals and shape memory alloys are limited to a few
percent, but, in recent years, the soft-matter community has demonstrated
several systems in which patterns of geometrically large strains can be
programmed into soft solids. Prominent examples include patterns of swelling
in gels Klein, Efrati, and Sharon (2007); Kim  _et al._ (2012); Na  _et al._
(2016); Gladman  _et al._ (2016), patterns of contraction in liquid crystal
elastomers de Haan  _et al._ (2012); Ware  _et al._ (2015); Aharoni  _et al._
(2018); Barnes and Verduzco (2019), and patterns of inflation in “baromorphs”
Siéfert  _et al._ (2019); Warner and Siéfert (2020). These large and
exquisitely programmable shape changes can appear to bring the matter to life
Camacho-Lopez  _et al._ (2004); White  _et al._ (2008); Gladman  _et al._
(2016), and this is no coincidence, as they strongly resemble the patterns of
muscular contraction that drive biological locomotion and the patterns of
growth that underpin biological development Thompson  _et al._ (1942).

This special issue focuses on the programming of such spontaneous shape
changes in liquid crystal elastomers (LCEs). These are rubbery networks of
rod-like mesogens which spontaneously align along a director, ${{\bm{n}}}$, to
form a nematic phase Warner and Terentjev (2007). On heating or illumination,
the nematic order can be disrupted (reflecting the nematic to isotropic
transition in conventional liquid crystals, fig. 3a) and, in LCEs, this
transition is accompanied by a dramatic and reversible contraction by a factor
of $\lambda_{s}\sim 0.5$ parallel to the alignment director ${{\bm{n}}}$
Warner and Terentjev (2007); Küpfer and Finkelmann (1991). LCEs are thus
promising artificial muscles and soft actuators De Gennes, Hébert, and Kant
(1997); Wermter and Finkelmann (2001).

Shape programming in LCEs is typically achieved by fabricating an elastomer in
which the director $\bm{n}$ is spatially varying, generating a corresponding
pattern of contraction on heating. Director programming can be implemented by
using surface-anchoring to pattern the director field in a nematic liquid
sheet, and then crosslinking to form an elastomer de Haan  _et al._ (2012);
Ware  _et al._ (2015). Alternatively, one may use an aligning stress field to
orient the director during crosslinking. This latter strategy was used to
create the original globally aligned monodomain LCEs Küpfer and Finkelmann
(1991), and, much more recently, has been deployed to generate patterns of
alignment during extrusion based 3-D printing Ambulo  _et al._ (2017);
Kotikian  _et al._ (2018); López-Valdeolivas  _et al._ (2018) and during the
direct shape programming of dual-network LCE sheets Barnes and Verduzco
(2019).

The significant majority of work on LCE shape programming has focused on 2D
sheets that bend and morph into curved surfaces on heating Warner (2020);
Aharoni  _et al._ (2018); Mostajeran  _et al._ (2016); Ware  _et al._ (2015).
However, recently, _Nocentini et al_ demonstrated an LCE fibre that was
twisted and stretched during cross-linking to imprint a helical director
fieldNocentini  _et al._ (2017). On heating, the resulting helical contraction
caused the fibre to spontaneously twist and contract. The designHaines  _et
al._ (2016); Aziz and Spinks (2020) and mechanicsCharles, Gazzola, and
Mahadevan (2019) of other torsional artificial muscles has been the focus of
much recent attention, as they offer a minaturizable version of a conventional
rotary engine. However, these previous torsional muscles are fabricated by
twisting multiple component fibres togetherForoughi  _et al._ (2011); Tondu
(2012); Zhang and Philen (2012); Haines  _et al._ (2014); Yuan  _et al._
(2019) (for example by twisting passive in-extensible fibres around an
inflating core Tondu (2012); Zhang and Philen (2012)) and their action relies
on slip within the resultant fibre bundle. In contrast, the LCE torsional
muscle is a monolithic cross-linked solid, and must be understood within the
framework of misfit elasticity.

In this manuscript, we seek to construct such an elastic theory to predict and
explain the spontaneous twist/stretch actuation of _Nocentini et al_ ’s LCE
fibres. At first sight, the natural starting point is Timishenko’s
paradigmatic calculation of the curvature of a bimetalic stripTimoshenko
(1925). Indeed, analogues of this small-strain analysis are frequently
deployed to model LCE bilayers, and other spontaneously bending sheets and
strips. Moreover, very recently, several authors have developed a
corresponding theories of spontaneous bend and twist in elastic rodsAharoni
_et al._ (2012); Kohn and O’Brien (2018); Moulton, Lessinnes, and Goriely
(2020); Bauer, Neukamm, and Schäffner (2020); Cicalese, Ruf, and Solombrino
(2017), and these certainly offer considerable insight into _Nocentini et
al’s_ LCE fibres. However, all such theories are only valid for small
spontaneous strains (or, more precisely, small incompatibilities of
spontaneous strains), leading to the complete decoupling of stretch, twist and
bend in the resultant elastic energies, and making such theories formally
inapplicable to the large strains generated in LCEs.

In contrast, in the field of rubber elasticity, there is a classic result,
discovered by PoyntingPoynting (1913) in 1913, that a rubber fibre that is
twisted substantially will also stretch in response. This large strain effect
cannot be captured by the small-strain approaches, but is captured by a simple
and exact large-deformation solution for a twisted and stretched rubber
cylinder Horgan and Murphy (2011); Zurlo  _et al._ (2020). Here, we derive an
analogue of this exact large deformation solution for a rubber cylinder
subject to an axisymetric pattern of spontaneous distortion. We find that
these spontaneous distortions introduce simple modifications to the energy, so
that it is minimized by an overall spontaneous twist and contraction. We then
compute these spontaneous twists and contractions for various LCE fibres with
different director fields, highlighting how a helical field with both
azimuthal and longitudinal components is required to produce spontaneous
twist. We validate our results with full 3-D finite element simulations.
Finally, we estimate the spontaneous twist and stretch expected in fibres
created by applying twist and stretch during cross-linking, as reported by
_Nocentini et. al._ Nocentini  _et al._ (2017). Counter-intuitively, we find
that maximum spontaneous twist is achieved by an optimal finite twist during
cross-linking, with both too little or too much twist yielding lower
performance. In conclusion, we discuss how these twisting fibres could be used
to create coiling artificial muscles, with greatly amplified stroke compared
to the intrinsic actuation of the LCE itself.

##  II Twist and stretch of a rubber cylinder

We start by recalling the classic large deformation solution for a twisted and
stretched passive rubber cylinderHorgan and Murphy (2011); Zurlo  _et al._
(2020). More precisely, we consider a long cylindrical rubber fibre with
undeformed radius $R_{0}$ and length $L\gg R_{0}$, that is subject to an
angular twist $\Delta\theta$ between the two ends and an overall extension by
a factor of $\lambda$. Working in cylindrical coordinates, if this deformation
maps the material point initially at ${{\bm{R}}}=(R,\Theta,Z)$ to the point
${{\bm{r}}}=(r,\theta,z)$, then the resultant deformation gradient is simply:

$\displaystyle{{\bm{F}}}=\frac{\partial{{\bm{r}}}}{\partial{{\bm{R}}}}=\left(\begin{array}[]{ccc}\frac{\partial
r}{\partial R}&\frac{1}{R}\frac{\partial r}{\partial\Theta}&\frac{\partial
r}{\partial Z}\\\ r\frac{\partial\theta}{\partial
R}&\frac{r}{R}\frac{\partial\theta}{\partial\Theta}&r\frac{\partial\theta}{\partial
Z}\\\ \frac{\partial z}{\partial R}&\frac{1}{R}\frac{\partial
z}{\partial\Theta}&\frac{\partial z}{\partial Z}\\\ \end{array}\right).$ (4)

In general, the hyper-elastic energy of a deformed solid may be written as

$\mathcal{E}=\int_{V}W({{\bm{F}}})dV$ (5)

where $V$ is the volume in the reference configuration and $W$ is a energy
density which depends on the deformation gradient. For an incompressible
rubber, the simplest (neo-Hookean) energy density is given by

$W({{\bm{F}}})=\frac{1}{2}\mu\,\text{Tr}\left({{\bm{F}}}.{{\bm{F}}}^{T}\right)+p\left(\det({{\bm{F}}})-1\right),$
(6)

where $\mu$ is the shear modulus and the second term arises to impose volume
conservation, $\det({{\bm{F}}})=1$, with $p$ being a spatially dependent
Lagrange multiplier describing the pressure field in the material.

Our challenge is to minimize this energy for the fibre subject to an overall
imposed twist $\tau=\Delta\theta/L$ and longitudinal stretch $\lambda$. Since
we are looking for states of uniform stretch and twist along the length of the
cylinder, the outer surface must have the form

$z(R_{0},\Theta,Z)=\lambda
Z\,\,\,\,\,\,,\,\,\,\,\,\,\,\theta(R_{0},\Theta,Z)=\Theta+\tau Z.$ (7)

These forms clarify that $\tau=\frac{d\theta}{dZ}$ corresponds to the angular
twist per unit _reference state_ length.

Figure 1: a) Schematics of the decomposition of
${{\bm{F}}}={{\bm{T}}}\cdot{{\bm{\lambda}}}$. b) Poynting effect: the
equilibrium stretch increases as more twist is imposed on the system.

Minimising the elastic energy with respect to variations in
${{\bm{r}}}({{\bm{R}}})$ leads to the traditional bulk equations of mechanical
equilibrium

$\nabla\cdot{{\bm{\Sigma}}}=0,$ (8)

and free boundary condition

${{\bm{\Sigma}}}\cdot\hat{{{\bm{R}}}}=0$ (9)

where the first Piola-Kirchhoff stress tensor is given by

${{\bm{\Sigma}}}\equiv{\frac{\partial
W}{\partial{{\bm{F}}}}}=\mu{{\bm{F}}}+p\det({{\bm{F}}}){{\bm{F}}}^{-T},$ (10)

Finally, minimising with respect to variations in $p$ returns the expected
bulk condition of impressibility

$\det({{\bm{F}}})=1.$ (11)

In the case of the twisted stretched fibre, the simplest possible fields are
those with homogeneous twist and (isochoric) stretch

$\displaystyle r(R,\Theta,Z)$ $\displaystyle=R/\sqrt{\lambda}$
$\displaystyle\theta(R,\Theta,Z)$ $\displaystyle=\Theta+\tau Z$ $\displaystyle
z(R,\Theta,Z)$ $\displaystyle=\lambda Z$ $\displaystyle p(R,\Theta,Z)$
$\displaystyle=-1/\lambda,$

which do indeed solve all the relevant bulk and boundary conditions. The
associated deformation gradient is simply

$\displaystyle{{\bm{F}}}=\left(\begin{array}[]{ccc}\lambda^{-1/2}&0&0\\\
0&\lambda^{-1/2}&R\tau\lambda^{-1/2}\\\ 0&0&\lambda\\\ \end{array}\right),$
(16)

and, upon substituting this back into the elastic energy, we obtain the total
energy of a twisted and stretched fibre as:

$\frac{\mathcal{E}}{\mu\pi
LR_{0}^{2}}=\frac{1}{2}\left(\lambda^{2}+\frac{2}{\lambda}+\frac{R_{0}^{2}\tau^{2}}{2\lambda}\right).$
(17)

We note that the first two terms in the parenthesis are simply the familiar
uni-axial stretching energy for a neo-Hookean rubber, while the final term
determines the energy cost of twisting the fibre. The energy is trivially
minimised by $\tau=0$ for all values of $\lambda$, but the minimum stretch is
given by

$\lambda_{m}=\sqrt[3]{\frac{(R_{0}\tau)^{2}+4}{4}}\geq 1,$ (18)

with $\lambda_{m}=1$ only when $\tau=0$, as illustrated in fig. 1b. This
asymmetric coupling between the two quantities is a direct manifestation of
the Poynting effect Billington (1986) in nonlinear elasticity. To clarify its
origin, we note that the deformation gradient can be achieved as a pure twist
followed by a pure stretch, as illustrated in figure 1a. Mathematically, this
corresponds to decomposing the deformation gradient as

${{\bm{F}}}={{\bm{\lambda}}}.{{\bm{T}}},$ (19)

where

${{\bm{\lambda}}}=\text{diag}(\lambda^{-1/2},\lambda^{-1/2},\lambda)\,\,\,\,\,\,\text{,}\,\,\,\,\,{{\bm{T}}}={{\bm{\delta}}}+R\tau\,\,\hat{{{\bm{e}}}}_{r}\hat{{{\bm{e}}}}_{\Theta},$
(20)

and ${{\bm{\delta}}}$ is the identity matrix. From the form of ${{\bm{T}}}$,
we see that a smaller radius implies a smaller deformation for the same twist.
Indeed, the twist energy $R_{0}\tau^{2}/\lambda$ can be written as
$R_{f}^{2}\tau^{2}$ where $R_{f}=R_{0}/\sqrt{\lambda}$ is the final radius of
the fibre. These consideration suggest that twist energy is partially relieved
by stretch as this reduces the radius via Poisson effects. The resultant
twist-stretch coupling is a paradigmatic example of the inherent geometric
non-linearity of large strains.

##  III Twist and stretch induced by a spontaneous deformation.

###  III.1 Spontaneous deformation field

In light of what we learned from the simple twisting case, we now turn our
attention to a cylindrical rubber fibre that undergoes a heterogeneous
spontaneous distortion, such that, locally, the energy minimizing deformation
is given by ${{\bm{F}}}={{\bm{G}}}(R,\Theta,Z)$. Since our ultimate motivation
is to understand the spontaneous twist and stretch of nematic LCE fibres, we
restrict consideration to ${{\bm{G}}}$ that are axisymetric, isochoric, and
independent of $Z$. However, the pattern of spontaneous distortion is allowed
to be incompatible, so that the cylinder cannot attain ${{\bm{F}}}={{\bm{G}}}$
throughout, but will instead relax to an internally stressed state that
minimizes the total elastic energy.

If the actual local deformation from the original state, prior to spontaneous
distortion is ${{\bm{F}}}$, then the elastic deformation from the local
relaxed state is simply ${{\bm{F}}}\cdot{{\bm{G}}}^{-1}$, where the second
term reverses the effect of spontaneous deformation, and the first applies the
actual deformation. The new elastic energy of the fibre, after spontaneous
distortions, is thus

$\mathcal{E}=\int_{V}\left[{\textstyle\frac{1}{2}}\mu{\text{Tr}{({{\bm{F}}}\cdot{{\bm{G}}}^{-1}\cdot{{\bm{G}}}^{-T}\cdot{{\bm{F}}}^{T})}}+p\,(\det({{\bm{F}}})-1)\right]dV.$
(21)

This "multiplicative decomposition” formDiCarlo and Quiligotti (2002) was
first introduced for elasto-plastic deformationsLee and Liu (1967) and now
pervades and unifies the study of solids with spontaneous deformations,
including growing tissues Dervaux and Amar (2008); Tallinen  _et al._ (2014),
swelling gels Dervaux and Amar (2011); Tallinen, Biggins, and Mahadevan
(2013), thermal expansion Vujošević and Lubarda (2002); Kohn and O’Brien
(2018) and, as we shall clarify later, nematic elastomers heated to the
isotropic state. Importantly, the spontaneous deformation only affects the
energy via the combination ${{\bm{g}}}={{\bm{G}}}^{-1}\cdot{{\bm{G}}}^{-T}$
(corresponding to the Finger tensor of ${{\bm{G}}}$) which, in our fibers,
takes the symmetric and axisymetric form

${{\bm{g}}}=\left(\begin{array}[]{ccc}g_{RR}(R)&g_{R\Theta}(R)&g_{RZ}(R)\\\
g_{R\Theta}(R)&g_{\Theta\Theta}(R)&g_{\Theta Z}(R)\\\ g_{RZ}(R)&g_{\Theta
Z}(R)&g_{ZZ}(R)\\\ \end{array}\right).$ (22)

The isochoric condition on all rubbers, including LCEs, requires
$\det{({{\bm{g}}})}=\det({{\bm{G}}})=1$.

###  III.2 Resultant deformation fields

Minimising (21) with respect to variations in ${{\bm{r}}}({{\bm{R}}})$ and
$p({{\bm{R}}})$ leads to the same bulk and boundary eqns as before ((8), (9),
(11)), but with the new Piola-Kirchhoff tensor given by

${{\bm{\Sigma}}}=\mu{{\bm{F}}}\cdot{{\bm{g}}}+p\det({{\bm{F}}}){{\bm{F}}}^{-T}.$
(23)

To solve for the actual deformation field, we first observe that, given
axisymmetry, $Z$ independence and incompresibility, the deformation must take
the simple form

$\displaystyle r(R,\Theta,Z)$ $\displaystyle=R/\sqrt{\lambda}$
$\displaystyle\theta(R,\Theta,Z)$
$\displaystyle=\Theta+\tau\left(Z+f_{\theta}(R)\right)$ $\displaystyle
z(R,\Theta,Z)$ $\displaystyle=\lambda(Z+f_{z}(R))$ $\displaystyle
p(R,\Theta,Z)$ $\displaystyle=-f_{p}(R)/\lambda.$ (24)

where $f_{\theta}$, $f_{z}$ and $f_{p}$ are as-yet unknown functions of $R$,
which modify the solution from the original one for a passive rubber cylinder.
The $f_{\theta}$ term introduces rotation as a function of radius, and allows
radii in the reference configuration to become curves in the final
configuration. Similarly, the $f_{z}$ term encodes warping of each cross
section into an identical surface of revolution, while $f_{p}$ allows the
pressure to vary with radius.

At this stage, one could methodically substitute these fields into the bulk
and boundary equations, and then solve for $f_{\theta}$, $f_{z}$ and $f_{p}$.
However, one can substantially simplify this process by noting that the
resultant deformation field can now be decomposed as

$\displaystyle{{\bm{F}}}$
$\displaystyle=\left(\begin{array}[]{ccc}\frac{1}{\sqrt{\lambda}}&0&0\\\
\frac{Rf_{\theta}^{\prime}}{\sqrt{\lambda}}&\frac{1}{\sqrt{\lambda}}&\frac{R\tau}{\sqrt{\lambda}}\\\
\lambda f_{z}^{\prime}&0&\lambda\\\
\end{array}\right)={{\bm{\lambda}}}\cdot{{\bm{T}}}\cdot{{\bm{\Phi}}},$ (28)

where ${{\bm{\lambda}}}$ and ${{\bm{T}}}$ are again pure stretch and twist
(eq. (20)) while

${{\bm{\Phi}}}=\left(\begin{array}[]{ccc}1&0&0\\\
R\tau\left(f_{\theta}^{\prime}-f_{z}^{\prime}\right)&1&0\\\
f_{z}^{\prime}&0&1\\\ \end{array}\right)$ (29)

describes the deformation of the cross-section in the absence of twist or
stretch, and is the only part that depends on $f_{\theta}$ and $f_{z}$. This
decomposition is illustrated in Fig. 2

Figure 2: a) Schematics of the decomposition of
${{\bm{F}}}={{\bm{T}}}\cdot{{\bm{\lambda}}}\cdot{{\bm{\Phi}}}$.

Furthermore, since only derivatives of the two fields $f_{\theta}$ and $f_{z}$
appear in ${{\bm{F}}}$, we can substitute this deformation into the energy,
and minimise directly with respect to variations in $f^{\prime}_{z}(R)$ and
$f^{\prime}_{\theta}(R)$, leading to the conditions:

$\displaystyle{\frac{\partial W}{\partial f^{\prime}_{\theta}}}$
$\displaystyle=\mu\mathrm{Tr}\left({{\bm{\lambda}}}\cdot{{\bm{T}}}\cdot{{\bm{\Phi}}}\cdot{{\bm{g}}}\cdot{\frac{\partial{{\bm{\Phi}}}^{T}}{\partial
f^{\prime}_{\theta}}}\cdot{{\bm{T}}}^{T}\cdot{{\bm{\lambda}}}^{T}\right)=0$
(30) $\displaystyle{\frac{\partial W}{\partial f^{\prime}_{z}}}$
$\displaystyle=\mu\mathrm{Tr}\left({{\bm{\lambda}}}\cdot{{\bm{T}}}\cdot{{\bm{\Phi}}}\cdot{{\bm{g}}}\cdot{\frac{\partial{{\bm{\Phi}}}^{T}}{\partial
f^{\prime}_{z}}}\cdot{{\bm{T}}}^{T}\cdot{{\bm{\lambda}}}^{T}\right)=0.$ (31)

Evaluating these yields the simple uncoupled differential equations

$\displaystyle g_{R\Theta}+R\tau g_{RR}f^{\prime}_{\theta}+R\tau g_{RZ}$
$\displaystyle=0,$ (32)
$\displaystyle\lambda^{2}\left(g_{RR}f^{\prime}_{z}+g_{RZ}\right)$
$\displaystyle=0,$ (33)

which can be directly integrated to obtain the fields

$\displaystyle f_{\theta}$
$\displaystyle=-\int_{R_{0}}^{R}\frac{g_{RZ}(u)}{g_{RR}(u)}du-\tau^{-1}\int_{R_{0}}^{R}\frac{1}{u}\left(\frac{g_{R\Theta}(u)}{\,g_{RR}(u)}\right)du,$
$\displaystyle f_{z}$
$\displaystyle=-\int_{R_{0}}^{R}\frac{g_{RZ}(u)}{g_{RR}(u)}du.$

These solutions then imply the full form of the deformation fields:

$\displaystyle r$ $\displaystyle=R/\sqrt{\lambda}$ $\displaystyle\theta$
$\displaystyle=\tau\,z/\lambda-\int_{R_{0}}^{R}\frac{1}{u}\left(\frac{g_{r\Theta}(u)}{\,g_{RR}(u)}\right)du,$
$\displaystyle z$
$\displaystyle=\lambda\left(Z-\int_{R_{0}}^{R}\frac{g_{RZ}(u)}{g_{RR}(u)}\right)du.$
(34)

Finally, if we now substitute these back into the original bulk and boundary
equations ((8) and (9)) we can confirm they are fully solved (e.g. in
mathematica) provided the pressure is taken as

$\displaystyle f_{p}(R)$
$\displaystyle=\tau\int_{{R0}}^{R}f_{\theta}^{\prime}(u)\left(\tau\,uf_{\theta}^{\prime}(u)g_{{RR}}(u)+2\tau
ug_{{RZ}}(u)+2g_{R\Theta}(u)\right)\,du$
$\displaystyle\,\,\,\,\,-g_{{RR}}(R)+\int_{{R0}}^{R}\left(g_{\Theta\Theta}(u)-g_{RR}(u)\right)u^{-1}\,du$
$\displaystyle\,\,\,\,\,+\int_{{R0}}^{R}\tau\left(2g_{{\Theta Z}}(u)+\tau
ug_{{ZZ}}(u)\right)\,du.$ (35)

###  III.3 Twist and stretch elastic energy

We now substitute the fields into the elastic energy (eqn. 21) to evaluate the
energy for a rubber rod with imposed twist $\tau$, stretch $\lambda$ and a
pattern of spontaneous deformation ${{\bm{g}}}$. To do this, it is convenient
to use the decomposition in (28). We then see immediately that the effects of
${{\bm{g}}}$ in the energy are entirely contained within the symmetric tensor

$\displaystyle\tilde{{{\bm{g}}}}={{\bm{\Phi}}}\cdot{{\bm{g}}}\cdot{{\bm{\Phi}}}^{T}=\left(\begin{array}[]{ccc}g_{RR}&0&0\\\
0&g_{\Theta\Theta}-\frac{g_{R\Theta}^{2}}{g_{{RR}}}&g_{{\theta
Z}}-\frac{g_{{R\Theta}}g_{{RZ}}}{g_{{RR}}}\\\ 0&g_{{\Theta
Z}}-\frac{g_{{R\Theta}}g_{{RZ}}}{g_{{RR}}}&g_{{ZZ}}-\frac{g_{{RZ}}^{2}}{g_{{RR}}}\\\
\end{array}\right),$ (39)

which captures the residual part of ${{\bm{g}}}$ after allowing the cross-
sections to relax via ${{\bm{\Phi}}}$, but in the absence of twist or stretch.
By inspection we may then multiply out all the terms in the energy to find:

$\displaystyle\frac{\mathcal{E}}{\mu\pi LR_{0}^{2}}$
$\displaystyle=\frac{a_{0}}{\lambda}+a_{1}\lambda^{2}+\frac{b}{\lambda}\tau+\frac{c}{\lambda}\tau^{2}$
(40)

with the coefficients given by:

$\displaystyle a_{0}$
$\displaystyle={\textstyle\frac{1}{2}}\int_{0}^{R_{0}}\left(\tilde{{{\bm{g}}}}_{RR}+\tilde{{{\bm{g}}}}_{\Theta\Theta}\right)RdR$
$\displaystyle a_{1}$
$\displaystyle={\textstyle\frac{1}{2}}\int_{0}^{R_{0}}\tilde{{{\bm{g}}}}_{ZZ}RdR$
$\displaystyle b$ $\displaystyle=\int_{0}^{R_{0}}\tilde{{{\bm{g}}}}_{\Theta
Z}R^{2}dR$ $\displaystyle c$
$\displaystyle={\textstyle\frac{1}{2}}\int_{0}^{R_{0}}\tilde{{{\bm{g}}}}_{ZZ}R^{3}dR.$
(41)

We note that these terms are related to the zeroth, first and second moments
of aspects of the spontaneous distortion, reminiscent of those derived via
Gamma convergence for linear elastic rods Kohn and O’Brien (2018).

As in the simple case explored in section II, twist and stretch are coupled
via the $\lambda^{-1}$ in the twist energy. However, with the addition of the
linear term in $\tau$, when $b\neq 0$ the energy is minimised by a non-zero
twist $\tau_{m}$ and a finite stretch $\lambda_{m}$ given by:

$\tau_{m}=-\frac{b}{2c}\,,\,\,\,\,\,\,\,\,\,\,\lambda_{m}=\sqrt[3]{\frac{4a_{0}c-b^{2}}{2a_{1}c}}.$
(42)

We can use these and $\bar{\lambda}=\lambda/\lambda_{m}$ to rewrite an energy
with the same structure as in eqn (17), but now including the source terms,

$\displaystyle\frac{\mathcal{E}}{\mu\pi LR_{0}^{2}}$
$\displaystyle=\lambda_{m}^{2}\left[a_{1}\left(\frac{2}{\bar{\lambda}}+\bar{\lambda}^{2}\right)+\frac{(c/\lambda_{m}^{3})}{\bar{\lambda}}(\tau-\tau_{m})^{2}\right].$
(43)

This is our main result and captures the emergent twist/stretch nature of
fibres with a spontaneous deformation field.

Thus, our calculation reveals that an axisymmetric spontaneous deformation can
induce 4 key effects in a fibre:
1\. spontaneous length change.
2\. spontaneous twist.
3\. cross-section warping, via the function $f_{z}$, so initial disks become
surfaces of revolution.
4\. cross-section twisting via $f_{\theta}$, so initial radii become curves.
These effects stem from different aspects of the spontaneous deformation, and
may or may not co-occur in a given case.

In the following sections, we concentrate on spontaneous deformation fields
arising from a LCE undergoing nematic-isotropic transition.

##  IV Nematic elastomer fibres

Figure 3: a) Schematics of isotropic-nematic transition and the resulting
alignment along ${{\bm{n}}}$. b) Example of a $R$-independent director field
with $\alpha=\pi/4$ and $\beta=\pi/2$, as shown on the top face. In red, we
mark a helical integral curve to highlight the chirality of the director
field. Figure 4:  a) 2D plot of twist as a function of both angles $\alpha$
and $\beta$. The twist is maximised when $\beta=\pi/2$ and $\alpha\sim\pi/4$.
b) – c) Variation of the twist and stretch as a function of $\alpha$ for
different values of $\lambda_{s}$.

We want to apply our theoretical machinery to nematic elastomer cylinders,
encoded with a spatially variable director pattern ${{\bm{n}}}$. On heating
through the nematic-isotropic transition schematically shown in fig 3 a), the
elastomer will spontaneously contract by a factor $\lambda_{s}$ along the
director and, to preserve volume, elongate by $1/\lambda_{s}$ in the two
perpendicular directions, corresponding to a spontaneous deformation
${{\bm{G}}}=\text{diag}(\lambda_{s}^{-1/2},\lambda_{s}^{-1/2},\lambda_{s})$ in
a frame aligned with the director Warner and Terentjev (2007). We can thus
write ${{\bm{g}}}$ for such spontaneous actuation as:

${{\bm{g}}}=\lambda_{s}^{-2}{{\bm{n}}}{{\bm{n}}}+\lambda_{s}({{\bm{\delta}}}-{{\bm{n}}}{{\bm{n}}}).$
(44)

As discussed in Appendix A, the resulting pre-strained neo-Hookean exactly
reproduces the “trace-formula” energy commonly encountered in the LCE
literatureWarner and Terentjev (2007) for an LCE cross-linked in the nematic
state then heated to the isotropic state. This form neglects stress induced
changes in the degree of nematic alignment, which is appropriate in the
isotropic phase, where an LCE is simply a classical rubber provided one is not
too close to the transition temperature. _Nocentini et al_ ’s experimental LCE
undergoes its transition in a 10K window around 373K, and is heated to 393K
for actuation, allowing for this approximation. We also note that the nematic
actuation strains during the transition, $\delta L/L\sim 0.3$, vastly exceed
other thermal effects such as (linear) thermal expansion which has
characteristic size $\delta L/L\sim 10^{-5}K^{-1}$.

We again limit attention to axisymmetric director patterns. In cylindrical
coordinates, the (unit) director field is then described by two angles
$\alpha(R)$ and $\beta(R)$ (fig 3 b) so that

${{\bm{n}}}=\cos\beta\sin\alpha{{\bm{\hat{e}}}}_{R}+\sin\beta\sin\alpha{{\bm{\hat{e}}}}_{\Theta}+\cos\alpha{{\bm{\hat{e}}}}_{Z}.$
(45)

In what follows, we first consider the case of $\alpha$ and $\beta$
independent of $R$, and show that spontaneous twist requires an oblique
director field between azimuthal and longitudinal directions, while warping of
the cross-sections requires a director that is oblique between the
longitudinal and radial directions. Finally, we consider some more realistic
$R$ dependent director fields for the fibres created by twisting and
stretching during cross-linkingNocentini  _et al._ (2017), and show that
spontaneous twisting is maximized by a finite optimum degree of twisting at
genesis.

###  IV.1 Example: R-independent Director field.

####  IV.1.1 Theoretical predictions

We start by considering director patterns in which the angles $\alpha$ and
$\beta$ are constants, independent of $R$. In this case, the integrals in our
solutions can be conducted analytically and yield the displacement fields

$\displaystyle r$ $\displaystyle=R/\sqrt{\lambda}$ $\displaystyle\theta$
$\displaystyle=\Theta+\tau\left(Z+\frac{g_{RZ}}{g_{RR}}(R_{0}-R)\right)-\log(R/R_{0})\frac{g_{R\Theta}}{g_{RR}}$
$\displaystyle z$
$\displaystyle=\lambda\left(Z+\frac{g_{RZ}}{g_{RR}}(R_{0}-R)\right).$

From the linear dependence on $R$ in $z$ we note that flat disc cross sections
turn into cones, while $\theta(R)$ shows that radii turn into spirals which
tend to conical-spirals at the centre of the fibre. Though such singular
structures may seem surprising at first glance, a more careful look reveals
that integral curves of the director-field ${{\bm{n}}}$ are also conical
spirals making a constant angle $\beta$ with the radial direction.
Furthermore, in planar LCEs encoded with constant angle +1 defects, the
integral curves also form planar log spirals (the planar projection of a
conical spiral) and it is well established that the resultant actuation
transforms the sheet into twisted conical surfaces in which the radii
transform into conical spiralsModes, Bhattacharya, and Warner (2011);
Mostajeran  _et al._ (2016); Warner and Mostajeran (2018).

In the simple case in which $\beta=\pi/2$, the integral curves become simple
helices, as shown in Fig. 3 b), and are helpful to understand how twist is
developed during activation. In the most simple sketch of the mechanics,
heat/illumination drives a contraction by $\lambda_{s}$ along the integral
curve as well as an increase in its radius by a factor of
$\sqrt{\lambda_{s}}$. The twisting thus occurs to reduce the length of the
integral curve, with the integral curve playing the same role of a sub-fibre
in twisted fibre bundle torsional musclesForoughi  _et al._ (2011); Tondu
(2012); Zhang and Philen (2012); Haines  _et al._ (2014); Yuan  _et al._
(2019).

Given constant $\alpha$ and $\beta$, the coefficients in our energy (eqn.
(43)) are:

$\displaystyle a_{0}=$
$\displaystyle{\textstyle\frac{1}{8}}\lambda_{s}^{-2}R_{0}^{2}d2\sin^{4}\alpha\sin^{2}2\beta\left(1\\!-\\!\lambda_{s}^{3}\right)$
$\displaystyle\,\,\,\,\,\,+{\textstyle\frac{1}{8}}\lambda_{s}^{-2}R_{0}^{2}\left(\\!3\lambda_{s}^{3}\\!+\\!1\\!+\cos
2\alpha(\lambda_{s}^{3}-1)\right)$ $\displaystyle a_{1}=$
$\displaystyle{\textstyle\frac{1}{4}}R_{0}^{2}\lambda_{s}d\left(3\\!+\\!\lambda_{s}\\!-\\!\left(\cos
2\alpha\\!+\\!2\cos
2\beta\sin^{2}\alpha\right)(\lambda_{s}^{3}\\!-\\!1)\right)$ $\displaystyle
b=$ $\displaystyle{\textstyle\frac{2}{3}}R_{0}^{3}d\sin
2\alpha\sin\beta\lambda_{s}\left(1-\lambda_{s}^{3}\right)$ $\displaystyle c=$
$\displaystyle{\textstyle\frac{1}{8}}R_{0}^{4}\lambda_{s}d\left(3\\!+\lambda_{s}\\!+\\!\left(\cos
2\alpha\\!+\\!2\cos
2\beta\sin^{2}\alpha\right)\left(1-\lambda_{s}^{3}\right)\right)$ (46)

where

$\displaystyle d=\left(2\sin^{2}\alpha\cos
2\beta\left(1-\lambda_{s}^{3}\right)+\cos
2\alpha\left(\lambda_{s}^{3}-1\right)+3\lambda_{s}^{3}+1\right)^{-1}.$

These can be used together with eqns. (44) and (42), to obtain the spontaneous
twist:

$\tau_{m}=\frac{-8\sin 2\alpha\,\sin\beta}{R_{0}\left(6\sin^{2}\alpha\cos
2\beta\\!+\\!3\cos
2\alpha\\!-\\!3\left(\lambda_{s}^{3}\\!+\\!3\right)/\left(\lambda_{s}^{3}\\!-\\!1\right)\right)}.$
(47)

The behaviour of the twist as a function of $\alpha$ and $\beta$ is shown in
Figure 4 a) and b) for typical values during a nematic-isotropic (heating)
transition ($1>\lambda_{s}>0.5$). Importantly, we see that the twist vanishes
when $\beta=0$ or $\alpha=0,\,\pi/2$, corresponding to a purely longitudinal
or azimuthathal-radial director-field where no shear $\Theta$-$Z$ is present.
On the other hand, $\tau_{m}$ is maximised when $\beta=\pi/2$ (independent of
$\alpha$ and $\lambda_{s}$), indicating a director with no radial component.
The largest possible twist is then given by choosing
$\alpha=\frac{1}{2}\cos^{-1}\left(\frac{\lambda_{s}^{3}-1}{\lambda_{s}^{3}+1}\right)\sim\pi/4+{\textstyle\frac{3}{4}}(\lambda_{s}-1)+O((\lambda_{s}-1)^{2})$,
showing that maximum twist is achieved at an oblique angle biased towards the
azimuthal direction for larger $\lambda_{s}$, as shown in Fig 4 b).

Similarly, the overall spontaneous stretch, $\lambda_{m}$, is given by eqn
(42). Although the full expression is too complicated to reproduce here, in
the simple and twist-maximising case of $\beta=\pi/2$, it reduces to

$\lambda_{m}\\!=\\!\frac{36\\!\cos\\!2\alpha\left(\\!1\\!-\\!\lambda_{s}^{3}\\!\right)\lambda_{s}^{3}\\!-\\!\cos\\!4\alpha\left(\\!\lambda_{s}^{3}\\!-\\!1\\!\right)^{2}\\!+\\!37\lambda_{s}^{6}\\!+\\!106\lambda_{s}^{3}\\!+\\!\\!1}{36\left(\cos
2\alpha\left(1-\lambda_{s}^{3}\right)+\lambda_{s}^{3}+1\right)^{2}}.$

In Fig. 4 c), we can see how this stretch behaves between the two extremes of
a longitudinally aligned ($\alpha=0$) director field, yielding a simple
contraction of $\lambda_{m}=\lambda_{s}$, and of a azimuthal director field
yielding
$\lambda_{m}=\left({\textstyle\frac{1}{2}}\left(1+\lambda_{s}\right)/\lambda_{s}^{3}\right)^{1/3}$.

Figure 5:  Comparison between theoretical predictions and finite element
simulation for twist (a) and stretch (b) as a function of the spontaneous
deformation parameter $\lambda_{s}$. For these calcuations, we used a fibre
with aspect ratio $R_{0}/L=1/15$ represented by $5400$ hex-8 elements.

####  IV.1.2 Finite element verification

To test our results, we use the open-source finite element software FEBio Maas
_et al._ (2012, 2016, 2018) to compute the spontaneous deformations of an LCE
fibre encoded with an $R$-independent director field. We used FEBio’s standard
prestrained neo-Hookean material on cylindrical fibres with an almost
incompressible Poisson ratio of $\nu=0.45$. The isochoric prestrain parameter
$1/\lambda_{s}$ was set to vary between values of $0.5$ to $2$. After applying
the pre-strain, the energy minimising deformation was found using a static
analysis.

Figure 6:  Comparison between Finite element simulations and theoretical data
for the case $\alpha=\pi/4$ and $\beta=\pi/4$. a) An example of a deformed
cross section directly from FE analysis. In red are theoretical lines for the
deformed radii while dots are node positions. b) Comparison between predicted
and theoretical coning. c) Comparison between theory and FE on the predicted
rotation angle difference as a function of radius. A fibre with aspect ratio
$R_{0}/L=1/5$, and $54000$ hex-8 elements was used for this simulation.

In figure 5, we compare the predicted spontaneous twist and stretch of a fibre
encoded with $\beta=\pi/2$ and $\alpha=\pi/4$ with finite element simulations.
The theory shows excellent agreement with the numerical simulations,
accurately capturing the non-linearities in both twist and stretch.

Fibres with $\beta=\pi/2$ have cross-sections that remain flat during
deformation, as coning is driven by the $g_{RZ}$ component of the spontaneous
deformation, which is only present if the director has an $RZ$ component.
Therefore, to confirm our predictions about cross-section warping, we also
computed the deformation of a fibre with $\beta=\pi/4$ and $\alpha=\pi/4$, as
shown in Fig. 6 a). Again, comparing our theoretical results with FE
simulations we obtain excellent agreement between the numerical and
theoretical coning (Fig. 6 (b)) and winding (Fig. 6 (c)) of the cross-section.
The logarithmic nature of the spirals implies a theoretically infinite number
of rotations at the centre of the cross section, although the stress, strain
and energy are all finite. The divergent rotation stems from the line of
director discontinuity (disinclination) along the central axis of the fibre.
Accordingly, in a real fibre rotation would be cut-off near the axis by a
regularisation of the director discontinuity within a defect core De Gennes
and Prost (1993), and in our finite elements it is cut-off by the element size
near the axis. However, the director discontinuity line and associated
infinite rotation are an artefact of patterns of constant $\alpha$ and
$\beta$, and, as discussed in the next section, are not expected in the
experimentally generated fibres.

###  IV.2 Nematic fibres produced by stretching and twisting during cross-
linking

Finally, we consider the twisting LCE fibres reported by _Nocentini et al_
Nocentini  _et al._ (2017). These fibres were produced by pulling a filament
out of a viscous LC monomer mixture while rotating the drawing end, and
simultaneously cross-linking with a UV light. The director alignment is
imprinted through the strains induced during this drawing process, shown in
Figure 7 a). Given there is both twisting and stretching during crosslinking,
and twisting strains are larger at larger radii, we expect this fabrication to
produce an $R$-dependent director field with azimuthal and longitudinal
components.

In reality, the imprinting of the director-field is a complex visco-elastic
process involving sticky polymers being cross-linked into a rubber. However,
here, we take a simple approximation, and assume the deformation is mainly
elastic, and the director aligns with the direction of maximum strain. The
elastic approximation is clearly appropriate once there is sufficient
crosslinking, but is probably also applicable to the initial visco-elastic
drawing as the strain rate is rather high.

To find the imprinted director pattern, we model the fibre during drawing as
an elastic cylinder that is stretched by a factor of $\lambda_{0}$ and twisted
by $\Delta\theta_{0}$ as cross-linking proceeds, resulting in a fibre of
length $L$ and radius $R_{0}$ and (final state) twist density
$\tau_{0}=\Delta\theta_{0}/L$. We use $(\tilde{R},\tilde{\Theta},\tilde{Z})$
as the reference state coordinates for this problem, so that we may use
$(R,\Theta,Z)$ for the final state coordinates, which then become the
reference state coordinates in our spontaneous deformation analysis. The
elastic deformation follows the simple treatment in section II, leading, in
our coordinate system, to the deformation fields

$\displaystyle R(\tilde{R},\tilde{\Theta},\tilde{Z})$
$\displaystyle=\tilde{R}/\sqrt{\lambda_{0}}$
$\displaystyle\Theta(\tilde{R},\tilde{\Theta},\tilde{Z})$
$\displaystyle=\tilde{\Theta}+\tau_{0}\lambda_{0}\tilde{Z}$ $\displaystyle
Z(\tilde{R},\tilde{\Theta},\tilde{Z})$ $\displaystyle=\lambda_{0}\tilde{Z}$

Figure 7:  a) Schematic of a fibre being drawn and twisted from a drop of LC
monomer while being cured with UV light. b) An example of the resulting
director field in a fibre of radius $r=2/\tau_{0}$. The two integral curves
highlight the change in azimuthal component as a function of the radius.

To obtain the direction of maximum strain in the final (post cross-linking)
configuration, we use the left-Cauchy deformation tensor

${{\bm{b}}}={{\bm{F}}}\cdot{{\bm{F}}}^{T}=\left(\begin{array}[]{ccc}\frac{1}{\lambda_{0}}&0&0\\\
0&\frac{1}{\lambda_{0}}+\lambda_{0}^{2}R^{2}\tau_{0}^{2}&\lambda_{0}^{2}R\tau_{0}\\\
0&\lambda_{0}^{2}R\tau_{0}&\lambda_{0}^{2}\\\ \end{array}\right).$ (49)

The largest eigenvalue of ${{\bm{b}}}$ identifies the largest component of
stretch while its corresponding eigenvector (which is a target state object)
is its direction. Since this is the direction along which the director will
orient, it can be used to express the values of the angles $\alpha$ and
$\beta$ in the fibre. We trivially obtain that $\beta=\pi/2$, since the
twisting during the manufacturing of the fibre induces no coupling of the
$R$-$Z$ components. For $\alpha$, one obtains

$\tan\alpha\\!=\\!\frac{\left(\bar{\tau}_{0}^{2}-1\right)\\!+\\!\sqrt{\left(\bar{\tau}_{0}^{2}+1\right){}^{2}+2\lambda_{0}^{-3}\left(\bar{\tau}_{0}^{2}-1\right)+\lambda_{0}^{-3}}+\lambda_{0}^{-3}}{2\bar{\tau}_{0}}$
(50)

where $\bar{\tau}_{0}=R\tau_{0}$. We note that, when the fibre is drawn from a
drop, $\lambda_{0}\gg 1$, which simplifies the angle to:

$\alpha=\tan^{-1}\left(R\,\tau_{0}\right).$ (51)

This implies that the director points along the $Z$ direction in the centre of
the fibre and tilts in the $\Theta$-$Z$ plane as one moves outwards. This is a
reflection of the fact that, during formation, $\Theta$-$Z$ shears grow like
$R$ as the filament is drawn and twisted, thus inducing a greater azimuthal
component further from the centre as shown in figure 7 b).

###  IV.3 Comparison between twist during cross-linking and twist during
activation

Finally, we can use the form of $\alpha$ and $\beta$ to obtain ${{\bm{g}}}$.
We then use equations (41) and (42) to find the twist and stretch capability
of a fibre given the twist imposed at its genesis. The results for values of
$\lambda_{s}<1$ are shown in figure 8. Remarkably, the output twist does not
monotonically grow as a function of $\tau_{0}$, but reaches a maximum and
decays to zero thereafter. Recall that, in the $R$-independent field, we
discussed how the twist is maximised when $\alpha\sim\pi/4$. In this case, a
small $\tau_{0}R_{0}$ implies the director is on average mainly longitudinally
aligned ($<\alpha>\sim 0$), inducing mostly a contraction by a factor
$\lambda_{s}$. On the other hand, a large $\tau_{0}R_{0}$ leads to a
dominantly azimuthal director orientation on the cross-section
($<\alpha>\sim\pi/2$), inducing mainly stretch by a factor of
$\left({\textstyle\frac{1}{2}}(\lambda_{s}^{3}+1)/\lambda_{s}^{3}\right)^{1/3}$.
The optimum twist output $\tau$ is maximised in between the two, when the
coupling between azimuthal and longitudinal component is greatest.

Figure 8:  Relationship between output twist a) and stretch b) as a function
of the twist imposed at genesis. Shaded in grey, the region for known
experimental results Nocentini  _et al._ (2017).

The fibres produced by _Nocentini et al_ Nocentini  _et al._ (2017) were made
from LCEs capable of a maximum spontaneous contraction of $\lambda_{s}\sim
0.71$ during heating. The twist imposed at genesis on the fibres was of about
$10$ turns with their diameter and length varying between $50$ and
$300\mathrm{\mu m}$ and $1$ to $5\mathrm{cm}$ respectively. This suggests
their fibres were fabricated with a (dimensionless) genesis twist of around
$0.16<R_{0}\tau_{0}<1.8$, and the corresponding experimentally-explored region
is shaded on Fig. 8. The authors only reported the output twist and
contraction for one fibre (of unknown length and radius), when activated in a
$\sim 0.5\mathrm{cm}$ long region via light. This fibre generated an overall
contraction $\sim 0.84$, suggesting that in the activated region
$0.71<\lambda_{m}<0.84$. On the other hand, activation induced a rotation of
about $\Delta\theta\sim 460^{o}\sim 8.2\,\mathrm{rad}$ corresponding to an
output twist of $0.04<R_{0}\tau<0.5$. Both these twist and stretch values fall
in the shaded region of Fig. 8, consistent with our theoretical results. We
highlight how, in general, these experimental fibres appear to have been
generated with too little twist at genesis, yielding a sub-optimal output
twist. This could perhaps be improved by increasing the number of turns during
fabrication.

##  V conclusion and discussion

In conclusion, we have obtained an exact elastic solution for a twisted and
stretched cylindrical neo-Hookean fibre subject to an axisymmetric isochoric
spontaneous deformation field. The solution yields a full non-linear elastic
energy for such a fibre, which is minimized by a spontaneous twist and
stretch. The energy also highlights a large-deformation coupling between twist
and stretch, as familiar from the classical Poynting effect. Finally, the
elastic fields also capture the large deformations of the fibre’s cross-
section, which is predicted to warp into a surface of revolution and twist
such that radii become curves during activation.

When applied to LCEs, our results show that a helical director-field oblique
in the longitudinal-azimuthal plane is required to induce twist. The twist
output depends on the spontaneous elongation $\lambda_{s}$ and is maximised
when no radial director component is present ($\beta=\pi/2$) as well as when
the azimuthal and longitudinal components are coupled through an angle
$\alpha\sim\pi/4+{\textstyle\frac{3}{4}}(\lambda_{s}-1)$.

It is instructive to compare our theory with recent work on the spontaneous
bending and twisting of rods via incompatible (aka misfit) spontaneous
distortionsAharoni  _et al._ (2012); Cicalese, Ruf, and Solombrino (2017);
Bauer, Neukamm, and Schäffner (2020); Kohn and O’Brien (2018); Moulton,
Lessinnes, and Goriely (2020). These treatments go beyond ours in that the
spontaneous strain is not assumed to be axisymetric, and the rod is allowed to
bend so that the center line no-longer remains straight. However, these
treatments do assume small spontaneous strains and high aspect-ratio rods,
allowing a linear elastic treatment similar to the original Kirchhoff model.
The small-strain thin-rod regime limits the theories to stretch free
deformations, and leads to simple bend twist energies of the form

$\frac{\mathcal{E}}{\mu\pi
LR_{0}^{2}}=R_{0}^{2}\left({\textstyle\frac{3}{8}}({{\bm{\kappa}}}-{{\bm{\kappa}}}_{m})^{2}+{\textstyle\frac{1}{4}}(\tau-\tau_{m})^{2}\right)$
(52)

where ${{\bm{\kappa}}}$ and $\tau$ are the bend (curvature vector) and twist
of the rod, while ${{\bm{\kappa}}}_{m}$ and $\tau_{m}$ are their minimising
values. These minimising values were first estimated by linearising the
spontaneous deformation in a Taylor series about the rod’s central axis
Aharoni  _et al._ (2012). More recently, Gamma convergenceCicalese, Ruf, and
Solombrino (2017); Bauer, Neukamm, and Schäffner (2020); Kohn and O’Brien
(2018) and 3D energy minimizationMoulton, Lessinnes, and Goriely (2020) have
been used to derive rigorous forms, yielding averages and moments of various
terms of the spontaneous deformation over the rod’s cross-section. Our large-
deformation axisymmetric treatment reduces to the twisting portion of these
rod-theories in the limit of small spontaneous distortions and little imposed
stretch. Indeed, if we Taylor expand our expressions for $\lambda_{m}$ and
$\tau_{m}$ (eqn. (42)) in the limit of small spontaneous deformations
${{\bm{g}}}={{\bm{\delta}}}+\epsilon{{\bm{\delta g}}}(R)$, we find that

$\displaystyle\lambda_{m}$
$\displaystyle=1+\epsilon\frac{1}{R_{0}^{2}}\int_{0}^{R_{0}}\left(\delta\\!g_{RR}+\delta\\!g_{\Theta\Theta}\right)R\,dR$
(53) $\displaystyle\tau_{m}$
$\displaystyle=-\epsilon\frac{4}{R_{0}^{4}}\int_{0}^{R_{0}}\delta\\!g_{\Theta
Z}\,R^{2}dR.$ (54)

If we then also assume the applied twist and stretch are small,
$\lambda-\lambda_{m}\sim\epsilon$ , $\tau\sim\epsilon$, we may expand our
energy (Eqn. (43)) to $\epsilon^{2}$ to obtain

$\frac{\mathcal{E}}{\mu\pi
LR_{0}^{2}}={\textstyle\frac{3}{2}}\left(\lambda-\lambda_{m}\right)^{2}+{\textstyle\frac{1}{4}}R_{0}^{2}\left(\tau-\tau_{m}\right)^{2}+\mathcal{O}(\epsilon^{3}).$
(55)

The twisting portion of this energy agrees with that in eqn. (52), and the
linearized form of $\tau_{m}$ above agrees with that given by _Kohn and
O’Brien_ Kohn and O’Brien (2018). Interestingly, although such small-strain
and high-aspect ratio assumptions appear necessary to resolve bending, our
treatment demonstrates they can be avoided entirely when only treating
twisting and stretching. The resulting non-linear formulation generates the
highly non-linear form $\tau_{m}(\lambda_{s})$ and $\lambda_{m}(\lambda_{s})$,
as seen in fig. 5, and is clearly required for accurate predictions in large
strain systems such as LCEs.

Remarkably, even a small strain rod theory can describe large bend and twist
displacements in a suitably long rod. This introduces fascinating and rich
geometric coupling between twist and bend deformations Van Der Heijden and
Thompson (2000); Charles, Gazzola, and Mahadevan (2019); Fuller (1971);
Purohit (2008). For example, if one straightens a wound headphone wire, it
becomes highly twisted. Similarly, if one twists a fibre then brings the ends
together, it will spontaneously untwist into a lower energy spring-like coiled
state. In general, torsional bend and twist can be exchanged in a rod or
fibre, without rotating the ends, provided the total number of turns is
conserved. This twist-torsion coupling is an example of a geometric phase and
is key in the winding of DNA molecules Marko (1998); Olsen and Bohr (2011);
Fuller (1978), and the coiling of plant tendrilsGerbode  _et al._ (2012).

Therefore, although our treatment includes no mention of bend, we can infer
from this coupling effect that a spontaneously twisting LCE fibre will bend
into a coiled spring-like configuration if it activates under boundary
conditions that prevent the ends from twisting relative to each other. Given
the large-strain nature of LCE actuation, a formal treatment of this torsional
effect appears to be a formidable challenge. However, as shown in Fig. 9, we
were easily able to observe it in finite elements simulations of a twisting
LCE fibre actuated under a constant longitudinal force, but with the
constraint that the ends may not rotate. This coiling mechanics is commonly
deployed in other twisting artificial muscles to generate linear actuation
Aziz and Spinks (2020); Charles, Gazzola, and Mahadevan (2019). In contrast,
LCEs are intrinsically contractile actuators, and simple linear contraction
can be trivially achieved with a monodomain strip Küpfer and Finkelmann
(1991); De Gennes, Hébert, and Kant (1997); Wermter and Finkelmann (2001).
However, the coiling mechanism would allow an LCE actuator with much higher
stroke and lower stiffness, forming switchable soft-springs that, like plant
tendrilsGerbode  _et al._ (2012), could be used to gently anchor and position
an object in 3-D space.

Figure 9:  Finite element simulation of a fibre with Young modulus $E=1$,
Poisson ratio $\nu=0.45$, 4800 elements, aspect ratio of $R_{0}/L=0.07$ and a
$R$ independent director-field with $\alpha=\pi/4$ and $\beta=\pi/2$. The
fibre is fixed at one end and subject to a tension $T=0.2\pi\mu R_{0}^{2}$
while preventing rotation at the other end. As the value of $\lambda_{s}$ is
decreased, the fibre wants to shrink and twist but cannot do the latter as the
ends are not allowed to rotate. The fibre thus coils to release some of the
twist energy, resulting in a greater stroke amplitude then that of classical
linear actuation.

###### Acknowledgements.

A.G. thanks J.L. for the support and the EPSRC for funding, project 2108804.
J.S.B. is supported by a UKRI Future Leaders Fellowship MR/S017186/1.

DATA AVAILABILITY The data that supports the findings of this study are
available within the article.

##  Appendix A Modelling a nematic elastomer as a neo-Hookean with spontaneous
distortion

The free energy density of a nematic elastomer Warner and Terentjev (2007) has
its roots in statistical physics, and contains two contributions: an elastomer
energy from the polymer network, and a nematic energy from the rods,

$W=W_{pol}+W_{rod}.$

The rod energy depends on the scalar order parameter $Q$ of the nematic field
(i.e. the degree of alignment) but not its direction (unit vector
${{\bm{n}}}$). Appropriate forms for $W_{rod}(Q)$ are provided by the Landau-
de-Gennes (phenomenological) or Maier Saupe (microscopic) theories of liquid
nematics De Gennes and Prost (1993). Either way, $W_{rod}(Q)$ has
characteristic size of $k_{B}T$ per rod, and transitions from having a minimum
at $Q=0$ (isotropic) to a finite $Q$ (nematic) below a critical temperature
$T^{*}$.

The polymer free energy is dominated by conformational entropy, like in the
statistical theory of conventional rubber. However, in the presence of a
nematic field, the polymer random walks are not isotropic, but biased along
the director ${{\bm{n}}}$ by an amount $r(Q)$ determined by the degree of
alignment, as encoded in a a step-length tensor
${{\bm{\ell}}}\propto\delta+(r-1){{\bm{n}}}{{\bm{n}}}$. The resultant polymer
energy is described by the “trace formula” Bladon, Terentjev, and Warner
(1993); Warner and Terentjev (2007)

$W_{pol}={\textstyle\frac{1}{2}}n_{s}k_{B}T\,{\text{Tr}{({{\bm{\ell}}}_{0}\cdot{{\bm{F}}}^{T}\cdot{{\bm{\ell}}}^{-1}\cdot{{\bm{F}}})}}$

where $n_{s}$ is the density of polymer strands, ${{\bm{F}}}$ is the
deformation from the cross-linking state to the final state, $\ell_{0}$ is the
step length tensor at cross-linking (which depends on the nematic variables at
cross-linking, $Q_{0}$ ${{\bm{n_{0}}}}$) and ${{\bm{\ell}}}$ is the step-
length tensor in the final state (which depends on the final state nematic
variables at $Q$, ${{\bm{n}}}$). The full behaviour of the nematic elastomer
is now given by minimizing the sum of both energies over elastic deformations
(${{\bm{F}}}$), final state order parameter ($Q$), and final state director
(${{\bm{n}}}$). In general, this minimisation gives a two way coupling between
the nematic order and the LCE deformation. However, the characteristic size of
$W_{rod}$ is $k_{B}T$ per rod, while the characteristic size of $W_{pol}$ is
$k_{B}T$ per polymer strand. Since in an elastomer there are typically more
than ten rods per strand, the nematic energy dominates the elastic one during
minimization over $Q$. Therefore, the degree of nematic alignment is
essentially that which minimizes the nematic energy alone, and is only
modestly affected by the polymer network. For example, the shift of the
nematic-isotropic transition caused by the presence of the network is
typically only a few Kelvin Zubarev, Talroze, and Plate (1997) (compared to a
transition temperature of 350K), and other mechanically induced changes in $Q$
are similarly small Warner and Terentjev (2007); Kaufhold, Finkelmann, and
Brand (1991). Such effects are more pronounced in some modern LCE compositions
which have more crosslinks per rod Jampani  _et al._ (2019); Mistry  _et al._
(2018), but, for simplicity, we focus on the traditional case here.

We may thus consider the simpler problem of minimizing the polymer energy,
with the magnitude of nematic order, $Q$, effectively fixed as a constraint by
the rods energy. In general, when we minimize the polymer energy in the
nematic state, we must still do it over final state director ${{\bm{n}}}$;
indeed the nematic director can be observed to rotate within an elastomer in
response to stretch Warner and Terentjev (2007); Finkelmann  _et al._ (1997).
However, if we further assume (as here) that the final state is isotropic,
$Q=0$, then $\ell\propto{{\bm{\delta}}}$ must also be isotropic, and there is
no final state director to minimize over. In this case, the elastomer energy
is simply

$W={\textstyle\frac{1}{2}}n_{s}k_{B}T\,{\text{Tr}{({{\bm{F}}}\cdot{{\bm{\ell}}}_{0}\cdot{{\bm{F}}}^{T})}}$

which corresponds to the standard neo-Hookean energy with spontaneous
distortion ${{\bm{G}}}=\ell_{0}^{1/2}$ encoded by the nematic field in the
cross-linking state and shear modulus $\mu=n_{s}k_{B}T$.

## References

  * Klein, Efrati, and Sharon (2007) Y. Klein, E. Efrati, and E. Sharon, “Shaping of elastic sheets by prescription of non-euclidean metrics,” Science 315, 1116–1120 (2007).
  * Kim  _et al._ (2012) J. Kim, J. A. Hanna, M. Byun, C. D. Santangelo, and R. C. Hayward, “Designing responsive buckled surfaces by halftone gel lithography,” Science 335, 1201–1205 (2012).
  * Na  _et al._ (2016) J.-H. Na, N. P. Bende, J. Bae, C. D. Santangelo, and R. C. Hayward, “Grayscale gel lithography for programmed buckling of non-euclidean hydrogel plates,” Soft Matter 12, 4985–4990 (2016).
  * Gladman  _et al._ (2016) A. S. Gladman, E. A. Matsumoto, R. G. Nuzzo, L. Mahadevan, and J. A. Lewis, “Biomimetic 4d printing,” Nature materials 15, 413–418 (2016).
  * de Haan  _et al._ (2012) L. T. de Haan, C. Sánchez-Somolinos, C. M. Bastiaansen, A. P. Schenning, and D. J. Broer, “Engineering of complex order and the macroscopic deformation of liquid crystal polymer networks,” Angewandte Chemie International Edition 51, 12469–12472 (2012).
  * Ware  _et al._ (2015) T. H. Ware, M. E. McConney, J. J. Wie, V. P. Tondiglia, and T. J. White, “Voxelated liquid crystal elastomers,” Science 347, 982–984 (2015).
  * Aharoni  _et al._ (2018) H. Aharoni, Y. Xia, X. Zhang, R. D. Kamien, and S. Yang, “Universal inverse design of surfaces with thin nematic elastomer sheets,” Proceedings of the National Academy of Sciences 115, 7206–7211 (2018).
  * Barnes and Verduzco (2019) M. Barnes and R. Verduzco, “Direct shape programming of liquid crystal elastomers,” Soft matter 15, 870–879 (2019).
  * Siéfert  _et al._ (2019) E. Siéfert, E. Reyssat, J. Bico, and B. Roman, “Bio-inspired pneumatic shape-morphing elastomers,” Nature materials 18, 24–28 (2019).
  * Warner and Siéfert (2020) M. Warner and E. Siéfert, “Inflationary routes to gaussian curved topography,” Proc. R. Soc. A 476, 20200047 (2020).
  * Camacho-Lopez  _et al._ (2004) M. Camacho-Lopez, H. Finkelmann, P. Palffy-Muhoray, and M. Shelley, “Fast liquid-crystal elastomer swims into the dark,” Nature materials 3, 307–310 (2004).
  * White  _et al._ (2008) T. J. White, N. V. Tabiryan, S. V. Serak, U. A. Hrozhyk, V. P. Tondiglia, H. Koerner, R. A. Vaia, and T. J. Bunning, “A high frequency photodriven polymer oscillator,” Soft Matter 4, 1796–1798 (2008).
  * Thompson  _et al._ (1942) D. W. Thompson _et al._ , “On growth and form.” On growth and form. (1942).
  * Warner and Terentjev (2007) M. Warner and E. M. Terentjev,  _Liquid crystal elastomers_ , Vol. 120 (Oxford university press, 2007).
  * Küpfer and Finkelmann (1991) J. Küpfer and H. Finkelmann, “Nematic liquid single crystal elastomers,” Die Makromolekulare Chemie, Rapid Communications 12, 717–726 (1991).
  * De Gennes, Hébert, and Kant (1997) P.-G. De Gennes, M. Hébert, and R. Kant, “Artificial muscles based on nematic gels,” in  _Macromolecular Symposia_ , Vol. 113 (Wiley Online Library, 1997) pp. 39–49.
  * Wermter and Finkelmann (2001) H. Wermter and H. Finkelmann, “Liquid crystalline elastomers as artificial muscles,” e-Polymers 1 (2001).
  * Ambulo  _et al._ (2017) C. P. Ambulo, J. J. Burroughs, J. M. Boothby, H. Kim, M. R. Shankar, and T. H. Ware, “Four-dimensional printing of liquid crystal elastomers,” ACS applied materials & interfaces 9, 37332–37339 (2017).
  * Kotikian  _et al._ (2018) A. Kotikian, R. L. Truby, J. W. Boley, T. J. White, and J. A. Lewis, “3d printing of liquid crystal elastomeric actuators with spatially programed nematic order,” Advanced Materials 30, 1706164 (2018).
  * López-Valdeolivas  _et al._ (2018) M. López-Valdeolivas, D. Liu, D. J. Broer, and C. Sánchez-Somolinos, “4d printed actuators with soft-robotic functions,” Macromolecular rapid communications 39, 1700710 (2018).
  * Warner (2020) M. Warner, “Topographic mechanics and applications of liquid crystalline solids,” Annual Review of Condensed Matter Physics 11, 125–145 (2020).
  * Mostajeran  _et al._ (2016) C. Mostajeran, M. Warner, T. H. Ware, and T. J. White, “Encoding gaussian curvature in glassy and elastomeric liquid crystal solids,” Proceedings of the Royal Society A: Mathematical, Physical and Engineering Sciences 472, 20160112 (2016).
  * Nocentini  _et al._ (2017) S. Nocentini, D. Martella, D. S. Wiersma, and C. Parmeggiani, “Beam steering by liquid crystal elastomer fibres,” Soft Matter 13, 8590–8596 (2017).
  * Haines  _et al._ (2016) C. S. Haines, N. Li, G. M. Spinks, A. E. Aliev, J. Di, and R. H. Baughman, “New twist on artificial muscles,” Proceedings of the National Academy of Sciences 113, 11709–11716 (2016).
  * Aziz and Spinks (2020) S. Aziz and G. M. Spinks, “Torsional artificial muscles,” Materials Horizons 7, 667–693 (2020).
  * Charles, Gazzola, and Mahadevan (2019) N. Charles, M. Gazzola, and L. Mahadevan, “Topology, Geometry, and Mechanics of Strongly Stretched and Twisted Filaments: Solenoids, Plectonemes, and Artificial Muscle Fibers,” Physical Review Letters 123 (2019), 10.1103/PhysRevLett.123.208003.
  * Foroughi  _et al._ (2011) J. Foroughi, G. M. Spinks, G. G. Wallace, J. Oh, M. E. Kozlov, S. Fang, T. Mirfakhrai, J. D. Madden, M. K. Shin, S. J. Kim, _et al._ , “Torsional carbon nanotube artificial muscles,” Science 334, 494–497 (2011).
  * Tondu (2012) B. Tondu, “Modelling of the McKibben artificial muscle: A review,” Journal of Intelligent Material Systems and Structures 23, 225–253 (2012).
  * Zhang and Philen (2012) Z. Zhang and M. Philen, “Pressurized artificial muscles,” Journal of Intelligent Material Systems and Structures 23, 255–268 (2012).
  * Haines  _et al._ (2014) C. S. Haines, M. D. Lima, N. Li, G. M. Spinks, J. Foroughi, J. D. Madden, S. H. Kim, S. Fang, M. J. De Andrade, F. Göktepe, Ö. Göktepe, S. M. Mirvakili, S. Naficy, X. Lepró, J. Oh, M. E. Kozlov, S. J. Kim, X. Xu, B. J. Swedlove, G. G. Wallace, and R. H. Baughman, “Artificial muscles from fishing line and sewing thread,” Science 343, 868–872 (2014).
  * Yuan  _et al._ (2019) J. Yuan, W. Neri, C. Zakri, P. Merzeau, K. Kratz, A. Lendlein, and P. Poulin, “Shape memory nanocomposite fibers for untethered high-energy microengines,” Science 365, 155–158 (2019).
  * Timoshenko (1925) S. Timoshenko, “Analysis of bi-metal thermostats,” Josa 11, 233–255 (1925).
  * Aharoni  _et al._ (2012) H. Aharoni, Y. Abraham, R. Elbaum, E. Sharon, and R. Kupferman, “Emergence of spontaneous twist and curvature in non-Euclidean rods: Application to Erodium plant cells,” Physical Review Letters 108, 1–5 (2012).
  * Kohn and O’Brien (2018) R. V. Kohn and E. O’Brien, “On the Bending and Twisting of Rods with Misfit,” Journal of Elasticity 130, 115–143 (2018).
  * Moulton, Lessinnes, and Goriely (2020) D. E. Moulton, T. Lessinnes, and A. Goriely, “Morphoelastic rods III: Differential growth and curvature generation in elastic filaments,” Journal of the Mechanics and Physics of Solids 142 (2020), 10.1016/j.jmps.2020.104022.
  * Bauer, Neukamm, and Schäffner (2020) R. Bauer, S. Neukamm, and M. Schäffner, “Derivation of a Homogenized Bending–Torsion Theory for Rods with Micro-Heterogeneous Prestrain,” Journal of Elasticity 141, 109–145 (2020), arXiv:arXiv:1903.08290v2 .
  * Cicalese, Ruf, and Solombrino (2017) M. Cicalese, M. Ruf, and F. Solombrino, “On global and local minimizers of prestrained thin elastic rods,” Calculus of Variations and Partial Differential Equations 56, 1–39 (2017), arXiv:1606.04524 .
  * Poynting (1913) J. Poynting, “The changes in length and volume of an indian-rubber cord when twisted,” India-Rubber Journal (1913).
  * Horgan and Murphy (2011) C. O. Horgan and J. G. Murphy, “Extension and torsion of incompressible non-linearly elastic solid circular cylinders,” Mathematics and Mechanics of Solids 16, 482–491 (2011).
  * Zurlo  _et al._ (2020) G. Zurlo, J. Blackwell, N. Colgan, and M. Destrade, “The poynting effect,” arXiv preprint arXiv:2004.09653 (2020).
  * Billington (1986) E. W. Billington, “The Poynting effect,” Acta Mechanica 58, 19–31 (1986), arXiv:2004.09653 .
  * DiCarlo and Quiligotti (2002) A. DiCarlo and S. Quiligotti, “Growth and balance,” Mechanics Research Communications 29, 449–456 (2002).
  * Lee and Liu (1967) E. Lee and D. Liu, “Finite-strain elastic—plastic theory with application to plane-wave analysis,” Journal of applied physics 38, 19–27 (1967).
  * Dervaux and Amar (2008) J. Dervaux and M. B. Amar, “Morphogenesis of growing soft tissues,” Physical review letters 101, 068101 (2008).
  * Tallinen  _et al._ (2014) T. Tallinen, J. Y. Chung, J. S. Biggins, and L. Mahadevan, “Gyrification from constrained cortical expansion,” Proceedings of the National Academy of Sciences 111, 12667–12672 (2014).
  * Dervaux and Amar (2011) J. Dervaux and M. B. Amar, “Buckling condensation in constrained growth,” Journal of the Mechanics and Physics of Solids 59, 538–560 (2011).
  * Tallinen, Biggins, and Mahadevan (2013) T. Tallinen, J. S. Biggins, and L. Mahadevan, “Surface sulci in squeezed soft solids,” Physical review letters 110, 024302 (2013).
  * Vujošević and Lubarda (2002) L. Vujošević and V. Lubarda, “Finite-strain thermoelasticity based on multiplicative decomposition of deformation gradient,” Theoretical and applied mechanics , 379–399 (2002).
  * Modes, Bhattacharya, and Warner (2011) C. D. Modes, K. Bhattacharya, and M. Warner, “Gaussian curvature from flat elastica sheets,” Proceedings of the Royal Society A: Mathematical, Physical and Engineering Sciences 467, 1121–1140 (2011), https://royalsocietypublishing.org/doi/pdf/10.1098/rspa.2010.0352 .
  * Warner and Mostajeran (2018) M. Warner and C. Mostajeran, “Nematic director fields and topographies of solid shells of revolution,” Proceedings of the Royal Society A: Mathematical, Physical and Engineering Sciences 474, 20170566 (2018).
  * Maas  _et al._ (2012) S. A. Maas, B. J. Ellis, G. A. Ateshian, and J. A. Weiss, “Febio: finite elements for biomechanics,” Journal of biomechanical engineering 134 (2012).
  * Maas  _et al._ (2016) S. A. Maas, A. Erdemir, J. P. Halloran, and J. A. Weiss, “A general framework for application of prestrain to computational models of biological materials,” journal of the mechanical behavior of biomedical materials 61, 499–510 (2016).
  * Maas  _et al._ (2018) S. A. Maas, S. A. LaBelle, G. A. Ateshian, and J. A. Weiss, “A plugin framework for extending the simulation capabilities of febio,” Biophysical journal 115, 1630–1637 (2018).
  * De Gennes and Prost (1993) P.-G. De Gennes and J. Prost,  _The physics of liquid crystals_ , Vol. 83 (Oxford university press, 1993).
  * Van Der Heijden and Thompson (2000) G. H. Van Der Heijden and J. M. Thompson, “Helical and localized buckling in twisted rods: A unified analysis of the symmetric case,” Nonlinear Dynamics 21, 71–99 (2000).
  * Fuller (1971) F. B. Fuller, “The writhing number of a space curve.” Proceedings of the National Academy of Sciences of the United States of America 68, 815–819 (1971).
  * Purohit (2008) P. K. Purohit, “Plectoneme formation in twisted fluctuating rods,” Journal of the Mechanics and Physics of Solids 56, 1715–1729 (2008).
  * Marko (1998) J. F. Marko, “DNA under high tension: Overstretching, undertwisting, and relaxation dynamics,” Physical Review E - Statistical Physics, Plasmas, Fluids, and Related Interdisciplinary Topics 57, 2134–2149 (1998).
  * Olsen and Bohr (2011) K. Olsen and J. Bohr, “The geometrical origin of the strain-twist coupling in double helices,” AIP Advances 1 (2011), 10.1063/1.3560851.
  * Fuller (1978) F. B. Fuller, “Decomposition of the linking number of a closed ribbon: A problem from molecular biology,” Proceedings of the National Academy of Sciences of the United States of America 75, 3557–3561 (1978).
  * Gerbode  _et al._ (2012) S. J. Gerbode, J. R. Puzey, A. G. McCormick, and L. Mahadevan, “How the cucumber tendril coils and overwinds,” Science 337, 1087–1091 (2012).
  * Bladon, Terentjev, and Warner (1993) P. Bladon, E. Terentjev, and M. Warner, “Transitions and instabilities in liquid crystal elastomers,” Physical Review E 47, R3838 (1993).
  * Zubarev, Talroze, and Plate (1997) E. Zubarev, R. Talroze, and N. Plate, “Phase behavior of nematic liquid-crystalline polymer networks,” Polymer science. Series A 39, 699–704 (1997).
  * Kaufhold, Finkelmann, and Brand (1991) W. Kaufhold, H. Finkelmann, and H. R. Brand, “Nematic elastomers, 1. effect of the spacer length on the mechanical coupling between network anisotropy and nematic order,” Die Makromolekulare Chemie: Macromolecular Chemistry and Physics 192, 2555–2579 (1991).
  * Jampani  _et al._ (2019) V. Jampani, R. Volpe, K. R. de Sousa, J. F. Machado, C. Yakacki, and J. Lagerwall, “Liquid crystal elastomer shell actuators with negative order parameter,” Science advances 5, eaaw2476 (2019).
  * Mistry  _et al._ (2018) D. Mistry, S. D. Connell, S. Mickthwaite, P. B. Morgan, J. H. Clamp, and H. F. Gleeson, “Coincident molecular auxeticity and negative order parameter in a liquid crystal elastomer,” Nature communications 9, 1–9 (2018).
  * Finkelmann  _et al._ (1997) H. Finkelmann, I. Kundler, E. Terentjev, and M. Warner, “Critical stripe-domain instability of nematic elastomers,” Journal de Physique II 7, 1059–1069 (1997).