A Simulation of an Elastic Filament Using Kirchhoff Model

This paper presents numerical simulations and comparisons between different approaches concerning elastic thin rods. Elastic rods are ideal for modeling the stretching, bending, and twisting deformations of such long and thin elastic materials. The static solution of Kirchhoff's equations [2] is produced using ODE45 solver where Kirchhoff and reference system equations are combined instantaneously. Solutions using formulations are based on Euler’s elastica theory [1] which determines the deformed centerline of the rod by solving a boundary-value problem, on the Discreet Elastic Rod method using Bishop frame (DER) [5,6] which is based on discrete differential geometry, it starts with a discrete energy formulation and obtains the forces and equations of motion by taking the derivative of energies. Instead of discretizing smooth equations, DER solves discrete equations and obeys geometrical exactness. Using DER we measure torsion as the difference of angles between the material and the Bishop frame of the rod so that no additional degree of freedom is needed to represent the torsional behavior. We found excellent agreement between our Kirchhoff-based solution and numerical results obtained by the other methods. In our numerical results, we include simulation of the rod under the action of the terminal moment and illustrations of the gravity effects.


Introduction
Elastic rods have many interesting applications, at large length scale they are used to study the mechanics of macro structural elements such as marine cables. At small length scales, these theories have been used to model biological threads, hair beams [15], DNA molecules [16] rubber bands and microfibers [17]. They are ideal for modeling the stretching, bending, and twisting deformations of such long, thin elastic materials. In a seminal work [2] Gustav Kirchhoff proposed a new theory in order to model the bending and torsion of a rod. In the 20th century, the Cosserat brothers introduced the directory theory which was an extension of Kirchhoff's rod theory in [3] and [4]. The rod theory is also considered as an example of a Cosserat rod theory. As far as we are aware that the analytic solutions of Kirchhoff equations for static rods are well studied. In [7] the most general analytical solution for a static Kirchhoff model is presented. In [8,9] a time-dependent perturbation scheme is employed to study the stability of equilibrium solutions of the Kirchhoff rod model. This technique has been extensively used, to study DNA rings [10,11]. Despite the success of the above examples, many contributions toward considering equilibrium solutions of the Kirchhoff model subject to some particular boundary conditions have been made, see [12][13][14] and references therein. This list is by no means exhaustive, and research in this direction is actively ongoing. The present contribution might be helpful in this endeavor.
The paper is organized as follows. In the next section of this paper, we present a short introduction to the general Kirchhoff rod model, we consider separately the geometric description of the rod (subsection 2.1), its deformation geometry (subsection 2.2), and lastly, we derive the equilibrium equations for a Kirchhoff rod. In Section 3 we describe our implementation procedure and present some numerical examples, including simulation of the rod under the action of the terminal moment, illustrations of the gravity force, we also consider examples of helices with different orientations, the first one is clamped vertically, while the other one is a horizontal helix. Finally, Section 4 contains our main results.

Kirchhoff Rod Model
In this study, the composition of the rod is assumed to be isotropic and linear elastic. According to the physics of Kirchhoff, the rod Γ is a deformable body in which the length, dominates the other two dimensions of the cross-section. We divide the rod into thin disks of length ds and cross-section ( ) A s . We apply the conservation law of the linear and angular momentum to every disk. We will study a segment of a rod with length ds . A segment of the rod can be evaluated by the arc-length parameter [0, ] s L ∈ . Let us denote with A the cross section of the rod.

Geometric Description
The configuration of the rod as shown in figure 1 is a parameterized space curve ( , ) s t r along with a parameterized orthonormal basis 3 ( , ), ( , ), ( , ). s t s t s t  ( , ) s t r and its derivative.

Geometry of Deformation
where u is Darboux vector.

Mechanics of the Rod
We now consider the balance laws that allow us to determine which configurations are equilibrium shapes.
We suppose the forces acting on the cross section [ ] The equation (4) is true for Since n + is a continuous quantity, from this condition We drop the superscripts Using some elementary calculation, we have We also assume a resultant contact torque

t a t s t s t ds s t A s s t ds t
Using elementary calculations leads to: Equations (5) and (7) are the equilibrium form for a special Cosserat rod. Substituting the equation (5) into the equation (7), we have  Equation (10), yields Substituting the equation (9) into the equation (10) Before we compare our method with the other methods, let us briefly remind the idea behind Discreet Elastic Rod method and Euler's classic theory. The discrete elastic rod method (DER) proposed by Bergou et al. [5,6], based on discrete differential geometry, starts with a discrete energy formulation and obtains the forces and equations of motion by taking the derivative of energies. Instead of discretizing smooth equations, DER solves discrete equations and obeys geometrical exactness. While the Euler's elastica theory determines the deformed centerline of the rod by solving a boundary-value problem. In the next section, we implement and perform simulations for different scenarios.

Implementation and Numerical Results
Let us consider the balance laws for a special Cosserat rod

t s t A s s t x x t s s s s s t s t s t s s s
The partial differential equations (11) enhanced by constitutive relations for the moment m we have Where 1 2 , EI EI are flexural rigidity, D is a torsional stiffness and E is Young's Modulus.
For a static rod, all derivative with respect to t vanish, that is ( ) ( )   we obtain a solution of the equation (14). But this solution is full understood if it is written with respect to the external and fixed system. It is known that there exists a transformation which takes vectors from the fixed system to the local one and vice versa (see figure 4). We write this transformation as follows where i e is a base vector of the fixed system and ( , ) Q s t is proper rotation matrix of the form 11 If we substitute equation (17) and (15) to (16), we obtain equation (2) Figure 6 shows a numerical simulation of the rod under the action of terminal moment. In Figure 7, we consider two helices with different orientations under the action of gravity, the first one is clamped vertically, while the other one is a horizontal helix.

Conclusion
We concluded this work with some numerical experiments. In this paper, a static solution for Kirchhoff's equations is obtained. We have integrated the static equations as an initial value problem, to integrate our equations the ODE45 scheme is used, we compare our solution with classical Euler theory of the Elastica and Discreet Elastic Rod model and see an excellent agreement. Using the proposed technique, solutions under the action of external forces can be easily obtained. In the present paper, we have only considered the gravity force as an external force, but we believe the technique presented here can be used in other classes of static solutions for Kirchhoff's equations, we can mention here, filament in a rotating system, where centrifugal acceleration plays the role of the external force and the electric force for a charged filament as well.