Buckling of chiral rods due to coupled axial and rotational growth

Abstract

We present a growth model for special Cosserat rods that allows for induced rotation of cross-sections. The growth law considers two controls, one for lengthwise growth and another for rotations. This is explored in greater detail for straight rods with helical and hemitropic material symmetries by introduction of symmetry-preserving growth to account for the microstructure. The example of a guided-guided rod possessing a chiral microstructure is considered to study its deformation due to growth. We show the occurrence of growth-induced out-of-plane buckling in such rods.

1. Introduction

Several theoretical models for elastic rods have been proposed. From Euler’s elastica to Kirchhoff rods, a very rich literature is available, including the general model developed by Green, Naghdi and their collaborators. The general rod theory proposed by Green and Naghdi subsumes classical theories like the Cosserat rod theory as special cases under appropriate constraints. A comprehensive description of different rod theories is provided by Antman [1] and O’Reilly [2].

Rod theories have been employed in many interesting applications in the last few decades, such as in DNA biophysics [3], marine cables [4], tendril perversion in plants [5, 6], surgical filaments [7], slender viscous jets [8], hair curls [9] and carbon nanotubes [10, 11].

Growing filamentary structures are ubiquitous in nature. Plant organs such as tendrils, roots and stem tend to twist while growing axially [12]. There are studies with helical growth models in which straight axial growth is accompanied by rotation of cross-sections [13, 14]. In this paper, we focus on this type of twisting growth, which can lead to non-planar configurations if the material of the rod exhibits some sort of twist–extension coupling.

The standard multiplicative decomposition [15] used to model biological growth has been specialized for one-dimensional growth by Moulton et al. [16]. A recent study by Moulton et al. [17] gives the reduction of three-dimensional energy for a tubular structure to a one-dimensional equivalent via minimization in cross-sections and subsequent averaging; it further demonstrates the generation of intrinsic twists and curvatures due to differential growth. A diverse account on biological growth is available from Goriely [18], containing both mathematical and biomechanical aspects.

Euler buckling of filaments evolving their shape under time-varying loads has been considered by Goldstein and Goriely [19]. Works like that of McMillen et al. [6] consider plant tendrils as Kirchhoff rods, straight in their initial states, which subsequently develop intrinsic curvatures in the grown equilibrium states. Another evolution law for intrinsic curvatures has been proposed by O’Reilly and Tresierras [20] with a focus on tip growth. Guillon et al. [21] modelled tree growth by considering the branch to be a special Cosserat rod growing in both length and diameter. They modelled the reference, relaxed and current configuration of the growing rod with separate base curves and director fields.

Several growing filaments in nature are known to have non-planar configurations [12, 13, 22]. Most existing works that study one-dimensional growth in filaments model them as isotropic rods. Such models usually rely on differential growth [17] or the presence of an external elastomeric matrix [23], multi-rod composites [24] or phototropism [25] to model the generation of curvature and torsion in non-planar configurations. Growth in chiral rods can be another way to obtain such non-planar deformations; this has not been explored in the literature. In this work, we show that growing chiral rods can buckle out of plane, simply with a boundary condition that arrests relative axial rotation at the ends.

There have been several attempts to understand material symmetries in rods. Different types of chiral material symmetries – such as hemitropy and helical symmetry – in initially straight rods with uniform circular cross-section have been investigated in great detail by Healey [26]. Other treatments of material symmetry in the context of rods include the works of Luo and O’Reilly [27] and Lauderdale and O’Reilly [28]. The latter authors [29] draw a few parallel comparisons with some results by Healey [26]. In this paper, we follow the definitions and ideas of material symmetries introduced by Healey [26, 30].

Energy representations for helical symmetry and hemitropy have been derived by Healey [26]. Multi-fold helical symmetry is useful in modelling rods whose micro-structure mimics the symmetries of a rope made up of entwined helices. Hemitropic rods possess the centre-line rotational symmetry of an isotropic rod but lack the reflection symmetries with respect to the longitudinal planes. Energy functions for rods with such chiral symmetries are typically characterized by coupled stretch, twist, shear and curvature terms. These couplings physically manifest as different types of non-traditional Poisson effects [31]. Moreover, the conventional quadratic energy densities associated with linear elasticity are incapable of distinguishing between different orders of helical symmetries and hemitropy.

Out-of-plane deformations are yet another feature of rods with such symmetries. Unshearable hemitropic rods can give rise to out-of-plane buckling when subjected to end displacements with fixed–fixed boundary conditions, but on the other hand an axial load applied to a fixed–free rod always results in a planar solution [32]. Similar bifurcation analysis has also been replicated for chiral rings with circular cross-sections under central loading [33]. Both in-plane and out-of-plane buckling of isotropic rods embedded in elastomeric matrix have been examined [23], revealing that non-planar configurations are obtained whenever the matrix is stiff enough, compared to the bending stiffness of the rod. Primary root growth of certain plants has been investigated [22], drawing analogies from mechanical buckling of a metal filament embedded in a matrix comprising two different gels whose interface is transverse to the filament.

The main focus of this work is to study growth-induced deformation in rods possessing chiral material symmetries – transverse hemitropy and dihedral helical symmetry. The growth law is also assumed to be chiral. Straight growth, where cross-sections do not rotate as they translate lengthwise, is not appropriate for modelling growth in rods with helical symmetry. If the chirality in material symmetry stems from some helical substructure associated with the rod’s microstructure or from some sort of helical fibre-reinforcement, then simple translational growth can alter the pitch of the helix. This, in turn, may modify the chiral constitutive quantities associated with the material law. Moreover, straight translational growth without any rotation can lead to unwinding or over-winding, thus inducing additional stresses; stress-free growth in such cases requires the consideration of a coupled axial and rotational growth. Modelling the virtual configuration obtained from stress-free growth as a special Cosserat rod allows us to consider growth-induced rotation of cross-sections. The exact relationship between the growth law and rod’s microstructure is not well established. In this work, we assume growth and constitutive laws to be independent in general. Additionally, for rods with helical symmetry we postulate the growth law to be symmetry preserving, so that any imaginary helix associated with the microstructure remains unaltered as the rod grows. Such a growth problem depends only on the microstructural pitch and the constitutive laws, keeping aside the boundary conditions and other external factors.

A rod constrained to grow (or decay) in a guided–guided environment is considered, with a chiral constitutive law that is applicable to both helical symmetry and transverse hemitropy. Out-of-plane buckling is observed to occur at certain growth (or atrophy) stages, corresponding to the bifurcation modes. We demonstrate that an exact reversal in chirality of these non-planar solutions requires us to mirror the chiral parameters in both growth and constitutive laws simultaneously. Comparisons are made for the end-to-end distance in the buckled configuration with that in the virtual state to see if the ends have come closer or moved apart, than what they would have been in the absence of the guides. We also show that total growth-induced extension in a rod does not depend monotonically on the degree of chirality – that is, total extension in an isotropic rod need not lie between the total extension of rods with opposite material chirality.

This paper is organized as follows. We begin with a theoretical background of material symmetries in the context of special Cosserat rods in Section 2. A twisting growth law with two control parameters is systematically derived using certain kinematic assumptions such as homogeneity in lengthwise growth and relative rotation of cross-sections in Section 3. In Section 4, we solve the problem of growth-induced out-of-plane bifurcation in a chiral rod with guided–guided boundary conditions to study the interplay between chiralities in growth and material laws. We present our conclusions in Section 5.

1.1. Notation

Throughout this text, the indices $i,j,k\in \{1,2,3\}$ and $\alpha ,\beta \in \{1,2\}$ , unless mentioned otherwise. We let $\{{\mathbf{e}}_1,{\mathbf{e}}_2,{\mathbf{e}}_3\}$ be a right-handed, fixed, orthonormal basis for the Euclidean space ${\mathbb{E}}^3$ . Boldface symbols are used to denote tensors, lowercase letters for first-order tensors (e.g. $\mathbf{v}$ ) and uppercase letters for second-order tensors (e.g. $\mathbf{T}$ ). Underlined symbols such as $\underset{¯}{\mathrm{v}}$ and $\underset{¯}{\mathrm{T}}$ denote matrix representation of tensors with respect to a basis.

2. Special Cosserat rod formulation

Consider a straight rod of unit length in its stress-free reference configuration as shown in Figure 1. Assumption of the special Cosserat rod requires the transverse cross-sections to stay rigid during the deformation. Let $s\in [-\frac{1}{2},\frac{1}{2}]$ denote a signed arc-length parameter of the centre-line in the reference configuration. Let $\mathbf{r}(s)$ define the centre-line of the deformed rod. Let $\mathbf{R}(s)\in \mathrm{SO}(3)$ be the rotation of transverse cross-sections in the reference configuration of the rod, mapping the fixed basis $\{{\mathbf{e}}_1,{\mathbf{e}}_2,{\mathbf{e}}_3\}$ to a triad of orthonormal directors given by

$${\mathbf{d}}_i(s)=\mathbf{R}(s){\mathbf{e}}_i.$$

(1)

Figure 1.

Kinematics of a special Cosserat rod depicting the deformed centre-curve and the triad of orthonormal directors.

The vector fields

$$\mathbf{\nu }:=\mathbf{r}\prime ,\mathit{\kappa }:=\mathrm{axial}(\mathbf{R}\prime {\mathbf{R}}^T)$$

(2)

define the convected coordinates $\mathbf{\nu }={\nu }_i{\mathbf{d}}_i$ and $\mathit{\kappa }={\kappa }_i{\mathbf{d}}_i$ with respect to the director frame field, along with the ordered triples $\underset{¯}{\mathrm{v}}:=({\nu }_1,{\nu }_2,{\nu }_3)$ and $\underset{¯}{\mathrm{k}}:=({\kappa }_1,{\kappa }_2,{\kappa }_3)$ . The strains ${\nu }_{\alpha }$ correspond to shear, ${\nu }_3$ corresponds to stretch, ${\kappa }_{\alpha }$ correspond to curvatures, and ${\kappa }_3$ corresponds to twist.

We further assume the rod to be hyperelastic with a differentiable energy density (per unit length) function $\mathrm{\Phi }(\mathbf{r}\prime ,\mathbf{R},\mathbf{R}\prime ,s)$ . Material objectivity allows for a simpler version of energy function in terms of strains [26], given by

$$\mathrm{\Phi }=W(\underset{¯}{\mathrm{v}},\underset{¯}{\mathrm{k}},s),$$

(3)

where W is another differentiable scalar valued function.

The internal force and moment on the transverse cross-section are denoted by $\mathbf{n}(s)=n_i{\mathbf{d}}_i$ and $\mathbf{m}(s)=m_i{\mathbf{d}}_i$ , respectively, along with the corresponding triples $\underset{¯}{\mathrm{n}}:=(n_1,n_2,n_3)$ and $\underset{¯}{\mathrm{m}}:=(m_1,m_2,m_3)$ . The components $n_{\alpha }$ are essentially the shear forces, $n_3$ is axial force, $m_{\alpha }$ are bending moments and $m_3$ is the torsional moment. These are related to the strain components as

$$\underset{¯}{\mathrm{n}}=\frac{\partial W}{\partial \underset{¯}{\mathrm{v}}},\underset{¯}{\mathrm{m}}=\frac{\partial W}{\partial \underset{¯}{\mathrm{k}}}.$$

(4)

To prevent self-penetration, we require

$${\nu }_3=\mathbf{r}\prime ·{\mathbf{d}}_3>0,$$

(5)

and the unshearability constraint is expressed as

$${\nu }_{\alpha }=\mathbf{r}\prime ·{\mathbf{d}}_{\alpha }=0.$$

(6)

2.1. Material symmetry in rods

In this section, we present a brief overview of certain classes of material symmetry for special Cosserat rods, as described by Healey [26, 30].

2.1.1. Helical symmetry

Consider a straight rod possessing helical material symmetry [26]. A unique flip axis (or symmetry axis) is associated with every transverse cross-section that rotates as the section plane moves along the length of the rod (Figure 2(a)).

Figure 2.

A depiction of symmetry and associated kinematic parameters in chiral rods [26].

A 180-degree rotation (flip) about this axis renders the rod the same as before. We denote by $\mathcal{M}\ne 0$ its signed pitch, so defined that $\mathcal{M}>0$ for right-handed helices and $|\mathcal{M}|$ is the least axial translation of the cross-section needed for the flip axis to complete a full rotation in the ${\mathbf{e}}_1$ – ${\mathbf{e}}_2$ plane (Figure 2(b)).

Unlike flips, reflections about a transverse plane do not result in a coincident helix, nor do the reflections through longitudinal planes. In fact, these reflections change the sign of $\mathcal{M}$ , keeping its magnitude, the same.

We introduce a rotating basis field

$$\{{\mathbf{e}}_1^{*}\left(\frac{s}{\mathcal{M}}\right){\mathbf{e}}_2^{*}\left(\frac{s}{\mathcal{M}}\right),{\mathbf{e}}_3^{*}\left(\frac{s}{\mathcal{M}}\right)={\mathbf{e}}_3\}$$

(7)

and a corresponding triad of director fields given by

$${\mathbf{e}}_i^{*}(\varphi )={\mathbf{\Theta }}_{\varphi }{\mathbf{e}}_i,0\le \varphi <2\pi$$

(8)

$${\mathbf{d}}_{\alpha }^{*}(s)=\mathbf{R}(s){\mathbf{e}}_{\alpha }^{*}\left(\frac{s}{\mathcal{M}}\right),$$

(9)

where ${\mathbf{\Theta }}_{\varphi }$ is a proper orthogonal tensor with matrix representation

$${\underset{¯}{\Theta }}_{\varphi }=\left[\begin{array}{ccc}\cos \varphi \hfill & -\sin \varphi \hfill & 0\hfill \\ \multicolumn{1}{c}{\sin \varphi }& \cos \varphi \hfill & 0\hfill \\ \multicolumn{1}{c}{0}& 0\hfill & 1\hfill \end{array}\right],$$

(10)

in the fixed basis.

Assuming ${\mathbf{e}}_1^{*}(\varphi )$ to be the rotating flip axis, we denote by ${\mathbf{H}}_{\varphi }^{\pi }$ the flip about ${\mathbf{e}}_1^{*}(\varphi )$ , so that

$${\mathbf{H}}_{\varphi }^{\pi }{\mathbf{e}}_1^{*}(\varphi )={\mathbf{e}}_1^{*}(\varphi ),{\mathbf{H}}_{\varphi }^{\pi }{\mathbf{e}}_2^{*}(\varphi )=-{\mathbf{e}}_2^{*}(\varphi ),{\mathbf{H}}_{\varphi }^{\pi }{\mathbf{e}}_3^{*}(\varphi )=-{\mathbf{e}}_3.$$

(11)

Material properties with respect to the symmetry axis ${\mathbf{e}}_1^{*}\left(\frac{s}{\mathcal{M}}\right)$ are assumed not to change as the cross-section s moves along the rod. This motivates the definition of a symmetry-adapted energy function [26] independent of s, given by

$$W(\underset{¯}{\mathrm{v}},\underset{¯}{\mathrm{k}},s)=\mathrm{\Phi }=W^{*}({\underset{¯}{\mathrm{v}}}^{*},{\underset{¯}{\mathrm{k}}}^{*}),$$

(12)

where ${\underset{¯}{\mathrm{v}}}^{*}=({\nu }_1^{*},{\nu }_2^{*},{\nu }_3)$ and ${\underset{¯}{\mathrm{k}}}^{*}=({\kappa }_1^{*},{\kappa }_2^{*},{\kappa }_3)$ are a result of the change of coordinates

$$\mathit{\kappa }={\kappa }_{\alpha }^{*}{\mathbf{d}}_{\alpha }^{*}+{\kappa }_3{\mathbf{d}}_3,\mathbf{\nu }={\nu }_{\alpha }^{*}{\mathbf{d}}_{\alpha }^{*}+{\nu }_3{\mathbf{d}}_3.$$

(13)

Helical symmetry is characterized by

$$W^{*}({\nu }_1^{*},{\nu }_2^{*},{\nu }_3,{\kappa }_1^{*},{\kappa }_2^{*},{\kappa }_3)=W^{*}(-{\nu }_1^{*},{\nu }_2^{*},{\nu }_3,-{\kappa }_1^{*},{\kappa }_2^{*},{\kappa }_3),$$

(14)

in terms of the new energy function without s as an argument.

2.1.2. n-fold helical symmetry

Consider a rod with a symmetry analogous to $n\ge 2$ helices entwined together, such that each cross-section at s has n equally spaced flip axes. A 180-degree rotation about each of these gives a symmetry (Figure 2(b)). Such a rod is said to have n-fold dihedral helical symmetry, which is characterized by the condition

$$W^{*}(-{\underset{¯}{H}}_{\frac{2\pi }{n}}^{\pi *}{\underset{¯}{\mathrm{v}}}^{*},-{\underset{¯}{H}}_{\frac{2\pi }{n}}^{\pi *}{\underset{¯}{\mathrm{k}}}^{*})=W^{*}({\underset{¯}{\mathrm{v}}}^{*},{\underset{¯}{\mathrm{k}}}^{*}),$$

(15)

in addition to equation (14), where ${\underset{¯}{H}}_{\frac{2\pi }{n}}^{\pi *}$ is the matrix of ${\mathbf{H}}_{\frac{2\pi }{n}}^{\pi }$ with respect to the rotating basis (equation (7)):

$${\underset{¯}{H}}_{\frac{2\pi }{n}}^{\pi *}=\left[\begin{array}{ccc}\cos \left(\frac{2\pi }{n}\right)\hfill & \sin \left(\frac{2\pi }{n}\right)\hfill & 0\hfill \\ \multicolumn{1}{c}{\sin \left(\frac{2\pi }{n}\right)}& -\cos \left(\frac{2\pi }{n}\right)\hfill & 0\hfill \\ \multicolumn{1}{c}{0}& 0\hfill & -1\hfill \end{array}\right].$$

(16)

2.1.3. Continuous helical symmetry

For $n>>1$ , a straight rod with n-fold dihedral helical symmetry approaches so-called continuous helical symmetry. In this type of symmetry all vectors of the cross-section act as symmetry axes, or equivalently any fixed flip axis, say ${\mathbf{e}}_1$ , acts as a symmetry axis for all cross-sections. Continuous helical symmetry can be characterized by

$$W(-{\underset{¯}{H}}_{\varphi }^{\pi }\underset{¯}{\mathrm{v}},-{\underset{¯}{H}}_{\varphi }^{\pi }\underset{¯}{\mathrm{k}})=W(\underset{¯}{\mathrm{v}},\underset{¯}{\mathrm{k}}),\forall \varphi \in [0,\pi ).$$

(17)

2.1.4. Transverse hemitropy and isotropy

Let $\mathbf{E}$ denote the reflection tensor with matrix

$$\underset{¯}{E}=\left[\begin{array}{ccc}1\hfill & 0\hfill & 0\hfill \\ \multicolumn{1}{c}{0}& -1\hfill & 0\hfill \\ \multicolumn{1}{c}{0}& 0\hfill & 1\hfill \end{array}\right],$$

(18)

written in the fixed basis. A homogeneous hyperelastic straight rod with energy function $W(\underset{¯}{\mathrm{v}},\underset{¯}{\mathrm{k}})$ is transversely hemitropic if

$$W({\underset{¯}{\Theta }}_{\varphi }\underset{¯}{\mathrm{v}},{\underset{¯}{\Theta }}_{\varphi }\underset{¯}{\mathrm{k}})=W(\underset{¯}{\mathrm{v}},\underset{¯}{\mathrm{k}})\forall \varphi \in [0,2\pi ),$$

(19)

and flip-symmetric if

$$W(\underset{¯}{E\mathrm{v}},\underset{¯}{E}\underset{¯}{\mathrm{k}})=W(\underset{¯}{\mathrm{v}},\underset{¯}{\mathrm{k}}).$$

(20)

Note that flip-symmetry does not belong to the class of transverse symmetry, defined by Healey [26]. A straight rod is transversely isotropic if in addition to equation (19), it also satisfies

$$W(\underset{¯}{\mathrm{v}},\underset{¯}{\mathrm{k}})=W(\underset{¯}{E\mathrm{v}},-\underset{¯}{E}\underset{¯}{\mathrm{k}}).$$

(21)

Flip-symmetric hemitropy is equivalent to continuous helical symmetry [26]. Another way to obtain flip-symmetric hemitropy is to consider a rod with helical symmetry and take the limit $\mathcal{M}\to 0$ [30].

2.2. Energy function

The energy density per unit length of unshearable hemitropic rods can be expressed as [1, 26]

$$W=\mathrm{\Upsilon }({\kappa }_{\alpha }{\kappa }_{\alpha },{\nu }_3,{\kappa }_3),$$

(22)

where $\mathrm{\Upsilon }$ is a scalar valued function. This representation is also valid for flip-symmetry. For calculations in this paper, we adopt a model considered by Papadopoulos and Healey [31, 32], defined as

$$\mathrm{\Upsilon }=\frac{1}{2}[\mathrm{\Psi }({\nu }_3)+2A[{\nu }_3-1]{\kappa }_3+B{\kappa }_3^2+C{\kappa }_{\alpha }{\kappa }_{\alpha }],$$

(23)

where $\mathrm{\Psi }:(0,\infty )\to \mathbb{R}$ is a function such that $g:=\frac{1}{2}\mathrm{\Psi }\prime$ obeys $g({\nu }_3)\to -\infty$ as ${\nu }_3\to 0$ . The function $g(·)$ allows us to modify the axial force response of the model, and it must satisfy $g(1)=0$ . The constant C corresponds to bending stiffness, $B-\frac{A^2}{g\prime (1)}$ is equivalent to torsional rigidity and $g\prime (1)-\frac{A^2}{B}$ to axial stiffness, where A is the degree of hemitropy. We assume $B>0$ , $C>0$ and $\mathrm{Bg}\prime ({\nu }_3)>A^2$ for all ${\nu }_3$ to ensure convexity. This in turn implies that $g(·)$ should be monotonic and hence invertible. For example, a response function satisfying all our criteria can be chosen as [31]

$$g({\nu }_3)=F\ln ({\nu }_3)+\frac{A^2}{B}[{\nu }_3-1],$$

(24)

where $F>0$ is a constant. This energy allows for infinite compressive axial force $n_3\to -\infty$ whenever an unrealistically extreme strain ${\nu }_3\to 0$ is present.

As demonstrated by Healey [26], quadratic energy functions are incapable of distinguishing between different types of n-fold helical symmetry ( $n\ge 3$ ) and hemitropy. Along similar lines, the energy function (23) can be shown to be applicable to n-fold helical symmetry.

3. Growth formulation

Growth in elastic bodies is typically modelled by introduction of a multiplicative decomposition of the deformation gradient into pure growth and pure elastic deformation parts [15, 34]. This decomposition assumes a virtual stress-free incompatible configuration. For one-dimensional structures where growth manifests as increase in overall length, first the stress-free rod isolated from its environment and boundary conditions can be allowed to grow free into a virtual state, and then the boundary and environmental factors can be forcibly imposed [18, 35].

One-dimensional growth models, where cross-sections simply translate during free growth, are not suitable for several classes of chiral rods. Chiral rods usually have a physical winding bias intrinsic to the microstructure [26]. Length-wise growth with no cross-sectional rotation can modify this microstructure. For example, a rod with helical symmetry made to grow axially will have to change its inherent pitch if the cross-sections are not allowed to rotate during growth; and as a result the constitutive parameters controlling material chirality must change accordingly. In such examples, to be able to look at growth that does not alter the microstructure, or restricts the microstructure to modify itself in a particular manner, it is essential that we look at rod growth in a more general setup.

Similarly, in chiral rods where material symmetry arises from fibre reinforcement [36], rod growth is a result of individual fibre growth and it is the growth pattern of these fibres that dictates whether the rod’s cross-sections must rotate, as they are translated during growth. Consider a rod that is composed of fibres twisted helically in the unstressed reference state; if the cross-sections are not allowed to rotate during growth, it would have an unwinding or over-winding effect on the fibres, thus generating stresses. In such cases, for growth to take place without the generation of any stress, the cross-sections must rotate.

This is why we choose to individually treat all three configurations (reference, virtual and current) as special Cosserat rods, and then analyse the relative rotations.

3.1. General framework for growing rods

Let ${\mathcal{R}}_o$ denote the initial stress-free reference configuration of the rod, occupying $\{S{\mathbf{e}}_3:-\frac{1}{2}\le S\le \frac{1}{2}\}$ and denote by S a signed arc-length parameter of the centre-line in ${\mathcal{R}}_o$ . Let $\tilde{\mathbf{r}}(S)$ be the curve taken by the centre-line in the virtual grown configuration $\tilde{\mathcal{R}}$ , such that the point $S{\mathbf{e}}_3$ in ${\mathcal{R}}_o$ gets mapped to $\tilde{\mathbf{r}}(S)$ in $\tilde{\mathcal{R}}$ (Figure 3). The virtual configuration is assumed to be stress-free. We define a signed arc-length $s(S)$ in $\tilde{\mathcal{R}}$ by

$$s(S):={\int }_0^S\parallel \tilde{\mathbf{r}}\prime (\tau )\parallel d\tau ,$$

(25)

where ∥·∥ denotes the Euclidean vector norm.

Figure 3.

Kinematics of an initially straight rod growing from origin $S_o$ , depicting the configurations: reference ${\mathcal{R}}_o$ , virtual $\tilde{\mathcal{R}}$ and current $\mathcal{R}$ , along with the multiplicative decomposition $\mathbf{Q}=\mathbf{RW}$ .

We denote the transverse cross-section at S in ${\mathcal{R}}_o$ by ${\mathrm{\Gamma }}_o(S)$ and let it get mapped to $\tilde{\mathrm{\Gamma }}(S)$ in the virtual configuration $\tilde{\mathcal{R}}$ . We define $\mathbf{W}(S)\in \mathrm{SO}(3)$ to be the rotation of $\tilde{\mathrm{\Gamma }}(S)$ with respect to ${\mathrm{\Gamma }}_o(S)$ , and let it map the fixed basis $\{{\mathbf{e}}_1,{\mathbf{e}}_2,{\mathbf{e}}_3\}$ to a virtual director field given by

$${\tilde{\mathbf{e}}}_i(S)=\mathbf{W}(S){\mathbf{e}}_i.$$

(26)

When the boundary conditions and environmental factors are imposed, let the centre-line take the curve $\mathbf{r}(S)$ in the current configuration $\mathcal{R}$ , and the cross-section $\tilde{\mathrm{\Gamma }}(S)$ in $\tilde{\mathcal{R}}$ be mapped to $\mathrm{\Gamma }(S)$ in $\mathcal{R}$ . We define $\mathbf{R}(S)\in \mathrm{SO}(3)$ to be the rotation of $\mathrm{\Gamma }(S)$ with respect to $\tilde{\mathrm{\Gamma }}(S)$ and $\mathbf{Q}(S)\in \mathrm{SO}(3)$ to be the net rotation of $\mathrm{\Gamma }(S)$ with respect to ${\mathrm{\Gamma }}_o(S)$ , so that

$$\mathbf{Q}(S)=\mathbf{R}(S)\mathbf{W}(S).$$

(27)

The virtual director field is transformed into another director field in the current configuration, given by

$${\mathbf{d}}_i(S)=\mathbf{R}(S){\tilde{\mathbf{e}}}_i(S)=\mathbf{Q}(S){\mathbf{e}}_i.$$

(28)

All the maps we have introduced are assumed to be smooth for the sake of convenience. Analogous to $\mathbf{r}:[-\frac{1}{2},\frac{1}{2}]\to {\mathbb{E}}^3$ , we define another map $\hat{\mathbf{r}}:[s(-\frac{1}{2}),s\left(\frac{1}{2}\right)]\to {\mathbb{E}}^3$ to denote the same curve via the reparametrization,

$$\mathbf{r}(S)=(\hat{\mathbf{r}}°s)(S).$$

(29)

This implies

$$\mathbf{r}\prime (S)=\parallel \tilde{\mathbf{r}}\prime (S)\parallel \frac{\partial \hat{\mathbf{r}}}{\partial s},$$

(30)

where $(·)\prime :=\frac{\partial }{\partial S}(·)$ , taking the magnitude of which gives the one-dimensional multiplicative decomposition

$$\Vert \frac{\partial \mathbf{r}}{\partial S}\Vert =\Vert \frac{\partial \tilde{\mathbf{r}}}{\partial S}\Vert \Vert \frac{\partial \hat{\mathbf{r}}}{\partial s}\Vert .$$

(31)

Similarly, we define $\hat{\mathbf{R}}:[s(-\frac{1}{2}),s\left(\frac{1}{2}\right)]\to \mathrm{SO}(3)$ by

$$\mathbf{R}(S)=(\hat{\mathbf{R}}°s)(S).$$

(32)

We assume the transverse cross-sections to remain orthogonal to the centre-line in both virtual and current configurations, hence the conditions

$$\tilde{\mathbf{r}}\prime ·{\tilde{\mathbf{e}}}_{\alpha }=0$$

(33)

$$\mathrm{and}\mathbf{r}\prime ·{\mathbf{d}}_{\alpha }=0$$

(34)

must hold, where equation (34) is equivalent to the unshearability constraint (6). The symbols and notations introduced in this section are pictorially represented in Figure 3.

3.1.1. Homogeneous growth kinematics

We consider the growth to be homogeneous throughout the rod. This assumption leads to the following constraints:

• The length-wise growth parameter denoted by $\gamma :=||\tilde{\mathbf{r}}\prime (S)||$ is a constant – that is, it is independent of S.

• Let $h\in \mathbb{R}$ be such that $0<|h|<1$ . Consider the relative rotation of cross-section $\tilde{\mathrm{\Gamma }}(S+h)$ with respect to $\tilde{\mathrm{\Gamma }}(S)$ :

$${\tilde{\mathbf{e}}}_i(S+h)=\mathbf{W}(S+h)\mathbf{W}(S{)}^{-1}{\tilde{\mathbf{e}}}_i(S).$$

(35)

For all permissible h, the relative rotation $\mathbf{W}(S+h)\mathbf{W}(S{)}^{-1}$ is assumed to be independent of S, and hence can be denoted as a function of h only:

$$\mathbf{W}(S+h)\mathbf{W}(S{)}^{-1}=:\mathbf{\Pi }(h).$$

(36)

This gives us the decomposition

$$\mathbf{W}(S+h)=\mathbf{\Pi }(h)\mathbf{W}(S).$$

(37)

Choosing $\hslash \ne 0$ such that all the tensor fields appearing in the following calculation make sense, we have

$$[\mathbf{\Pi }(h+\hslash )-\mathbf{\Pi }(h)]\mathbf{W}(S)=\mathbf{W}(S+h+\hslash )-\mathbf{W}(S+h)$$

(38a)

$$=\mathbf{\Pi }(h)[\mathbf{W}(S+\hslash )-\mathbf{W}(S)].$$

(38b)

Dividing by $\hslash$ and taking the limit $\hslash \to 0$ yields

$$\frac{\partial \mathbf{\Pi }(h)}{\partial h}\mathbf{W}(S)=\mathbf{\Pi }(h)\frac{\partial \mathbf{W}(S)}{\partial S}.$$

(39)

We define the following tensor fields for our convenience:

$$\mathbf{\Lambda }(S):=\mathbf{W}(S{)}^T\frac{\partial \mathbf{W}(S)}{\partial S},\mathrm{and}\mathbf{\Omega }(S):=\frac{\partial \mathbf{W}(S)}{\partial S}\mathbf{W}(S{)}^T.$$

(40)

Now equation (39) implies

$$\mathbf{\Pi }(h{)}^T\frac{\partial \mathbf{\Pi }(h)}{\partial h}=\mathbf{\Omega }(S).$$

(41)

Since h and S can be chosen arbitrarily, independent of each other, we conclude that $\mathbf{\Omega }(S)$ is constant. Another way to interpret equation (36) is to set

$$\frac{\partial }{\partial S}[\mathbf{W}(S+h)\mathbf{W}(S{)}^{-1}]=\mathbf{O}.$$

(42)

We expand the derivative to get

$$\frac{\partial }{\partial S}[\mathbf{W}(S+h)\mathbf{W}(S{)}^{-1}]=\frac{\partial }{\partial S}[\mathbf{W}(S+h)]\mathbf{W}(S{)}^T-\mathbf{W}(S+h)\mathbf{W}(S{)}^T\frac{\partial \mathbf{W}(S)}{\partial S}\mathbf{W}(S{)}^T$$

(43a)

$$=\mathbf{W}(S+h)[\mathbf{\Lambda }(S+h)-\mathbf{\Lambda }(S)]\mathbf{W}(S{)}^T.$$

(43b)

This implies $\mathbf{\Lambda }(S+h)=\mathbf{\Lambda }(S)$ for all choices of S and h, chosen independent of each other, which means $\mathbf{\Lambda }$ is constant.

Moreover, differentiating the condition $\mathbf{W}{\mathbf{W}}^T=\mathbf{I}={\mathbf{W}}^T\mathbf{W}$ of orthogonality with respect to S yields

$${\mathbf{\Lambda }}^T=-\mathbf{\Lambda }\mathrm{and}{\mathbf{\Omega }}^T=-\mathbf{\Omega };$$

(44)

that is, $\mathbf{\Lambda }$ and $\mathbf{\Omega }$ are skew-symmetric.

• We fix a point on the centre-line that gets mapped to itself under the growth transformation, along with its corresponding cross-section. Thus, we assume the existence of a point $S_o\in [-\frac{1}{2},\frac{1}{2}]$ satisfying

$$\tilde{\mathbf{r}}(S_o)=S_o{\mathbf{e}}_3\mathrm{and}\mathbf{W}(S_o)=\mathbf{I}.$$

(45)

This can also be interpreted as if the rod is allowed to grow while being held at $S_o$ (origin of growth). It is held in such a way that no incompatibility or stress is caused due to growth. We define vectors $\mathit{a}:=\mathrm{axial}(\mathbf{\Lambda })$ and $\mathit{\omega }:=\mathrm{axial}(\mathbf{\Omega })$ ; these are actually constant vectors and can be related by

$$\mathit{\omega }=\mathbf{W}(S)\mathit{a}.$$

(46)

Since this is also satisfied for the specific point $S=S_o$ , we imply $\mathit{a}=\mathit{\omega }$ and $\mathbf{\Lambda }=\mathbf{\Omega }$ . This also means that $\mathrm{axis}(\mathbf{W}(S))=\mathit{a}$ for all S.

Thus one can solve the system

$${\mathbf{W}}^T\frac{\partial \mathbf{W}}{\partial S}=\mathbf{\Lambda }\mathrm{with}\mathbf{W}(S_o)=\mathbf{I},$$

(47)

for $\mathbf{W}$ to obtain (Appendix A)

$$\mathbf{W}(S)=e^{(S-S_o)\mathbf{\Lambda }},$$

(48)

where tensor exponential is defined by the series expansion of exponential function.

3.1.2. Extension to a general growing curve

Consider a general scenario in which the initial configuration ${\mathcal{R}}_o$ is a special Cosserat rod. Let $\overline{\mathbf{r}}:[-\frac{1}{2},\frac{1}{2}]\to {\mathbb{E}}^3$ be its centre-curve, where $\overline{\mathbf{r}}(S)$ is arc-length parametrized. Let $\overline{\mathbf{W}}(S)\in \mathrm{SO}(3)$ denote the orientation of ${\mathrm{\Gamma }}_o(S)$ with respect to the fixed basis, mapping those to an orthonormal director field ${\overline{\mathbf{e}}}_i(S):=\overline{\mathbf{W}}(S){\mathbf{e}}_i$ associated with initial configuration. Since $\mathbf{W}(S)\in \mathrm{SO}(3)$ maps ${\mathrm{\Gamma }}_o(S)$ to $\tilde{\mathrm{\Gamma }}(S)$ , the virtual director field should be given by ${\tilde{\mathbf{e}}}_i(S)=\mathbf{W}(S)\overline{\mathbf{W}}(S){\mathbf{e}}_i$ , so that equation (35) is modified as

$${\tilde{\mathbf{e}}}_i(S+h)=\mathbf{W}(S+h)\overline{\mathbf{W}}(S+h)\overline{\mathbf{W}}(S{)}^{-1}\mathbf{W}(S{)}^{-1}{\tilde{\mathbf{e}}}_i(S).$$

(49)

Homogeneous growth law still requires $\gamma$ to be constant. The tensor $\mathbf{W}(S+h)\mathbf{W}(S{)}^{-1}$ is again independent of S. In addition, the rod is assumed to be held at $S_o\in [-\frac{1}{2},\frac{1}{2}]$ while growing, so that we have

$$\tilde{\mathbf{r}}(S_o)=\overline{\mathbf{r}}(S_o)\mathrm{and}\mathbf{W}(S_o)=\mathbf{I}.$$

(50)

This assumption, along with the kind of homogeneity used in induced rotations, gives such a $\mathbf{W}(S)$ that makes all the cross-sections rotate about the particular axis a . Moreover, the solution is given by equation (48), which in turn implies

$${\tilde{\mathbf{e}}}_i(S)=e^{(S-S_o)\mathbf{\Lambda }}\overline{\mathbf{W}}(S){\mathbf{e}}_i.$$

(51)

In fact, the constant vector $\mathit{\omega }=a$ can be treated as the growth parameter controlling relative rotation of cross-sections while $\gamma$ controls the lengthwise growth as in the former case. Whenever the centre-curves are normal to the cross-sections, throughout ${\mathcal{R}}_o$ and $\tilde{\mathcal{R}}$ , we deduce

$$\tilde{\mathbf{r}}(S)=\overline{\mathbf{r}}(S_o)+\gamma {\int }_{S_o}^S{e}^{(\tau -S_o)\mathbf{\Lambda }}\overline{\mathbf{r}}\prime (\tau )d\tau .$$

(52)

We emphasise that equations (51) and (48) do not assume the respective centre-curves to be normal to the cross-sections in ${\mathcal{R}}_o$ and $\tilde{\mathcal{R}}$ .

3.2. Growth in straight rods

Consider a straight rod with flip-symmetric hemitropy in its reference configuration. A straight virtual configuration condenses to

$$\tilde{\mathbf{r}}(S)=[S_o+\gamma [S-S_o]]{\mathbf{e}}_3,$$

(53)

which, with the aid of equation (33), results in

$$\mathbf{W}(S){\mathbf{e}}_3={\mathbf{e}}_3\forall S\in [-\frac{1}{2},\frac{1}{2}].$$

(54)

This indicates that $\mathit{\omega }$ is along ${\mathbf{e}}_3$ . We introduce another growth parameter $\omega$ , defined by

$$\mathit{\omega }=\omega {\mathbf{e}}_3,$$

(55)

so that its corresponding skew tensor is

$$\mathbf{\Omega }=\omega \mathbf{A},\mathrm{with}\mathbf{A}={\mathbf{e}}_2\otimes {\mathbf{e}}_1-{\mathbf{e}}_1\otimes {\mathbf{e}}_2.$$

(56)

Since the rotation tensor can also be expressed as

$${\mathbf{\Theta }}_{\varphi }=e^{\varphi \mathbf{A}},$$

(57)

we get

$$\mathbf{W}(S)={\mathbf{\Theta }}_{(S-S_o)\omega }.$$

(58)

The parameters $\gamma$ and $\omega$ capture all the necessary information regarding growth. It is evident that $\gamma >1$ denotes growth while $\gamma <1$ denotes atrophy. Similarly, $\omega$ and $-\omega$ signify two opposite cross-sectional rotations caused by growth while $\omega =0$ indicates no growth-induced rotation.

This type of growth is helical in nature. Consider any line in the bulk of the rod parallel to its axis, but not the axis itself. As the rod grows this line transforms into a helix of pitch $\frac{\gamma }{\omega }$ . Similarly, any helix in the initial configuration transforms into another coaxial helix, not necessarily with the same pitch. This is a reflection of the fact that the Darboux vector of a helix is a constant vector aligned along its axis.

3.2.1. Growth law

The growth law adopted here considers rotation of cross-sections with respect to each other in the due course of growth. Consider a rotating basis field $\{{\mathbf{e}}_1^{*}(S),{\mathbf{e}}_2^{*}(S),{\mathbf{e}}_3^{*}(S)={\mathbf{e}}_3\}$ given by

$${\mathbf{e}}_i^{*}(S)={\mathbf{\Theta }}_{\frac{S}{\mathcal{M}}}{\mathbf{e}}_i,$$

(59)

representing a helix embedded in the initial configuration of a rod. As the rod grows this transforms into $\mathbf{W}(S){\mathbf{e}}_i^{*}(S)$ in the virtual configuration. Let us denote this by a basis field $\{{\mathbf{f}}_1^{*}(s),{\mathbf{f}}_2^{*}(s),{\mathbf{f}}_3^{*}(s)\}$ defined on the virtual arc-length parameter by

$$\mathbf{W}(S){\mathbf{e}}_i^{*}(S)=:({\mathbf{f}}_i^{*}°s)(S).$$

(60)

This is equivalent to

$${\mathbf{f}}_i^{*}(s)={\mathbf{\Theta }}_{\frac{s}{\gamma \mathcal{M}}+[\frac{s}{\gamma }-S_o]\omega }{\mathbf{e}}_i.$$

(61)

Let $h\ne 0$ be such that ${\mathbf{e}}_i^{*}(S+h)$ and ${\mathbf{f}}_i^{*}(s+h)$ are well defined, then we obtain

$${\mathbf{e}}_i^{*}(S+h)={\mathbf{\Theta }}_{\frac{h}{\mathcal{M}}}{\mathbf{e}}_i^{*}(S),$$

(62)

$${\mathbf{f}}_i^{*}(s+h)={\mathbf{\Theta }}_{\frac{h}{\gamma }[\frac{1}{\mathcal{M}}+\omega ]}{\mathbf{f}}_i^{*}(s).$$

(63)

This shows that our chosen growth map transforms the initial helix with pitch $\mathcal{M}$ into another helix with pitch, say $\mu$ , which can be expressed as

$$\mu =\frac{\gamma \mathcal{M}}{1+\omega \mathcal{M}}.$$

(64)

This motivates us to define a symmetry-preserving growth law for rods possessing helical symmetry.

Rods with helical symmetry. Consider a rod which, due to its microstructure, possesses simple helical symmetry or n-fold helical symmetry. Let $\mathcal{M}$ be the pitch associated with its microstructure. Once growth parameters $\gamma$ and $\omega$ are known, equation (64) serves as an evolution law for the pitch of its microstructure.

We introduce the idea of symmetry-preserving growth – wherein the growth map fixes all helices with pitch the same as that of the microstructure ( $\mu =\mathcal{M}$ ). Thus, for rods with a pitch associated with their microstructure we have the following helical growth law

$$\gamma =1+\omega \mathcal{M},$$

(65)

where $\gamma$ is the only growth parameter and $\mathcal{M}$ comes from the material symmetry. For rods having helical symmetry, this assumption of symmetry-preserving growth provides a rationale for relative rotation of cross-sections during growth.

Hemitropic rods. Although there are different versions [26, 30] of how helical symmetry can be used to arrive at hemitropy, there is no pitch directly associated with transverse hemitropy (equation (19)). So, for hemitropic rods (and even isotropic), one may use the same helical growth law (65) without any notion of microstructural pitch, in which case both $\gamma$ and $\mathcal{M}$ are independent growth parameters. For such a growth law, all helices with pitch $\mathcal{M}$ remain unaltered under the growth map, so we denote it as the characteristic pitch of growth.

3.2.2. Calculation of strains

The grown configuration is obtained by imposing environmental and boundary effects on the virtual stress-free configuration. Hence the strain energy is a function of $\frac{\partial \hat{\mathbf{r}}}{\partial s}$ , $\hat{\mathbf{R}}$ and $\frac{\partial \hat{\mathbf{R}}}{\partial s}$ . We define the vector fields

$$\hat{\mathbf{\nu }}=\frac{\partial \hat{\mathbf{r}}}{\partial s}\mathrm{and}\hat{\mathit{\kappa }}=\mathrm{axial}\left(\frac{\partial \hat{\mathbf{R}}}{\partial s}{\hat{\mathbf{R}}}^T\right).$$

(66)

Let their components be $\hat{\mathbf{\nu }}={\hat{\nu }}_i{\mathbf{d}}_i$ and $\hat{\mathit{\kappa }}={\hat{\kappa }}_i{\mathbf{d}}_i$ with respect to the director frame in the current configuration. Consider the derivative

$$\frac{\partial {\mathbf{d}}_i}{\partial S}=\frac{\partial \mathbf{Q}}{\partial S}{\mathbf{Q}}^{-1}{\mathbf{d}}_i=[\frac{\partial \mathbf{R}}{\partial S}\mathbf{W}+\mathbf{R}\frac{\partial \mathbf{W}}{\partial S}]{\mathbf{W}}^{-1}{\mathbf{R}}^{-1}{\mathbf{d}}_i$$

(67a)

$$=\frac{\partial \mathbf{R}}{\partial s}\frac{\partial s}{\partial S}{\mathbf{R}}^{-1}{\mathbf{d}}_i+\mathbf{R}\frac{\partial \mathbf{W}}{\partial S}{\mathbf{W}}^{-1}{\tilde{\mathbf{e}}}_i$$

(67b)

$$=\gamma \hat{\mathit{\kappa }}\times {\mathbf{d}}_i+\mathbf{R}[\mathit{\omega }\times {\tilde{\mathbf{e}}}_i]$$

(67c)

$$=[\gamma \hat{\mathit{\kappa }}+\mathbf{R}\mathit{\omega }]\times {\mathbf{d}}_i.$$

(67d)

Now define the axial vector $\mathit{\beta }:=\mathrm{axial}\left(\frac{\partial \mathbf{Q}}{\partial S}{\mathbf{Q}}^{-1}\right)$ , which, along with the straight growth assumption, implies

$$\mathit{\beta }=\gamma \hat{\mathit{\kappa }}+\omega {\mathbf{d}}_3.$$

(68)

Given the growth parameters, this relation will be used in retracting the actual strains from the apparent curvature $\mathit{\beta }$ . Corresponding to $\hat{\mathbf{\nu }}$ and $\hat{\mathit{\kappa }}$ we define

$$\mathbf{\nu }=\mathbf{r}\prime \mathrm{and}\mathit{\kappa }=\mathrm{axial}(\mathbf{R}\prime {\mathbf{R}}^T),$$

(69)

along with their convected components $\mathbf{\nu }={\nu }_i{\mathbf{d}}_i$ and $\mathit{\kappa }={\kappa }_i{\mathbf{d}}_i$ . These speeds and curvatures can be related to the actual strains by

$${\nu }_i=\gamma {\hat{\nu }}_i\mathrm{and}{\kappa }_i=\gamma {\hat{\kappa }}_i.$$

(70)

With the unshearability constraint in place, we have ${\nu }_{\alpha }={\hat{\nu }}_{\alpha }=0$ , while ${\nu }_3=\gamma {\hat{\nu }}_3$ represents the multiplicative decomposition for lengthwise growth. Using the energy density function (23), internal force $\mathbf{n}(S)=n_i(S){\mathbf{d}}_i(S)$ and moment $\mathbf{m}(S)=m_i(S){\mathbf{d}}_i(S)$ in the current configuration can be related to the strains as follows:

$$n_3=g({\hat{\nu }}_3)+A{\hat{\kappa }}_3,$$

(71)

$$m_3=A[{\hat{\nu }}_3-1]+B{\hat{\kappa }}_3,$$

(72)

$$m_{\alpha }=C{\hat{\kappa }}_{\alpha }.$$

(73)

3.2.3. Equilibrium equations

The local linear and angular momentum balance equations for static equilibrium [18, 35] are as follows:

$$\frac{\partial \mathbf{n}}{\partial s}+\mathbf{f}=0,$$

(74)

$$\frac{\partial \mathbf{m}}{\partial s}+\frac{\partial \mathbf{r}}{\partial s}\times \mathbf{n}+\mathbf{l}=0,$$

(75)

where $\mathbf{f}$ and $\mathbf{l}$ respectively denote the body force and body moment per unit virtual arc-length. The change of variable $\frac{\partial }{\partial s}(·)=\frac{1}{\gamma }\frac{\partial }{\partial S}(·)$ to reference coordinates results in

$$\mathbf{n}\prime +\gamma \mathbf{f}=0,$$

(76)

$$\mathbf{m}\prime +\mathbf{r}\prime \times \mathbf{n}+\gamma \mathbf{l}=0.$$

(77)

4. Growing rod with guided–guided ends

A fixed–fixed rod subject to axial displacement or load is constrained both axially and rotationally, and is known to buckle out-of-plane with a transversely hemitropic constitutive law. But in a fixed-free hemitropic rod subject to an axial load, material chirality does not lead to any chiral deformation and the solution is always planar [32]. For a growing rod, there can be another intermediate boundary condition pair – with guided ends – which is rotationally constrained, but axially free at both ends. A guided boundary condition is equivalent to fixing the end of the rod to a block constrained by a slot to translate only along the rod’s axis (Figure 4). In this section we show that a growing rod with guided ends can give rise to non-planar chiral solutions by itself, without any additional load. We use the energy function (23) and the growth law (65) to model the rod. Even though all the calculations would be similar, the results can be discussed separately for two different problems – first, a hemitropic rod, and second, a rod with n-fold helical symmetry.

Figure 4.

Schematic of guided–guided boundary condition. The guides arrest all degrees of freedom at the ends except for axial translation.

The linear and angular momentum balance equations are

$$\frac{d}{\mathrm{ds}}[n_{\alpha }\mathbf{Q}{\mathbf{e}}_{\alpha }+[g({\hat{\nu }}_3)+A{\hat{\kappa }}_3]\mathbf{Q}{\mathbf{e}}_3]=0,$$

(78)

$$\begin{array}{c}\frac{d}{\mathrm{ds}}[C{\hat{\kappa }}_{\alpha }\mathbf{Q}{\mathbf{e}}_{\alpha }+[A[{\hat{\nu }}_3-1]+B{\hat{\kappa }}_3]\mathbf{Q}{\mathbf{e}}_3]\hfill \\ +\hat{\mathbf{r}}\prime \times [n_{\alpha }\mathbf{Q}{\mathbf{e}}_{\alpha }+[g({\hat{\nu }}_3)+A{\hat{\kappa }}_3]\mathbf{Q}{\mathbf{e}}_3]=0,\hfill \end{array}$$

(79)

along with the boundary conditions

$$\mathbf{n}(\pm \frac{1}{2})·{\mathbf{e}}_3=0,$$

(80)

$$\mathbf{r}(\pm \frac{1}{2})·{\mathbf{e}}_{\alpha }=0$$

(81)

$$\mathrm{and}\mathbf{Q}(\pm \frac{1}{2})=\mathbf{I}.$$

(82)

The unshearability constraint (34) results in

$$\mathbf{r}\prime ·\mathbf{Q}{\mathbf{e}}_{\alpha }=0.$$

(83)

Equations (78)–(83) comprise our boundary value problem to be solved for the fields $\mathbf{r}$ , $\mathbf{R}$ and $n_{\alpha }$ . Since we have not imposed any sort of axial constraint, with these sets of boundary conditions we will get a family of solutions differing by a scalar multiple of ${\mathbf{e}}_3$ .

The rod is assumed to be of unit length; thus, all the kinematic quantities are dimensionless by default. The components of internal force, internal moment, material constants $A,B$ and the response function $g(·)$ can all be non-dimensionalized against C by either dividing the concerned quantities in equations (71)–(73) by C, or equivalently setting $C=1$ in the boundary value problem (94)–(101). We follow the bifurcation analysis methodology presented by Healey and Papadopoulos [32] and Smith and Healey [37], wherein first a primary solution is determined which is then perturbed and the boundary value problem is re-derived in terms of the perturbations to get linearized equations.

4.1. The straight solution

Let us consider a solution in which the rod always remains straight while growing, given by

$$\mathbf{r}(S)=\lambda S{\mathbf{e}}_3,\mathbf{Q}(S)=\mathbf{I},n_{\alpha }(S)=0,$$

(84)

where $S\in [-\frac{1}{2},+\frac{1}{2}]$ . This solution has its local force, moment and strain fields as follows:

$$\hat{\mathbf{\nu }}(s)=\frac{\lambda }{\gamma }{\mathbf{e}}_3,$$

(85)

$$\hat{\mathit{\kappa }}(s)=-\frac{\omega }{\gamma }{\mathbf{e}}_3,$$

(86)

$$\mathbf{n}(S)=[g\left(\frac{\lambda }{\gamma }\right)-A\frac{\omega }{\gamma }]{\mathbf{e}}_3,$$

(87)

$$\mathbf{m}(S)=[A[\frac{\lambda }{\gamma }-1]-B\frac{\omega }{\gamma }]{\mathbf{e}}_3.$$

(88)

For such a solution to comply with the force boundary condition (80), we require $\lambda$ to satisfy

$$g\left(\frac{\lambda }{\gamma }\right)=A\frac{\omega }{\gamma },$$

(89)

where A denotes the degree of hemitropy and $g(·)$ is the axial force response function. As we approach a ‘no-growth’ stage, the strain ${\hat{\nu }}_3=\frac{\lambda }{\gamma }\to 1$ and the ratio

$$\frac{m_3}{{\hat{\kappa }}_3}=B-\frac{A^2[\frac{\lambda }{\gamma }-1]}{g\left(\frac{\lambda }{\gamma }\right)}$$

(90)

approaches $B-\frac{A^2}{g\prime (1)}$ , which represents the torsional rigidity [31].

4.2. Perturbed solution

Consider a first-order perturbation of the straight solution (with $0<\epsilon <<1$ ) given by

$$\mathbf{r}(S)=\lambda S{\mathbf{e}}_3+\epsilon \mathit{\rho }(S),$$

(91)

$$\mathbf{Q}(S)=e^{\epsilon \mathbf{\Psi }(S)},$$

(92)

$$n_{\alpha }(S)=\epsilon {\eta }_{\alpha }(S),$$

(93)

where $\mathbf{\Psi }(S)$ is skew-symmetric with $\mathrm{axial}(\mathbf{\Psi })=:\mathit{\psi }$ . We require these perturbed fields to satisfy our boundary value problem. Substituting the perturbations (91)–(93) into our boundary value problem (78)–(83) results in the following linearized problem:

$${\eta \prime }_{\alpha }{\mathbf{e}}_{\alpha }=0,$$

(94)

$$\left[\mathit{\psi }″+\omega {\mathbf{e}}_3\times \mathit{\psi }\prime \right]·{\mathbf{e}}_{\alpha }{\mathbf{e}}_{\alpha }+[A[\lambda -\gamma ]-B\omega ]\mathit{\psi }\prime \times {\mathbf{e}}_3+\gamma \lambda {\mathbf{e}}_3\times {\eta }_{\alpha }{\mathbf{e}}_{\alpha }=0,$$

(95)

$$[g\prime \left(\frac{\lambda }{\gamma }\right)\mathit{\rho }″+A\mathit{\psi }″]·{\mathbf{e}}_3=0,$$

(96)

$$[A\mathit{\rho }″+B\mathit{\psi }″]·{\mathbf{e}}_3=0,$$

(97)

$$[\mathit{\rho }\prime -\lambda \mathit{\psi }\times {\mathbf{e}}_3]·{\mathbf{e}}_{\alpha }=0,$$

(98)

$$\mathit{\psi }(\pm \frac{1}{2})=0,$$

(99)

$$\mathit{\rho }(\pm \frac{1}{2})·{\mathbf{e}}_{\alpha }=0,$$

(100)

$$[g\prime \left(\frac{\lambda }{\gamma }\right)\mathit{\rho }\prime (\pm \frac{1}{2})+A\mathit{\psi }\prime (\pm \frac{1}{2})]·{\mathbf{e}}_3=0,$$

(101)

with details provided in Appendix B. Since $\mathrm{Bg}\prime \left(\frac{\lambda }{\gamma }\right)-A^2$ is non-zero (assumed to be positive), equations (96) and (97) imply

$$\mathit{\rho }″·{\mathbf{e}}_3=0\mathrm{and}\mathit{\psi }″·{\mathbf{e}}_3=0.$$

(102)

Boundary condition (99) requires $\mathit{\psi }(S)\in \mathrm{span}\{{\mathbf{e}}_1,{\mathbf{e}}_2\}$ , which motivates the introduction of the decomposition

$$\mathit{\rho }(S)={\mathit{\rho }}_t(S)+{\mathit{\rho }}_a(S),$$

(103)

where ${\mathit{\rho }}_t(S)\in \mathrm{span}\{{\mathbf{e}}_1,{\mathbf{e}}_2\}$ and ${\mathit{\rho }}_a(S)\in \mathrm{span}\{{\mathbf{e}}_3\}$ .

Equations (94)–(101) can now be reduced to the following (details in Appendix C):

$$\mathit{\psi }″+\zeta \mathit{\psi }\prime \times {\mathbf{e}}_3=\mathit{\psi }\prime (+\frac{1}{2})-\mathit{\psi }\prime (-\frac{1}{2}),$$

(104)

$${\mathit{\rho }\prime }_t=\lambda \mathit{\psi }\times {\mathbf{e}}_3,$$

(105)

$${\mathit{\rho }″}_a=0,$$

(106)

accompanied by the boundary conditions

$${\mathit{\rho }}_t(\pm \frac{1}{2})=0,$$

(107)

$$\mathit{\rho }{\prime }_a(\pm \frac{1}{2})=0.$$

(108)

The new parameter $\zeta$ appearing in equation (104) is defined as

$$\zeta :=A[\lambda -\gamma ]-[B+1]\omega .$$

(109)

It is clear that ${\mathit{\rho }}_a(S)=C_o{\mathbf{e}}_3$ for all S, where $C_o$ is a constant that appears because we have put no physical constraint in the axial direction. As the rod can slide axially without causing any strain, we can fix $C_o=0$ .

For $\zeta =0$ , the problem admits only trivial solutions (Appendix C). Now assuming $\zeta \ne 0$ , the differential equations (104) and (105) admit general solutions of the form

$$\underset{¯}{\psi }(S)=\frac{C_1}{\zeta }\left[\begin{array}{c}\sin (\zeta S)\hfill \\ \multicolumn{1}{c}{-\cos (\zeta S)+2S\sin \frac{\zeta }{2}}\\ \multicolumn{1}{c}{0}\end{array}\right]+\frac{C_2}{\zeta }\left[\begin{array}{c}-\cos (\zeta S)-2S\sin \frac{\zeta }{2}\hfill \\ \multicolumn{1}{c}{-\sin (\zeta S)}\\ \multicolumn{1}{c}{0}\end{array}\right]+\left[\begin{array}{c}{C}_3\hfill \\ \multicolumn{1}{c}{C_4}\\ \multicolumn{1}{c}{0}\end{array}\right],$$

(110)

$$\begin{array}{c}\underset{¯}{{\rho }_t}(S)=C_1\frac{\lambda }{{\zeta }^2}\left[\begin{array}{c}-\sin (\zeta S)+\zeta S^2\sin \frac{\zeta }{2}\hfill \\ \multicolumn{1}{c}{\cos (\zeta S)}\\ \multicolumn{1}{c}{0}\end{array}\right]+C_2\frac{\lambda }{{\zeta }^2}\left[\begin{array}{c}\cos (\zeta S)\hfill \\ \multicolumn{1}{c}{\sin (\zeta S)+\zeta S^2\sin \frac{\zeta }{2}}\\ \multicolumn{1}{c}{0}\end{array}\right]\\ +\lambda \left[\begin{array}{c}{C}_5+C_{4}S\hfill \\ \multicolumn{1}{c}{C_6-C_{3}S}\\ \multicolumn{1}{c}{0}\end{array}\right],\end{array}$$

(111)

where $C_1$ , $C_2$ , ⋯, $C_6$ are generic integration constants in $\mathbb{R}$ . The representations $\underset{¯}{\psi }$ and $\underset{¯}{{\rho }_t}$ are with respect to the fixed basis. The boundary conditions (99) when invoked into equation (110) leads to

$$[C_1-C_2]\sin \frac{\zeta }{2}=0,$$

(112)

simultaneously giving

$$C_3=\frac{C_2}{\zeta }\cos \frac{\zeta }{2},C_4=\frac{C_1}{\zeta }\cos \frac{\zeta }{2}.$$

(113)

The values of $\zeta \ne 0$ for which $\sin \frac{\zeta }{2}=0$ eventually leads to the trivial solution (Appendix C). Therefore, we assume $C_1=C_2$ , which when substituted in the general solution (111) and forced to satisfy equation (107), leads to the condition

$$\frac{1}{\zeta }\sin \frac{\zeta }{2}-\frac{1}{2}\cos \frac{\zeta }{2}=0.$$

(114)

It simultaneously requires

$$C_5=-\frac{C_1}{{\zeta }^2}[\frac{\zeta }{4}\sin \frac{\zeta }{2}+\cos \frac{\zeta }{2}]=C_6.$$

(115)

Hence we have an out-of-plane solution,

$$\begin{array}{c}\underset{¯}{{\rho }_t}(S)=C_1\frac{\lambda }{{\zeta }^2}[\cos (\zeta S)+[S^2-\frac{1}{4}]\zeta \sin \frac{\zeta }{2}-\cos \frac{\zeta }{2}]\left[\begin{array}{c}1\hfill \\ \multicolumn{1}{c}{1}\\ \multicolumn{1}{c}{0}\end{array}\right]\hfill \\ +C_1\frac{\lambda }{{\zeta }^2}[S\zeta \cos \frac{\zeta }{2}-\sin (\zeta S)]\left[\begin{array}{c}1\hfill \\ \multicolumn{1}{c}{-1}\\ \multicolumn{1}{c}{0}\end{array}\right],\hfill \end{array}$$

(116)

whose existence is subject to the condition that parameters $\gamma$ and $\lambda$ admit sensible solutions ( $\gamma >0$ and $\lambda >0$ ). A positive increasing sequence $(a_n{)}_{n=1}^{\infty }$ satisfying $\tan (a_n)=a_n$ can be defined. The values taken by $\zeta \in \{\pm 2a_n:n\in \mathbb{N}\}$ correspond to the discrete bifurcation modes.

4.3. Results and discussion

In view of the equivariance properties of our problem (Appendix D), any rotation of equation (116) about ${\mathbf{e}}_3$ is an acceptable solution. Hence the solution can be simplified to

$$\underset{¯}{\rho }(S)=C_1\frac{\lambda }{{\zeta }^2}\left[\begin{array}{c}\cos (\zeta S)+[S^2-\frac{1}{4}]\zeta \sin \frac{\zeta }{2}-\cos \frac{\zeta }{2}\hfill \\ \multicolumn{1}{c}{\sin (\zeta S)-S\zeta \cos \frac{\zeta }{2}}\\ \multicolumn{1}{c}{0}\end{array}\right],$$

(117)

$$\underset{¯}{\psi }(S)=C_1\frac{\lambda }{\zeta }\left[\begin{array}{c}\cos \frac{\zeta }{2}-\cos (\zeta S)\hfill \\ \multicolumn{1}{c}{2S\sin \frac{\zeta }{2}-\sin (\zeta S)}\\ \multicolumn{1}{c}{0}\end{array}\right],$$

(118)

$${\eta }_1(S)=0,$$

(119)

$${\eta }_2(S)=C_1\frac{\zeta }{\gamma }\cos \frac{\zeta }{2},$$

(120)

where representations (117) and (118) are with respect to the fixed basis. This solution is clearly flip-symmetric about ${\mathbf{e}}_1$ , thus suggesting that equation (116) is also flip-symmetric, but about an axis different from ${\mathbf{e}}_1$ . The deformed centre-line $\mathbf{r}(S)$ for this solution is

$$\underset{¯}{\mathrm{r}}(S)=\frac{\lambda }{{\zeta }^2}\left[\begin{array}{c}\cos (\zeta S)+[S^2-\frac{1}{4}]\zeta \sin \frac{\zeta }{2}-\cos \frac{\zeta }{2}\hfill \\ \multicolumn{1}{c}{\sin (\zeta S)-S\zeta \cos \frac{\zeta }{2}}\\ \multicolumn{1}{c}{{\zeta }^{2}S}\end{array}\right],$$

(121)

represented with respect to the fixed basis, wherein $\epsilon C_1=1$ is set for the sake of simplicity. For a particular $\zeta \in \{\pm 2a_n:n\in \mathbb{N}\}$ , the end-to-end distance $\lambda$ and growth stage $\gamma$ can be found by solving the system

$$g\left(\frac{\lambda }{\gamma }\right)=\frac{A}{\mathcal{M}}[1-\frac{1}{\gamma }],$$

(122)

$$\zeta =A[\lambda -\gamma ]-\frac{B+1}{\mathcal{M}}[\gamma -1],$$

(123)

simultaneously (Table 1, Appendix E). Equations (122)–(123) couple the axial force response of the rod with the bifurcation mode caused due to growth, via the kinematic constraint of symmetry-preserving growth (equation (65)). Whenever this system does not admit a solution $\gamma >0$ and $\lambda >0$ , the perturbation chosen gives only trivial solutions, indicating that out-of-plane buckling is not guaranteed. Note that, linear stability analysis reveals only the shape of the buckled state, without any information on its amplitude. In the ensuing discussion, all results are described in terms of $\gamma$ , $\lambda$ and $\zeta$ – these are the fundamental kinematic quantities that help us understand the interplay between growth and material chirality.

An inspection of equation (121) reveals that the sign change $\zeta \mapsto -\zeta$ reverses the chirality of the solution curve, reflecting it about the ${\mathbf{e}}_1$ – ${\mathbf{e}}_3$ plane (Figure 5). Moreover, since our solution is flip-symmetric, this is equivalent to the reflection in the ${\mathbf{e}}_1$ – ${\mathbf{e}}_2$ plane. Reflections in the ${\mathbf{e}}_2$ – ${\mathbf{e}}_3$ plane also give solutions of opposite chirality, but need to be rotated by ${180}^{°}$ about the ${\mathbf{e}}_3$ -axis to coincide with the reflections in ${\mathbf{e}}_1$ – ${\mathbf{e}}_3$ plane. These centre-line solutions with handedness are similar to those obtained by Healey and Papadopoulos [32] for a fixed–fixed rod under axial compression.

Figure 5.

Out-of-plane bifurcated solution for the case $\mathcal{M}=-0.1$ , $A=-8$ , $B=1.2$ and $F={10}^5$ . The spatial curves shown correspond to the shape of the rod’s centre-line in the first three buckled modes. These solutions are flip-symmetric about ${\mathbf{e}}_1$ (Appendix D). The values of end-to-end distance $\lambda$ and growth stage $\gamma$ at buckling are calculated using equations (122)–(123). The magnitude of $\zeta$ denotes the mode of bifurcated solution, while its sign controls the chirality. Reflections about the ${\mathbf{e}}_1$ – ${\mathbf{e}}_3$ plane correspond to the reversal in chirality $\zeta \to -\zeta$ .

Internal chirality of the rod is taken care of by the constants $\mathcal{M}$ and A. In the case of hemitropic rods, A captures chirality in load response of the rod while $\mathcal{M}$ contains information regarding the chiral growth law. For rods with n-fold helical symmetry, A denotes the same thing, but with the assumption of symmetry-preserving growth in place, $\mathcal{M}$ captures chirality in microstructure.

Consider two rods with opposite internal chirality with all other material properties the same. Let one of them with chiral constants A, $\mathcal{M}$ have a solution with bifurcation mode $\zeta$ , end-to-end distance $\lambda$ and growth stage $\gamma$ . Naturally the second rod with opposite internal chirality is expected to give rise to a reflected solution with bifurcation mode $-\zeta$ , while end-to-end distance and growth stage are still the same. Thus, equations (122) and (123) imply that the chiral constants associated with the second rod are $-A$ and $-\mathcal{M}$ . We infer that the complete reversal of internal chirality in rods requires the transformations $\mathcal{M}\mapsto -\mathcal{M}$ and $A\mapsto -A$ to be taken simultaneously. In addition, the $\zeta$ solution of a rod with internal chirality $\mathcal{M},A$ and the $-\zeta$ solution of a rod with opposite internal chirality $-\mathcal{M},-A$ are mirror images with respect to ${\mathbf{e}}_1$ – ${\mathbf{e}}_2$ and ${\mathbf{e}}_1$ – ${\mathbf{e}}_3$ planes.

In the absence of guides, the end-to-end distance of the rod would have been the same as its growth stage $\gamma$ . Therefore, in order to understand the influence of guides on this end-to-end distance, we compare the growth stage $\gamma$ at which the rod buckles with the corresponding value of $\lambda$ . Assuming that equations (122) and (123) admit an acceptable solution, the monotonicity of $g(·)$ and the condition $g(1)=0$ reveal the following observations:

Growth $\gamma >1$

• A and $\mathcal{M}$ are of the same sign if and only if $\lambda >\gamma$ , signifying that the ends in the current configuration have moved away from each other, as compared to both initial and virtual configurations.

• A and $\mathcal{M}$ are of opposite sign if and only if $\lambda <\gamma$ , signifying that the ends in the current configuration have come closer as compared to the virtual configuration, but no guaranteed comparison can be made with the initial configuration.

Atrophy $\gamma <1$

• A and $\mathcal{M}$ are of opposite sign if and only if $\lambda >\gamma$ , signifying that the ends in the current configuration have moved apart as compared to the virtual configuration, but no guaranteed comparison can be made with the initial configuration.

• A and $\mathcal{M}$ are of the same sign if and only if $\lambda <\gamma$ , signifying that the ends in the current configuration have come closer as compared to both initial and virtual configurations.

For a rod with n-fold helical symmetry with the growth law assumed to be symmetry-preserving, these results reveal an interesting interplay between chiralities in microstructure and load response of the rod. But for a hemitropic rod, the growth law allowing the cross-section to rotate makes the guided–guided problem similar to a non-growing rod subject to a axial twist at one end while the other end is free to move axially. The results above directly reflect the twist–extension type Poisson effect expected in hemitropic rods.

Case of isotropy A = 0. In this case, the solutions have $n_3=0$ with

$$\gamma =\lambda =1-\frac{\zeta \mathcal{M}}{B+1}.$$

(124)

A growing isotropic rod has an out-of-plane solution with sign of $\zeta$ opposite to that of $\mathcal{M}$ . But for a decaying isotropic rod, equation (124) guarantees an out-of-plane solution only if $|\mathcal{M}|<\frac{B+1}{2a_1}$ , and hence such solutions exist only up to the first few modes (for the chosen perturbation), with sign of $\zeta$ the same as that of $\mathcal{M}$ .

For small $A\ne 0$ , the solution is close (in terms of $\gamma$ and $\lambda$ ) to that of the isotropic case with B, $g(·)$ , $\mathcal{M}$ and $\zeta$ kept the same. In addition, the chirality of these solutions is the same as that of the corresponding isotropic case. With $A\ne 0$ , growing rods admit $\zeta \mathcal{M}>0$ and atrophying rods admit $\zeta \mathcal{M}<0$ only if A is taken to be very large, which in turn may be unrealistic.

Consider two rods with degrees of hemitropy $A^{+}>0$ and $A^{-}<0$ , such that $A^{+}+A^{-}=0$ , everything else being kept the same. One of these cases gives a solution where ends come closer, while the ends move apart in the other case (comparisons made here are with respect to the virtual configuration). Let ${\lambda }^{+}$ and ${\lambda }^{-}$ denote the respective solutions for $A^{+}$ and $A^{-}$ , whereas ${\lambda }^o$ denotes the same for the isotropic case. While ${\lambda }^o$ may lie between ${\lambda }^{+}$ and ${\lambda }^{-}$ , it is also a possibility that both ${\lambda }^{+}$ and ${\lambda }^{-}$ might lie on the same side of ${\lambda }^o$ (Figure 6), thus suggesting that no definitive comment can be made on this.

Figure 6.

Variation of end-to-end distance $\lambda$ and growth stage $\gamma$ at buckling with the degree of hemitropy A for first mode solutions ( $\zeta =8.986$ ).

5. Conclusion

In this work we study the growth of slender elastic rods with chiral material symmetries – transverse hemitropy and multi-fold dihedral helical symmetry. Based on the intuitive notion that rods with helical symmetry should twist during growth, we propose a homogeneous growth law that allows for relative rotation of cross-sections. A guided–guided rod setup is considered to illustrate the occurrence of out-of-plane buckling at certain stages of growth (or atrophy). These solutions obtained are flip-symmetric and chiral in nature. A complete mirroring of the rod, including both growth and constitutive properties, gives a solution with opposite chirality, under the same deformation. We show that the end-to-end distance at bifurcation modes for the isotropic case need not lie between those for rods of opposite material chiralities, with the rest of the elastic and growth properties kept the same. The end-to-end distances for different combinations of growth (atrophy) and material chiralities have also been examined to understand the effect of twisting growth on the constitutive twist–extension coupling.

Embedding our biologically active (growth or atrophy) chiral rod setup in an elastomeric matrix and introducing inhomogeneities similar to that by Almet et al. [38] could be an interesting direction to explore. One can also consider a ply of biologically active rods, like growing bi-rods [24], to study the effect of growth and material chiralities of individual rods on the total deformation.