Deformation induced loss of ellipticity in an anisotropic circular cylindrical tube

When a transversely isotropic circular cylindrical tube is subject to axial extension and inflation, the governing equations of equilibrium can lose ellipticity under certain combinations of deformation and direction of transverse isotropy. In this paper, it is shown how the inclusion of an axial shear deformation moderates the loss of ellipticity condition. In particular, this condition is analysed for a material model consisting of an isotropic neo-Hookean matrix within which are embedded fibres whose properties are characterized by the addition to the strain-energy function of a reinforcing model depending on the local fibre direction.

mention below only the limited list of works concerned with the analysis of ellipticity and associated discontinuous solutions for various specializations of the tube problem considered here, and its variations.
Fosdick and MacSithigh [2] studied the helical shear problem for an incompressible isotropic elastic material with an emphasis on the structure of the energy function, with particular reference to its non-convexity and the related emergence of equilibrium states with discontinuous deformation gradients. The azimuthal shear problem for an incompressible isotropic elastic material was studied by Abeyaratne [3] with a focus on loss of ellipticity and the emergence of discontinuous solutions, while for the anti-plane shear problem Silling [4] considered numerically the passage from ellipticity to hyperbolicity of the governing equations resulting from deformation of an incompressible isotropic material containing a crack or a screw dislocation. For an incompressible transversely isotropic elastic material associated with a single family of fibre directions, the problem of loss of strong ellipticity for the azimuthal shear problem was first studied by Kassianidis et al. [5], who examined, in particular, the emergence and disappearance of non-uniqueness of solution. This was extended to the case of two symmetrically arranged fibre families by Dorfmann et al. [6] and El Hamdaoui and Merodio [7].
In the present work, while we adopt the neo-Hookean model for the matrix material, in contrast to the analysis in [1], we consider a different reinforcing model for the description of the anisotropy, which leads to quite different results in general compared with those in [1]. The following sections contain a general set-up of the problem in Sect. 2, including the description of the considered geometry and deformation, the constitutive equations and loss of ellipticity condition. Then, in Sect. 3, the constitutive law is specialized and the effect of axial shear on the loss of ellipticity is analysed in some detail for the separate cases of radial and axial transverse isotropy. Finally, some concluding remarks are provided in Sect. 4.

Geometry and kinematics
In an undeformed and stress-free reference configuration, the geometry of a circular cylindrical tube may be defined in terms of cylindrical polar coordinates (R, Θ, Z ) by where A, B, respectively, are its internal and external radii and L its length. Material points in the reference configuration are labelled by their position vector X, which is defined by are the unit basis vectors associated with (R, Z ). The cylinder is subject to a uniform axial extension and purely radial deformation, superimposed on which is a helical shear composed of a combination of azimuthal shear and axial shear. Described in terms of cylindrical polar coordinates (r, θ, z), this results in a deformation given by the equations where λ z is the (constant) axial stretch of the cylinder and r (R), g(R) and w(R) are unknown functions corresponding to radial deformation, azimuthal rotation and axial displacement, respectively. In principle, these can be determined from the solution of the equilibrium equations and boundary conditions. The corresponding position vector x in the deformed configuration is given by x = r e r + ze z , where (e r , e z ) are the unit basis vectors associated with (r, z). In the following we use the notations for the inner and outer deformed radii.
The gradient of the deformation from X to x is denoted by F and is the second-order tensor given by ∇ X x, where ∇ X is the gradient operator with respect to X. Here we adopt the incompressibility constraint det F = 1. (4) With respect to the chosen axes, the matrix of components of F is denoted by F and given by where λ r = r (R) is the radial stretch, λ θ = r/R is the azimuthal stretch and the shear components γ θ and γ z are defined by the prime signifying differentiation with respect to R. Then, by (4), so r (R) is determined explicitly and the two functions g(R) and w(R) remain to be determined. With T signifying the transpose of a second-order tensor, the Cauchy-Green deformation tensors C = F T F (right) and B = FF T (left) have matrices of components, denoted by C and B, respectively, given by The only principal invariants of either C or B when the incompressibility constraint is in force are denoted by I 1 and I 2 and defined by noting that, by incompressibility, I 3 = det C = 1. Separately from these (isotropic) invariants, transverse isotropy requires two independent invariants associated with the direction of transverse isotropy. This direction is referred to as a preferred direction, defined in the reference configuration and denoted here by the unit vector A, in general dependent on X. The (transversely isotropic) invariants are denoted here by I 4 and I 5 and defined by where a = FA is the push forward of A from the reference to the deformed configuration. Note that in the undeformed configuration I 1 = I 2 = 3 and I 4 = I 5 = 1.

Constitutive equations and equilibrium
The properties of an incompressible elastic material are described in terms of a strain-energy function W (F) (per unit volume), associated with which are the nominal and Cauchy stress tensors, denoted by S and σ , respectively, and given by where p is a Lagrange multiplier associated with the constraint (4), I is the identity tensor and We now specialize these equations by considering a transversely isotropic material with the direction of transverse isotropy identified locally by A. Then W depends on F and A through the combined invariants of C and A given in (9) and (10). Thus, W = W (I 1 , I 2 , I 4 , I 5 ) and, from (11) 2 , it follows that and then, by using standard expressions for the derivatives ∂ I m /∂F, it expands fully into the form where W m = ∂ W/∂ I m , m = 1, 2, 4, 5. The corresponding expression for the nominal stress tensor can be obtained from the connection S = F −1 σ , but is not needed here.
Since the reference configuration is here considered to be stress-free it follows from (13) that p 0 being the appropriate value of p. These restrictions were originally given in [8].
When no body forces are present the nominal and Cauchy stress tensors, respectively, satisfy the equilibrium equations where ∇ X and ∇ x are the divergence operators with respect to X and x, respectively. The first of these will be considered in the following subsection in connection with the concept of ellipticity, while, for the deformation defined in Sect. 2.1, the second has radial, azimuthal and axial components that specialize to Once the constitutive equation and boundary conditions are made explicit, the radial equation can be used to determine the stress component σ rr , and hence σ θθ and p. The azimuthal and axial equations are integrated immediately to give where τ θ and τ z are the azimuthal and axial shear stress components on the deformed boundary r = b. Explicit expressions for σ r θ and σ r z obtained from a given constitutive equation, coupled with (17), can in principle be used to obtain two algebraic formulas (in general implicit and coupled) for γ θ and γ z in terms of the stretches and other parameters of the problem. In general, the solutions for γ θ and γ z may not be unique, as exemplified in the azimuthal shear problem discussed in [5] and [9].

Ellipticity
Returning to the equilibrium equation (15) 1 , we note that in terms of rectangular Cartesian components it can be expressed as where Greek and Roman indices relate to the reference and deformed configurations, respectively, so that X has components X α and x has components We have also used the convention that an index following a comma indicates differentiation with respect to the relevant coordinate, while the summation convention applies to repeated indices. Also, defines the components of the elasticity tensor. Associated with A αiβ j is the acoustic tensor Q(n), whose components are defined by The ellipticity status of the equilibrium equation depends on the properties of Q(n), and, in particular, it is said to be strongly elliptic if [Q(n)m] · m > 0 for all non-zero m and n such that m · n = 0.
We have a special interest in situations where strong ellipticity just fails, i.e.
[Q(n)m] · m ≥ 0 for all non-zero m and n such that m · n = 0, with equality holding for at least one non-zero pair of vectors m and n. Configurations in which (22) holds are those for which (strong) ellipticity is lost and which are identified by the (zero) eigenvalue problem Q * m = 0, where Q * is defined by Q * = Q(n) − n ⊗ Q(n)n and is in general non-symmetric. Since the material is incompressible, the latter forms a two-dimensional problem in the plane normal to n and we therefore project Q * onto this plane asQ =ĪQ * Ī =ĪQĪ, which is symmetric, by means of the projection operator I = I − n ⊗ n, leading to the two-dimensional eigenvalue problem in the plane normal to n. In the context of incompressible isotropic elasticity, the projected acoustic tensor was introduced by Scott and Hayes [10]. The equation determines possible values of the vector n, which, without loss of generality, may be taken as a unit vector. Surfaces with local unit normal n are considered as surfaces of discontinuity, across which the deformation is continuous and the deformation gradient is discontinuous (in the case of a so-called strong discontinuity) or the deformation gradient is continuous but its gradient is discontinuous (a weak discontinuity). The possible emergence of such a surface as the deformation proceeds is associated with the loss of (strong) ellipticity.

Application to I 5 reinforcement
We now specialize the form of strain-energy function in order to illustrate the loss of ellipticity. In a previous paper [1] we examined the influence of an energy function based on the use of the invariant I 4 and independent of I 5 . The effect of I 5 , however, is in general quite different and hence the focus of the remainder of this paper concerns the inclusion of I 5 in the strain-energy function in the absence of I 4 . We then specialize the constitutive law so that the strain energy depends only on I 1 and I 5 , i.e. W = W (I 1 , I 5 ). This difference in W is significant since, while I 4 is the square of the stretch in the preferred direction, the invariant I 5 involves not only the stretch but also shears in directions perpendicular to the preferred direction. The expression (13) for the Cauchy stress tensor then specializes to and the restrictions (14) in the reference configuration reduce accordingly to For the analysis of ellipticity, we require the second derivative of W with respect to the deformation gradient F in (19). This is given by where A particular context that this theory relates to is that of an isotropic elastic material reinforced by a single family of fibres aligned locally in the direction A. To be more specific, we now specialize further and consider W (I 1 , I 5 ) to have the simple form for which (26) is satisfied with p 0 = μ. This consists of an isotropic neo-Hookean base (matrix) material, with positive shear modulus μ, combined with a so-called standard reinforcing model that accounts for both fibre extension/contraction and shearing via I 5 . The parameter ρ > 0 is used to regulate the degree of anisotropy associated with the fibres. Then (30) simplifies accordingly with We note in passing that the counterpart of (30) in [1] (Eq. (33) therein), for which W = W (I 1 , I 4 ), was specialized prematurely with W 1 = constant.

Radial fibres
We now consider the fibres to be purely radial, so that A = E R . Figure 1 illustrates the deformed configuration of a single fibre when the deformation (2) creates fibre contraction for the situation in which m and n lie in the (e r , e z ) plane.
In this case, the invariants I 4 and I 5 are given by and their dependence on λ z for different values of λ θ , γ θ , γ z is illustrated in Fig. 2. The values of each of I 4 and I 5 are increased by the inclusion of one or both of γ θ and γ z , so it is clear that when loss of ellipticity is associated with I 4 < 1, for example, in the absence of shear, then shear will delay the loss of ellipticity. The details of this effect are studied in this and the following subsection. The components of σ in (25) are then easily obtained, and we note, in particular, that In the absence of shear, these two stress components vanish. In this case, with χ = π/2 (n = −e r ),Q has only diagonal components

Fig. 2
Representative plots of the invariants a in I 4 and b in I 5 as functions of λ z for λ θ = 0.85, 1.0, 1.15, the short dashed, continuous and long dashed curves, respectively, and the indicated values of (γ θ , γ z ). In particular, (b) illustrates the non-monotonicity of I 5 . Note that fibre contraction (I 4 < 1) is only possible for relatively small values of the shear components γ θ and γ z , and that I 5 < 1 implies I 4 < 1 and the loss of ellipticity condition (24) is satisfied withQ zz = 0 and m = e z in (23). Expressed in terms of the two stretches λ θ and λ z , this is This loss of ellipticity is associated with a discontinuity surface locally perpendicular to the radial direction, and the relevant failure mechanism is fibre kinking.
In a similar analysis for χ = 0 (n = e z ),Q has only components andQ zz = μλ 2 z , in which case failure of ellipticity corresponds to and is associated with a discontinuity surface locally parallel to the radial direction, the failure mechanism then being interpreted as fibre splitting.
In each case, the dependence of λ z on ρ is illustrated in Fig. 3 for λ θ = 1. For the considered simple strainenergy function, the deformation admits loss of ellipticity for all ρ > 0 in the case of χ = π/2 (Fig. 3a) and for all ρ 1.299 in the case of χ = 0 (Fig. 3b). The strong ellipticity region is below and to the left of the curve in each plot. Note that for a given value of ρ within the considered range, ellipticity is always lost in the case of χ = π/2 prior to that for χ = 0 as λ z increases its initial value 1.
Next, we analyse how the loss of ellipticity surface changes under an applied axial shear. For this purpose, we also illustrate this surface in the absence of any shear deformation (γ z = γ θ = 0) for specific values of the reinforcing parameter ρ, as a reference for comparing the solutions obtained when the shear deformation γ z is involved.
Without specializing the value of χ , the loss of ellipticity condition is obtained from (24) for the deformation (2), the strain-energy function (31) and A = E R . This yields the equation 2ρ λ 4 r − 1 2λ 4 r sin 2 χ 3 cos 2 χ + sin 2 χ + λ 2 r λ 2 z cos 2 (2χ) +16ρλ 8 r sin 2 χ cos 2 χ + λ 2 r sin 2 χ + λ 2 z cos 2 χ = 0. The results with χ plotted against λ z and λ r are illustrated in Fig. 4 for λ θ = 1, so that, by the incompressibility condition (7) 1 , λ r λ z = 1 and the deformation is homogeneous. Results are similar for other values of λ θ , but then the loss of ellipticity condition depends on R. Some comments on this non-homogeneous deformation are reserved for the concluding remarks in Sect. 4. Note that (36) and (38) are special cases of (39). In Fig. 4a, the value of χ obtained from (39) is plotted against λ z for three values of ρ. For each ρ, as λ z is increased from 1 ellipticity is lost at a value of λ z very close to 1 and for χ = π/2, so that the discontinuity surface is perpendicular to the fibre, this being interpreted as fibre kinking. As λ z increases there could in principle be two (symmetric) discontinuity surfaces in the (r, z) plane, one with χ > π/2 and the other with χ < π/2. Also, a discontinuity surface parallel to the fibre with χ = 0 (or π ) is possible, which is interpreted as corresponding to fibre splitting, but, as already mentioned, loss of ellipticity with χ = π/2 always occurs before that with χ = 0. Note that Fig. 4b is equivalent to Fig. 4a for λ r ∈ [1/3, 1] and shows χ plotted against λ r = λ −1 z instead of λ z . Figure 5 illustrates how non-zero shear γ z affects the onset of loss of ellipticity for ρ = 1, 2, 3. As is clear from Fig. 5, the main effect of introducing γ z = 0 is to delay the onset of loss of ellipticity with respect to λ z and to change the angle χ corresponding to loss of ellipticity from χ = π/2 to a pair of values symmetrically disposed about χ = π/2. For larger values of |γ z |, loss of ellipticity is totally ruled out. The loss of ellipticity curves for γ z = 0 are shown for reference. Increasing values of γ z also have the effect of counteracting the influence of increasing values of ρ shown in Fig. 4 which advance the onset of loss of ellipticity. Figure 6 provides a selected alternative view of the results in Fig. 5 for λ z = 1.2, and illustrates the symmetry (γ z , χ) ↔ (−γ z , π − χ). For sufficiently large |γ z | the deformation is strongly elliptic, but as the value of |γ z | is reduced ellipticity is lost for two values of χ , one for positive γ z and one for the negative value of equal magnitude, and the associated values of χ are symmetric with respect to χ = π/2. It is clear that a non-zero γ z delays the onset of ellipticity loss, with the fibre splitting mode (χ = 0 or χ = π ) possible at a smaller value of |γ z | than fibre kinking (χ = π/2). This effect is shown in Fig. 5. Note that we have set γ θ = 0 throughout. Indeed, with the choice we have made that m and n lie in the (r, z) plane, it can be shown that γ θ = 0 is a necessary consequence of (23) except for very special values of χ and/or the deformation. This can be shown by expressingQ in the form where α, β and γ are lengthy expressions depending on the deformation, the reinforcing parameter ρ and the angle χ . For convenience, the expressions for α, β and γ are given in the appendix, although γ is not needed here. We mention here, though, that β has the form β = μργ θ f (χ , λ θ , λ z , γ θ , γ z ), where f is a known function of its arguments. Using (40), Eqs. (23) and (24) yield From the first of these, α = β = 0, and the second is then automatically satisfied. Vanishing of β requires either γ θ = 0 or that χ is independent of the reinforcing parameter ρ. The latter is inconsistent with the loss of ellipticity condition α = 0, from which χ depends on ρ, so it is therefore necessary to set γ θ = 0, which has been done in Figs. 5 and 6.

Axial fibres
In this subsection, we consider the fibres to be purely axial, so that A = E Z . Figure 7 illustrates the deformed configuration of a single fibre when the deformation (2) creates fibre contraction, again with m and n lying in the (e r , e z ) plane. In this case, the σ r θ and σ r z components of σ in (25) simplify to with the invariants I 4 and I 5 given by In contrast to the previous section, the deformed fibre length √ I 4 is not affected by either the axial or azimuthal shear, and the expressions for I 4 and I 5 are simpler.
In the absence of shear deformation (γ z = γ θ = 0), the loss of ellipticity condition can be obtained from (39) by interchanging the roles of λ r and λ z , and, bearing in mind the different ranges of values of the angle χ , plots equivalent to those in Figs. 3 and 4 can be obtained but are not included here. The explicit loss of ellipticity condition α = 0 is independent of γ θ and yields (with the factor μ omitted) (λ r sin χ + γ z cos χ) 2 + λ 2 z cos 2 χ + 2ρλ 2 z (I 5 − 1) λ 2 r cos 2 2χ Fibre contraction in this case requires λ z < 1, and the next two figures illustrate the change in the loss of ellipticity curves under the applied axial shear deformation γ z for specific values of the reinforcing parameter ρ based on (45). In Fig. 8, for ρ = 1, with λ z decreasing from 1, only the fibre kinking failure mode arises, and non-zero γ z delays the onset of loss of ellipticity, which disappears altogether for γ z of sufficiently large magnitude. A similar trend can be seen for ρ = 2 and ρ = 3, with both kinking and splitting modes being possible for small values of γ z for ρ = 3. In Fig. 9, some plots of χ against γ z corresponding to loss of ellipticity are shown for specific values of ρ, analogously to those of Fig. 6, and the interpretation is similar to that in Fig. 6.

Concluding remarks
The results in Sect. 3 serve to illustrate how axial shear affects the onset of loss of ellipticity when superimposed on axial extension and radial deformation in a fibre-reinforced elastic solid for the situations in which the fibre direction is either radial or axial and there is no azimuthal shear. Very similar results can be obtained for the case of azimuthal shear in the absence of axial shear with the vectors m and n in the (r, θ) plane instead of the (r, z) plane. For other fibre orientations, the resulting loss of ellipticity condition is quite complicated, and its analysis requires a separate and purely numerical approach.
The model for the fibre reinforcement adopted here was taken to depend on the invariant I 5 and leads to results that are significantly different from those in [1], where an I 4 -based model was used. In particular, while the loss of ellipticity associated with the I 4 model admits only discontinuities related to fibre kinking, the I 5 model may be related to both fibre kinking and fibre splitting, and, in the case of kinking, the discontinuities here develop differently from those in [1]. The main difference between the I 4 and I 5 invariants themselves, which influences the loss of ellipticity results through the strain-energy function, is that the former involves only fibre stretch while the latter relates to both fibre stretch and, in general, shearing, although it may be misleading to suggest that fibre splitting is directly related to shearing alone. Consideration of a strain-energy function which combines I 4 and I 5 would provide an obvious extension of the analysis considered here and in [1], although necessarily more complicated, and it would be of interest to investigate this further.
In the numerical examples in Sect. 3 we have assumed that λ θ = 1, i.e. r = R, which means that, in the absence of shear, the deformation of the tube is homogeneous. When an axial shear is in place it depends on r and, depending on the constitutive law, it may result in the loss of ellipticity surface located at different values of r ∈ (a, b), with γ z having different values on either side and part of the tube in the strongly elliptic regime. This effect is modified to some extent when λ θ = 1 and the deformation is inhomogeneous in the absence of shear.