# Helical post-buckling of a rod in a cylinder: with applications to drill-strings

J. M. T. Thompson, M. Silveira, G. H. M. van der Heijden, M. Wiercigroch

## Abstract

The helical buckling and post-buckling of an elastic rod within a cylindrical casing arises in many disciplines, but is particularly important in the petroleum industry. Here, a drill-string, subjected to an end twisting moment combined with axial tension or compression, is particularly prone to buckling within its bore-hole—with potentially serious results. In this paper, we make a detailed theoretical study of this type of instability, deriving precise new results for the advanced post-buckling stage when the rod is in continuous contact with the cylinder. Results, including rigorous stability analyses and contact pressure assessments, are presented as equilibrium surfaces to facilitate comparisons with experimental results. Two approximate solutions give insight, universal graphs and parameters, for the practically relevant case of small angles, and highlight the existence of a critical cylinder diameter. Excellent agreement with experiments is achieved.

## 1. Introduction

The helical buckling of an elastic rod in contact with a surrounding cylindrical casing of larger diameter arises in several areas of science and engineering. It is of major interest to petroleum engineers, whose long thin drill-strings are simultaneously twisted and axially loaded. Near the surface the axial load is tensile (to reduce the self-weight load on the drill-bit) but lower down (and significantly in horizontal drilling) the axial load becomes compressive. We review the drill-string literature in §6. A second area of interest arises in the deployment of a stent when treating patients with heart disease. Here, the guide-wire, of length about 1 m and diameter of perhaps 0.4 mm, is inserted into a human artery whose internal diameter is ca 4 mm. During its use, the surgeon relies on the longitudinal and twisting movement on the input end to control the leading-end movement (Schneider 2003). Helical buckling and ply formation also play a key role in the functioning of the DNA molecule (Thompson et al. 2002; Calladine et al. 2004; Travers & Thompson 2004; Thompson 2008).

Quite apart from applications, the helix in a cylinder has become an archetypal problem in advanced continuum mechanics, because it exhibits (perhaps in their simplest form) complex nonlinear phenomena of spatial localization and chaos. This is presented in the papers by van der Heijden and co-workers (van der Heijden 2001; van der Heijden et al. 2002), which derived from earlier papers on complexity in rod mechanics (Champneys & Thompson 1996; Champneys et al. 1997; van der Heijden et al. 1998, 1999). In a similar vein, a recent paper by Chen & Li (2011) studies in fine detail a thin elastic rod constrained within a cylinder under the action of an end twisting moment. Akin to earlier work by Coleman & Swigon (2000) on ply formation, it is found that, with an increase of the moment, the rod progresses from one-point contact to two- and then three-point contact, and finally to line contact. The authors are able to get these results by making the assumption of small deformations (small slopes). This renders some of the equations linear, which helps in formulating a shooting problem suitable for numerical solution. The authors also perform experiments on a metal alloy rod inside a glass tube and find reasonable agreement with the theoretical predictions.

In the light of all this interest and relevance, we make here a careful theoretical study of the deformations and stability of a rod in continuous contact with a cylinder on the inside (with what we would call positive pressure) or the outside (which in our formulation would correspond to negative pressure). We study the rod under end twisting moment and an axial load that can be either tensile or compressive, and display the results in graphical forms that give a high level of understanding and applicability.

## 2. Theoretical approach

This paper is devoted to the helical buckling of a long elastic rod of length S within a wider cylindrical casing, owing to an applied tension T and an end twisting moment M. The rod is taken to have circular cross section, either solid or hollow. The bending stiffness about any axis is denoted by B. Note that T can be negative, implying a compressive load P=−T. Lying initially along the central axis of the cylinder, buckling and post-buckling of the free rod can eventually generate contact with the casing. We use a Cosserat energy formulation for the rod, and, restricting attention to uniform helical states (whether free or in contact with the casing), the problem can be completely cast in terms of three degrees of freedom.

The corresponding deflection of the tension, T, is the extension written as Sd, and that of the moment, M, is the end rotation Sr, and our aim is to examine in detail the equilibrium surfaces T(d,r) and M(d,r). In particular, we relate features of the surfaces to the stability changes that are encountered under dead, rigid and mixed loading conditions (Thompson 1961, 1979). We use curly brackets to indicate which two variables are controlled, the remaining two being then passive responses. The four possible cases are then totally dead loading {T,M}, totally rigid loading {d,r} and two forms of mixed loading {T,r} and {d,M}.

Defining the dimensionless load parameter, , the free rod (under tension) exhibits a subcritical bifurcation of helical post-buckling states at mc=2, from which the rod would jump dynamically under {T,M}—where it would eventually stabilize depends crucially on the dynamics and the gap size. Meanwhile under rigid loading, {d,r}, the rod is known to follow a spatially localizing path which also emanates from the bifurcation point and is energetically favourable (Thompson & Champneys 1996). This path eventually reaches a fold, from which a free rod would jump into a self-contacting loop or ply. The gap between the rod and casing will govern whether this jump is observed or not. A simple {d,r} experiment under fixed r and varying d in which there is no jump is illustrated in figure 1. The buckle amplitude grows steady under decreasing d, and a contacting helix gradually extends along the tube. This steady growth of a buckling deflection can be attributed to slight misalignment of the rod in the jaws of the drill-chucks that grip the loaded ends.

Figure 1.

Some basic experiments under rigid {d,r} loading, performed with the help of Charlie McRobie (Comberton Village College). Here all units are arbitrary.

In this paper, we analyse first the helical deformations of the free rod, by setting to zero derivatives of the total potential energy with respect to the internal torsion, τi, the helical angle, β, and the helical radius, ρ. When we turn to analyse the rod in contact with its casing the last condition is replaced by the constraint ρ=ρw, where ρw is the effective radius of the casing. Meanwhile, the sign of the ρ derivative informs us about the pressure between the rod and casing.

The situation is illustrated schematically in figure 2. This shows the variation of the total potential energy, V , of the free system at three values of m, two less than mc and one greater than mc. This energy is really a function of τi, β and ρ, but in this schematic diagram the variation with τi and β is ignored.

Figure 2.

Schematic of the form taken by the free total potential energy at three values of m. Free equilibrium solutions are illustrated by large balls, contacting solutions by the smaller balls.

Concerning the range and scope of the study, we should note that T can be positive or negative, a negative value corresponding to compression. Although we mainly concentrate on positive T, the compression results are perfectly valid, though overall Euler buckling is impossible because of our assumed helical form. We also focus on positive M, with negative M simply giving mirror-image results with the helical angle reversed. We allow the helical angle, β, to range from −π to +π. When β is positive and equal to +π/2 the rod forms a ring, and mathematically (though not physically owing to the inevitable self-contacting of the rod) this ring can pass through itself to give a straight rod again at β=π. Negative β, generated by negative M, changes the handedness of the helix, and a mirror-image sequence can be generated as β decreases to −π. Finally, we can contemplate either positive or negative pressure between the rod and its casing. Negative pressure could arise from magnetic or electrical attraction/repulsion, or more simply if the rod were located outside a cylinder.

The paper is organized as follows. Section 3 has the main definitions for internal and kinematic torsion, following the analysis of Love (1927). Some simple geometrical formulae are defined for the total angle turned by the end moment, arc-wavelength, spatial frequency and end shortening. The large amplitude equilibrium equations are derived from the total potential energy of the rod, following Thompson & Champneys (1996). Section 4 contains the results for the free rod before contact with the wall, leading to the (exact) equation of the post-buckling equilibrium path. Section 5 considers the rod in contact with the wall, with the helical radius fixed and equal to that of the constraining cylinder. An exact solution is derived which allows us to discriminate between solutions with positive and negative pressure on the casing, and the stability of the post-buckling solutions is assessed exactly, again with the help of M and T contours. Results for the rod under pure compressive load show that the origin is a singular point and that two real folds exist, implying a jump with increasing β. In addition, a first-order small angle approximation is pursued, noting that the solution passes around a fold bifurcation. A second approximation is also used which permits us to characterize wide and narrow bore-hole radii. Section 6 contains some qualitative results from a small scale experiment in which the rod dynamically assumes a helical configuration under stick–slip oscillations. Also, the theoretical results for pure compression are compared with the experimental work of Wu et al. (1993) with excellent agreement. Finally, conclusions are drawn in §7.

## 3. Cosserat formulation

### (a) The rod and its casing

Consider the buckling and post-buckling of a long rotationally symmetric elastic rod carrying an end tension, T, and an end twisting moment, M. It lies within a wider cylindrical casing as illustrated in figure 3.

Figure 3.

The cross section of the rod within its cylindrical casing. Here, for simplicity, the cross section of the rod is drawn as circular; strictly it will be elliptical owing to the helical deformation.

The rod is assumed to be axially inextensional, and the only elastic constants entering into the analysis are the bending stiffness, B, and the torsional stiffness, C. We assume that the rod deforms into a uniform helix, and therefore has three degrees of freedom represented by the internal torsion, τi, which is a generalized strain defined in §3b, the helical angle, β, and the helical radius, ρ. The rod in straight and helical configurations (without the casing) is illustrated in figure 4. In terms of a scaled total potential energy, v, the three equilibrium equations for the free rod are ∂v/∂τi=0, ∂v/∂β=0 and ∂v/∂ρ=0.

Figure 4.

The straight elastic rod and the helical deformation.

The rod of cross-sectional radius a lies at first centrally within a wider casing of radius A>a with which it can make continuous frictionless contact once the radius of the helix is large enough. Since, mathematically, the Cosserat rod has no cross-sectional dimension, the appropriate criterion for contact will be that ρ=ρw=Aa. Note that in offshore drilling the ratio of the bore-hole diameter to the drill-string diameter is usually about 3.5. Once contact is made, the last of the three equilibrium equations becomes invalid and is replaced by the constraint ρ=ρw.

For the free, non-contacting rod (under tension) the assumed helix is indeed the exact solution to the linearized equations for the buckling of the long (strictly infinite) straight trivial solution, and the analysis correspondingly predicts the correct buckling condition . Notice in passing that this result does not contain C, which does indeed play a rather insignificant role in twisted rod theory. As with the buckling condition, our three-mode theory provides an exact helical solution both before and during the contact with the wall. Notice finally that, as T tends to zero, the critical value of the moment, Mc, tends to zero. It is the tension that stabilizes the long free rod. There are, of course, no compressive solutions for the free rod in our analysis, since the Euler buckling solutions (known to be unique) violate our helical assumption.

### (b) Internal and kinematic torsion

Before starting our analysis of a Cosserat rod (namely, a mathematical line imbued with appropriate stiffnesses, but no cross-sectional area), we need to take note of Love's analysis of the kinematic twist rate, here written as torsion (Love 1927). The total twist per unit arc-length, τ, of the rod is decomposed into two parts, called the internal torsion, τi, and the kinematic torsion of the space curve defined by the rod, τs, as follows: 3.1

Note that this use of τ for a kinematic twist rate (a generalized strain) follows a long geometrical tradition; since we shall never be talking about stresses, there will be no cause for confusion with its use in continuum elasticity as a shear stress. For our linearly elastic rod with torsional stiffness, C, the twisting moment at any point is given by .

For the helix illustrated in figure 4, with pitch angle β and radius ρ, the geometrical theory of space curves gives the kinematic torsion, τs, and the total curvature, κ, as 3.2 and 3.3

The total angle in radians, R, turned by the end moment is given, using the topological concepts of twist, link and writhe as in van der Heijden & Thompson (2000), by 3.4 where S is the length of the (inextensional) rod. This expression is adequate for the present investigation, being strictly correct for a long uniform helix that is supported in a manner similar to that envisaged by Love (1927). Notice that r is not given simply by τ. If it were, we would find that under a prescribed end rotation, r, we would have a prescribed torsion, τ, and a buckling deformation into a helix would achieve no reduction in torsional strain energy (so there would be no buckling).

We must finally write down some simple geometrical formulae that we shall need in the following sections for the arc-wavelength, L, a non-dimensional spatial frequency parameter, σ, and the total axial end shortening of the deformed rod, D, as 3.5 3.6 and 3.7

Here, the word ‘arc’ emphasizes that we are defining these quantities in terms of arc-length along the curved rod, not along the straight central axis. Note that the dimension in figure 4 is commonly called the pitch of the helix.

### (c) Energy formulation

We aim to derive, following Thompson & Champneys (1996), the large amplitude equilibrium conditions for a rod deforming into a uniform helix which can make continuous frictionless contact with its outer casing. In particular, we want to determine the equilibrium value of the helical angle β (and hence the wavelength L) in terms of the tension in the rod, T, and the twisting moment, M.

For this analysis (but not necessarily in its physical interpretation), we can regard T and M as prescribed dead (generalized) loads which do work through the corresponding deflections −D and R. Ignoring for the time being the contact with the casing, the total potential energy, V , can then be treated as a function of the three variables (τi,β,ρ).

The potential energies of T and M are VT and VM, respectively, given by 3.8 and 3.9

The strain energy of bending, VB, and the strain energy of torsion, VC, are given by 3.10 and 3.11

The total potential energy of the system, rod plus dead loads, is now 3.12 and the necessary and sufficient conditions for equilibrium of the free non-contacting rod are that ∂v/∂τi=0, ∂v/∂β=0 and ∂v/∂ρ=0.

The first equation gives us immediately 3.13

Whatever the (unspecified) boundary conditions are, equation (3.13) is a natural overall requirement, demanding that the end of the rod carries about its own body axis the twisting moment M. The second equation gives, using (3.13), 3.14 and the third equation gives, using (3.13), 3.15

Here, we have included the full expression for ∂v/∂ρ, which we shall need subsequently to characterize the contact pressure between the rod and the casing (where the derivative is not equal to zero). Equations (3.14) and (3.15) are invariant under a simultaneous sign change of M and β, which is in line with the handedness discussion of §2.

## 4. The free rod before contact

In this section, we look briefly at the response of the free rod. Note that our analysis will supply no non-trivial solutions under compression because it is restricted to helical deformations. Although we prescribed dead T and M for our energy analysis, the equilibrium conditions that we have established do of course hold for either rigid or dead loading. Having used (3.13) to eliminate τ, the remaining equations (3.14) and (3.15) relate the loads (T,M) to the generalized helical coordinates (β,ρ) for the free rod without wall contact. Note that these remaining full nonlinear equilibrium equations do not contain the torsional stiffness C or the generalized coordinate τi. With a little manipulation they can be written as 4.1 and 4.2 confirming that T≥0.

These two equations, (4.1) and (4.2), give T and M in terms of β and ρ, and can be inverted to give β and ρ in terms of T and M.

Defining the dimensionless load parameter, 4.3 we eliminate ρ between (4.1) and (4.2) to give the equation of the post-buckling equilibrium path as , which simplifies to 4.4 and gives us the buckling load as mc=2. Equation (4.4) is an exact answer for one of the subcritical post-buckling paths of a long rod, because the assumed helix is indeed one of the buckling and post-buckling modes. For m>2, when we are beyond the subcritical buckling bifurcation, the only equilibrium solution is the unstable trivial state with β=0. Note, however, that there is also a spatially localized post-buckling path that is energetically favourable, and is always observed experimentally (Thompson & Champneys 1996). Substituting (4.4) back into (4.1) to eliminate β, we get the equation for ρ as 4.5

Results of this exact free-rod analysis are plotted in figure 5. We notice that reaches a maximum value of unity at m=1. This will clearly influence whether or not the rod makes contact with its casing.

Figure 5.

Universal plots of β and ρ versus the load parameter, m, for the free rod, using (4.4) and (4.5).

Notice that using (3.5) and (4.1), we can immediately show that under all free equilibrium conditions we have the universal result that 4.6 showing that this arc-wavelength depends only on the load T, and that L tends to infinity as T tends to zero. Had we chosen the axial wavelength, no such neat result would have been obtained. The second form of (4.6) makes an instructive comparison with planar Euler buckling because the right-hand side is the compressive critical load for a strut clamped at both ends. This can be compared, physically, with a single wave of the helix of axial length under the tension T. A useful non-dimensional parameter related to L, giving a measure of the spatial ‘arc-frequency’, is σ defined in (3.6), which we now find is given, in all free states, by 4.7

In addition to the final equilibrium equations, it is useful to write down the exact (nonlinear) results for the corresponding deflections of T and M. Equation (3.7) gives us immediately 4.8 and, after a little algebra, (3.4) gives us 4.9

Although eliminated from our equations at an early stage, we must remember our third generalized coordinate, τi, which if needed can be obtained using .

## 5. The rod in contact with the wall

### (a) Exact analysis

We aim now to make an exact analysis of the rod when it is in contact with its casing, so ρ is no longer an unknown and can be equated to ρw, which is the known effective radius of the casing. In particular, we want to determine T and M contours in the (d,r) space. Note that anywhere in this exact analysis we can set P=−T to generate results for the helix under axial compression, P.

For the displacements, we have d from (3.7) and r from (3.4) and (3.13) as 5.1 and 5.2 where ν is Poisson's ratio introduced to represent the ratio of B to C via the relationship B/C=1+ν, which is true for a solid circular elastic rod. When necessary, we take Poisson's ratio as one-third (typical for a metal), giving B/C=4/3.

It follows that for M contours we can immediately write down 5.3 where the + sign corresponds to solutions with β>0 and the − sign to solutions with β<0. For T contours, we first need to solve for M in terms of T and β using (3.14), after which we obtain 5.4

Looking ahead to the ninth figure, we note that there is an asymptote at d=3/2, separating small-d from large-d solutions. It means that at constant T no amount of applied moment M will make a solution path go past . There is a critical T contour (green) at the value of 2w/B=3/4, for which the numerator and denominator in (5.4) are simultaneously zero. This separates T contours that asymptote with from those that asymptote with , as we see in the ninth figure. We should remark here that solutions with β>π/2 are physically meaningful, but cannot be reached in an experiment that starts from β=0, since a real rod cannot pass through (or indeed near to) the ring solution at β=π/2.

Now we shall want to plot equilibrium paths in the form of a response, β (say), against a varying control, M (say), at constant T (as in figure 6). Alternatively, we might plot β against the varying control, T, at constant M (see, for example, the plot of P=−T against β in the tenth figure). These graphs are often called bifurcation diagrams, and on them we shall be particularly interested in fold bifurcations, where the varying control exhibits an extremum.

Figure 6.

(a) Bifurcation diagrams for fixed T (2w/B=−0.4 (dark green) to 1.2 (brown), in steps of 0.2, from left to right at (say) β=1) and (b) enlargement of the top right-hand quadrant with the bifurcation diagram for 2w/B=0.2.

A little thought shows that these fold points will be identified on our (d,r) contour diagrams as points of tangency, where a T contour is tangent to an M contour. There are two ways of confirming and illustrating this observation. Firstly, if we were (say) holding M constant and moving along an M contour until we encountered a tangency with a T contour, this would clearly signal an extremum (maximum or minimum) of T. So the constant-M equilibrium path being followed would display a fold on a graph of T versus a displacement measure such as β. Secondly, and more directly, we can argue that at a tangency there is in a linear sense an ‘adjacent position of equilibrium’, again implying a fold.

For a tangency of the M and T contours, and hence a fold in the M=const. or T=const. bifurcation diagram, we must have 5.5

Using (5.3), (5.4) and (3.14) to eliminate m and T, this yields the tangency curve given by 5.6

A double fold (or cusp) occurs where this tangency curve is in turn tangent to both the M and T contours (if it is tangent to the one it is also tangent to the other), i.e. where 5.7

The second equation (5.7) reduces to a quintic in d with solution d=0.172194, giving β=0.595612 or 34.13° (independent of ν). The first equation (5.7) then gives the critical value w/B=1.068666 and (5.4) gives the critical value . For loads larger than these critical values, there are no folds and hence no stable solutions at small β. Both these critical values are universal, i.e. independent of ν. The tangency curve has an asymptote at . It may also be verified that at d=3/2 the tangency curve has , which agrees with (5.4) after setting 2w/B=3/4, meaning that the tangency curve goes exactly through the critical saddle at d=3/2. These results can be seen later in the ninth figure (tangency curve in black).

We must next consider the pressure between the helical rod and the casing, and a measure of this pressure will be p:=−(ρ3w/B)∂v/∂ρ|ρ=ρw, as illustrated schematically in figure 2. For simplicity, we will in future just call measure p the pressure.

By (3.15), a zero-pressure solution (which must obviously be a solution of the free unconstrained rod) requires 5.8 (there is also the trivial solution β=0, for any M). This gives the curve 5.9

Note that for non-zero solutions β∈(−π,π) a sign change of the pressure requires M and β to have the same sign. For β>0, the pressure is positive below the (+) curve, while for β<0 it is positive above the (−) curve, where the sign refers to the sign in (5.9). The zero-pressure curve (orange in the ninth figure) intersects the tangency curve (black in the ninth figure) just before the latter's maximum (for β>0) or minimum (for β<0), implying that there is always a sign change of the pressure between the two (small-β) folds in either M=const. or T=const. bifurcation diagrams (confirmed in figure 6).

As we shall examine in the tenth figure, the T bifurcation diagram for purely compressive buckling (M=0) exhibits a jump instability at a fold. The corresponding critical load can be obtained from (3.14). Setting M=0 gives so that the condition for a fold, dT/=0, yields (54.74°). The critical compression is then found as .

In the T bifurcation diagram, there is a special point at β=2π/3 through which all M contours pass. It corresponds to the branch d=3/2 of the critical contour at 2w/B=3/4 in the ninth figure. All M contours intersect this branch.

Finally, to determine the stability characteristics of the equilibrium paths we compute the Hessian matrix, H, to check the second variation of the total potential energy 5.10 evaluated at a solution. Stability requires this matrix to be positive definite, which is the case only if both the trace Tr(H) and the determinant det(H) are positive. We find that stability changes at folds, as expected, and the stability regimes are indicated in our final figures. Note that our notion of stability is restricted within the space of helical solutions.

To examine the results of the exact analysis, we look first in figure 6 at some equilibrium paths in the bifurcation diagram. For a set of fixed T values, figure 6a shows the colour-coded equilibrium paths as the helical angle β against the twisting moment M. All equilibrium solutions in the β range of −π to +π are shown. The regimes in which the contact pressure is positive are shown in grey. These regimes are bounded by the p=0 curve in black. We remember that where the pressure is zero we are at a post-buckling solution of the free rod, with inevitably m<2. Figure 6b focuses on the red curves for 2w/B=0.2, and on the lower branch we see two folds bounding a stable region with a sign change of pressure between the folds at β=0.4636, w/B=0.8472, which correspond to β=27° and m=1.8944.

Turning to the contour diagrams, figure 7 shows T and M contours in the (d,r) space for a local area of particular interest. Since d is directly related to β by (3.7), we have added a β scale in brackets.

Figure 7.

Families of T (blue) and M (red) contours demonstrating an application of a stability theorem. The contour spacing is 0.05 for both T and M parameters. Zero p curve (σ=1) in green, with p positive below. Points (a) and (b) are on the 2w/B=0.2 contour. Point (a) is also on the w/B=0.85 contour.

Marked in black is the fold line, connecting tangencies between the T and M contours. Beneath this line, our stability investigation shows that the equilibrium states are stable, while above it they are unstable. Also in black are two curves showing where m=2, with the m<2 regime between them. The zero pressure line is shown in green, with the regime of positive pressure below it. Notice that the p=0 curve lies entirely in the m<2 regime.

To look independently at the stability of the equilibrium points highlighted by the small (yellow) circles on the 2w/B=0.2 contour, we invoke a recently established theorem (to appear) which checks the sign of the work done by T and M in virtual displacements of d and r to assess the change of the total potential energy. Since this check is within an equilibrium subspace, it can normally predict instabilities but not stabilities. Here, however, since the constrained rod has only two degrees of freedom, it is more comprehensive.

In its simplest graphical form, the theorem allows us to colour the circles shown in figure 7 into sectors where δ2v is guaranteed to be positive (green as if a field), guaranteed negative (blue, ocean) or unknown in sign (yellow, beach). Using this information, and working outwards from the ‘neutral’ fold, we are able to predict that under dead {T,M} loading the total potential energy has the local forms indicated, showing that state (a) is unstable, while state (b) is stable. This agrees with our main dead-load stability result described earlier. A change of this stability only occurs when crossing the fold line, so we see that, along the constant T arc through point (a), zero pressure occurs at a stable state, as we had encountered in figure 6.

A visual inspection of the local forms of δ2v and the associated circles allows us to establish the following results for the mixed loading conditions. Notice that, with only two degrees of freedom (DOF), there can be no virtual movement at all if both d and r are controlled. Approaching the fold from the stable state (b) in this two DOF system the second variation of the total potential energy, δ2v, will have the boat-shaped minimum displayed in figure 7. This becomes the horizontal parabolic valley at the fold, and then the geometrical saddle at the unstable state (a). So at an infinitesimal distance beyond the fold there is only one direction, defined by the critical eigenvector, along which the total potential energy will fall (owing to the nonlinear cubic term). Typically, this eigenvector will not coincide with either of the d or r directions, so if either of the generalized displacements is controlled it is guaranteed to block the instability. So close to the fold, we have:

View this table:

Meanwhile, to locate the stability transitions under the mixed loading conditions, we have just a simple generalization of the tangency rule. Under {T,r} loading, the stability of a path changes where a T contour is tangent to an r=const. line: such tangencies are easily seen in figure 7. Similarly, under {d,M} loading, the stability changes where an M contour is tangent to a d=const. line: this occurs only at the extreme (straight configuration) of d=2.

Finally, when looking at the comprehensive over-view of the helical behaviour, we must first take note of the two sheets of solutions, distinguished by the sign of β, sketched in figure 8. Keeping this three-dimensional structure in mind, we are in a position to understand the comprehensive contour diagram of figure 9, which contains all the various features that we have described above.

Figure 8.

Schematic of the equilibrium surface, with T contours on the upper (β>0) sheet (grey) and on the lower (β<0) sheet (beige). Both are shown in the projection onto the (d,r) base plane.

Figure 9.

Complete results (a) and the top sheet (b). w/B contours: −1.4 to 1.4, steps of 0.2 (increase with r) (red); 2w/B>0 contours: 0.2 to 1.4, steps of 0.2 (increase with r) (blue); 2w/B<0 contours: −1.4 to −0.2, steps of 0.2 (increase with r) (cyan); T=0 contour (purple); critical 2w/B=3/4 contour (green); p=0 curve, p>0 inside (orange); tangency (fold) curves (black).

### (b) Pure compression

Under a pure compressive load, P=−T, with no twisting moment (M=0), a stable equilibrium path emerges from the origin with zero slope as shown in figure 10. An identical mirror-image version (only partially displayed) exists in the negative β regime. This mirror image has the same stability characteristics as the positive β results, and we should note that this makes the origin an unusual singularity, not a fold (with a stability change) like the other two displayed in the figure.

Figure 10.

Solution path for the constrained rod under pure compression. Solid line, stable; dashed line, unstable; filled circles, folds; open circle, singularity.

In fact, the solution at the origin is really a limiting case of a helix with infinite wavelength (L), in which the rod lies straight along (what we might call) the bottom of a horizontal cylinder. As P is increased from zero, the sign of β will be decided by small imperfections in the system. Without loss of generality, we assume that positive β is chosen.

As P is increased from zero, the wavelength, L, steadily reduces as can be inferred from the rising β values. The path then loses its initial stability at the first fold where (54.74°) and the compression parameter is 5.11

From this fold, under dead compressive loading, the rod would jump as indicated. There being no other available solution, the rod (within the mathematical model) would probably pass dynamically through the ‘ring’ solution and end up as a straight rod under tension.

### (c) First approximation (β and small)

In many applications, both the helical angle, β, and will be small, and, in these cases, we can write down some useful approximate solutions which give easy design formulae and graphs. We present here a first-order approximation that is valid in the asymptotic limit where and β are both of order ϵ, where ϵ is a small parameter.

Returning to the basic equilibrium equations of §3, we approximate the second, (3.14), for small ϵ to obtain, after dividing through by β, 5.12 where ρw is the known effective radius of the casing. For given T and M, this single equilibrium equation can be solved for β, with the second remaining degree of freedom, τi, given by (3.1), (3.2) and (3.13). Now within our approximation to first order in ϵ, equation (3.5) becomes L=2πρw/β and for σ defined in (3.6) we find 5.13

Substituting into (5.12) using , we obtain the quadratic equation in σ, 5.14

Notice that this equation can be reinterpreted in the case of negative T if the dimensionless σ and m are defined with P=−T instead of T: the only change needed in the equation is to replace +2 with −2.

Solving for our arc-frequency parameter, σ, we have finally 5.15

This has real solutions for m2>32/9, namely for , so we will first solve it for the free critical buckling value of m=mc=2 to find 5.16

If we take the σ=1 solution, we can go back to physical variables and find the helical wavelength as , equal to the universal wavelength, LU, given for the free theory in (4.6). This solution is labelled as point C in figure 11, and according to our exact analysis it represents a stable state.

Figure 11.

Frequency parameter versus the load parameter m for the first approximation using (5.15). The grey curve shows, for comparison, the exact solution for 2w/B=0.2. Here, the solid line denotes a stable solution while the broken line denotes an unstable solution.

Conversely, if we take the σ=1/2 solution, we obtain L=2LU. As suggested by the passage around a fold bifurcation (and confirmed in our exact analysis) this solution is unstable. It is marked as point B in figure 11.

More generally, if we write the two numerical solutions (5.15) of the quadratic as σ=σ1,2(m) we can write the result for L as 5.17 and the contacting solution for β as 5.18

Numerical results of this first approximation are shown in black in figure 11, where we observe a fold bifurcation at A that corresponds to the vanishing of the square root in equation (5.15) and is given by , and . These results are compared with the exact solution for 2w/B=0.2 in grey. As we would expect, there is good agreement at low β, up to the fold (at A in the approximation). After this, the two solutions diverge spectacularly, with the exact solution exhibiting a second fold that is not picked up at all in the black curve.

The marked pressure change point (p=0) merits some discussion. As we noted in the exact analysis, the zero pressure solution must be one of the free post-buckling equilibrium states. In these states, equation (4.7) gives us σ=1, and figure 6 shows that non-trivial free states exist only for m<2. Now the vertical axis parameter of figure 11, namely , is equal to σ within the first approximation, so, within this, we would conclude that the pressure is zero at C where m=2. This is not possible, and points to a limitation of the first approximation. For the exact solution, the pressure change (marked by a grey diamond) is where the exact σ is equal to unity, and this, plus the shift in the curve, establishes that p=0 is indeed at m<2.

### (d) Second approximation (β small)

If only β is small, we retain one higher order term in β when approximating the T term in the second equilibrium equation (3.14) to obtain, after dividing through by β, 5.19

Introducing the dimensionless cylinder size parameter 5.20 this equation can be written as a quadratic in β as 5.21 with the solution 5.22

The criterion for two real solutions is 5.23

Setting the right-hand side of (5.23) equal to zero allows us to define a critical cylinder size defined by . For a wide casing with γw>γwC the right-hand side is negative, and any m gives us two solutions. Conversely, for a narrow casing, with γw<γwC, the right-hand side is a positive number, and m2 needs to be greater than this number to give solutions.

We should remark here that when considering what is a ‘large’ ρw we face the problem that there is no other physical length dimension with which to compare its magnitude. We can make some progress by interpreting the dimensionless cylinder size parameter, γw, defined in equation (5.20), as follows. Consider the radius of curvature, Rc, produced in the rod by a bending moment equal in magnitude to T acting on a lever arm of length ρw, which is given by Rc=B/(w). Substitution gives , showing that γw is a measure of how large ρw is compared with Rc. This helps a little, but depends on a knowledge of the magnitude of T.

For two numerical studies, we give γw half its critical value, , and a rather extreme figure of 1/30, as displayed in figure 12. These are compared with the exact results, for the corresponding values of T (as given by equation (5.20)), and once again agreement is good for low β.

Figure 12.

Solutions of the second approximation for two values of . The graph for γw=1/30 is indistinguishable by eye from the first approximation and from the exact solution shown in grey but is almost totally masked by the black. The extensive grey curve shows the exact solution for 2w/B=3, which corresponds to half the critical γwC. Here, the solid line denotes a stable solution while the broken line denotes an unstable solution.

The two approximate solutions for the contacting rod give a great deal of insight, and useful universal graphs and parameters, for the practically relevant case of small angles. The second approximation has the advantage of picking up the existence of a critical cylinder diameter.

### (a) Historical drill-string studies

For many years, there has been a continuing thread of interest in the petroleum industry on the helical buckling of drill-strings within their casing. Lubinski et al. (1962) made calculations for the buckling of vertical tubing subjected to compressive loading, without torsion, with and without packers (collars between the tubing and casing). This was re-studied by Mitchell (1982), who included the boundary conditions at the ends where packers enforce alignment with the central axis of the casing. Wu et al. (1993) and Wu & Juvkam-Wold (1993) studied theoretically and experimentally the buckling of a drill-string in extended well bores. He described how a heavy string lying straight on the lower side of a horizontal bore first buckles in compression into a sinusoidal shape along the lower side. This snaking arises because the self-weight of the rod, acting through the cylindrical curvature, provides the equivalent of a lateral elastic foundation. The sinusoidal amplitude increases with further loading, and the deformation eventually becomes helical around the casing. This transition is accomplished by the ‘raising up of every other half-wave to touch the upper side of the wellbore’. Tan & Digby (1993) used an energy approach to study buckling configurations in an inclined hole, including the effect of self-weight, while Miska & Cunha (1995) studied helical buckling in inclined well bores under axial and torsional loads.

Careful post-buckling experiments on stainless steel tubes compressed in a horizontal casing were described by Kuru et al. (1999, 2000) and Martinez et al. (2000). Load deflection curves show successive small peaks as successive helices are formed. The difference between the ‘top load’ and the ‘bottom load’ increases with the number of helices, implying an increase in the lateral contact force and hence an undesirable decrease in the transfer of axial force. Huang & Pattillo (2000) focused on the buckling of an elastic tube in an inclined well bore under the action of its own weight and a compressive force at its upper end. They estimate the load for the onset of helix formation using a Rayleigh approach. These early drill-string papers are usefully reviewed and summarized by Cunha (2004), who lists the various critical loads proposed.

### (b) Recent experiments at Aberdeen

At Aberdeen, UK, a small scale experiment was conducted by two of the present authors (M.S. and M.W.) to observe the behaviour of a rod rotating while experiencing stick–slip oscillations (Silveira 2012). The rig (figure 13) comprises a box of transparent acrylic with a long glass tube inside that simulates the bore-hole.

Figure 13.

Still picture of the small scale experiment with the constrained rod in a helical configuration. An electric motor inputs rotational motion from the right, and frictional torque is generated at the disc–rock interface at the left. The dynamical helical configuration is observed during the brief ‘stick phase’ (see the video in the electronic supplementary material).

Rotating inside the bore-hole is the rod (drill-string), which is a long and flexible steel cable. At the top end, the rod is attached to an electric DC motor (chosen to overcome the maximum resistive torque on the bottom end) with a built-in gearbox, permitting a maximum rotational velocity of the rod of 130 r.p.m. At the lower end, the rod is attached to a metallic disc (drill-bit), which rotates in contact with a rock sample. An axial compressive load is imposed through the rock, generating frictional torque at the disc–rock interface. The load applied to the rock (weight-on-bit, written as WOB) is taken by means of a load cell which withstands up to 70 N, and the velocities of the motor and the disc are measured by two high-precision rotary encoders. The transparent bore-hole enables clear visualization of the torsional buckling experienced by the rod when operating under stick–slip oscillations. The rod dynamically assumes a helical configuration as soon as the disc enters a ‘stick phase’, when torsional loads are higher, and rapidly returns to a straight configuration at the beginning of a ‘slip phase’. This illustrates one way in which a helical deformation can develop dynamically in a realistic drilling scenario.

### (c) Comparison with experiments

Only a few of the experimental results that we have been able to find in the literature allow a direct and relevant comparison with our solutions, which are exclusively focused on the advanced post-buckling state when the whole of the rod is helical and changes in end shortening are simply due to variation of β over the whole length in the absence of friction. Our results are also strictly restricted to conditions in which self-weight effects are small. Despite these considerations, we have thought it worthwhile to superimpose our results onto the experiment of Wu et al. (1993), as shown in figure 14.

Figure 14.

Comparison of our theoretical gradient with the experimental load versus end-shortening results of Wu et al. (1993).

This shows a sketch of Wu's results (the helical stage from his fig. 6) for a brass rod of length 3.048 m and diameter 2.4 mm in a clear horizontal plastic cylinder of inner diameter 25.7 mm. This rod is subjected to the pure compression that we studied in §5b. Before the sketched helical regime the rod had exhibited sinusoidal buckling on the floor of the cylinder, in the manner of a strut on an elastic foundation, followed by a rising towards the top of the cylinder. The curved black line then shows the load deflection response of the ensuing helical distortions.

The corresponding result from our exact theory is 6.1 where the right-hand approximation is amply justified by the very low experimental d values. Now the nature of Wu's experimental plot, including as it does the early non-helical deformations, precludes any sensible comparison with our equation as a whole. However, positioning this straight-line solution close to Wu's advanced helical deformations, we see that the gradients are in excellent agreement.

## 7. Concluding remarks

The buckling and post-buckling of an elastic rod deforming helically in contact with a surrounding cylindrical casing of larger diameter has taken its place as an archetypal problem in modern continuum mechanics. In particular, complex nonlinear phenomena of spatial localization and chaos have been presented in the papers van der Heijden (2001) and van der Heijden et al. (2002). In this paper, we have presented new and extensive Cosserat solutions for the advanced post-buckling stage under torsion and either tensile or compressive loading. Detailed stability and contact-pressure investigations have been made and the results are presented with equilibrium surfaces to add insight and facilitate comparisons with experimental results. One such experimental comparison, in figure 14, is extremely encouraging. Two approximate solutions give universal graphs, parameters and design formulae, for the practically relevant case of small angles, and highlight for the first time the existence of a critical cylinder diameter. The results will be of interest to researchers in many areas of science and technology, including petroleum engineers and surgeons working with stents.