Likely chirality of stochastic anisotropic hyperelastic tubes

When an elastic tube reinforced with helical ﬁbres is inﬂated, its ends rotate. In large deformations, the amount and chirality of rotation is highly non-trivial, as it depends on the choice of strain-energy density and the arrangements of the ﬁbres. For anisotropic hyperelastic tubes where the material parameters are single-valued constants, the problem has been satisfactorily addressed. However, in many systems, the material parameters are not precisely known, and it is therefore more appropriate to treat them as random variables. The problem is then to understand chirality in a probabilistic framework. Here, we develop a procedure for examining the elastic responses of a hyperelastic cylindrical tube of stochastic anisotropic material, where the material parameters are spatially-independent random variables deﬁned by probability density functions. The tube is subjected to uniform dead loading consisting of internal pressure, axial tension and torque. Assuming that the tube wall is thin and that the resulting deformation is the combined inﬂation, extension and torsion from the reference circular cylindrical conﬁguration to a deformed circular cylindrical state, we derive the probabilities of radial expansion or contraction, and of right-handed or left-handed torsion. We refer to these stochastic behaviours as ‘likely inﬂation’ and ‘likely chirality’, respectively.


Introduction
Biological, medical, and engineering applications of the extension-torsion-inflation couplings, in pressurised circular tubes with helical aeolotropy, range from plant and animal tissues [5,23], to cardiac and vascular systems [12,43], to soft actuators [7,8] and dielectric elastomers [17]. The mechanical effects of the stiffness and orientation of the fibres winding helically around the tube's axis were analysed in detail, within the framework of finite elasticity, in [14]. There it was shown that: (i) the tube can lengthen or increase in girth depending on the fibre orientation and the relative stiffness of the matrix and fibres; (ii) when the fibres are equal in stiffness and orientation, the tube may exhibit no rotation; (iii) a tube with a single pre-compressed right-handed fibre may rotate either clockwise or anti-clockwise depending on the material parameters, showing that transfer of chirality under external loads from the tube-wall micro-structure to the tubular macro-structure relies both on geometry and material properties.
In addition, for many solid materials, a crucial part in assessing their physical properties is to quantify the uncertainties in their mechanical responses, which cannot be ignored [18,40,42]. In finite elasticity, the use of the variability in observational data and information about uncertainties [11,55] were proposed in [50][51][52] and [35], where stochastic isotropic hyperelastic models are described by a strain-energy function, with the parameters defined as random variables characterised by probability distributions, while anisotropic stochastic models, where the parameters are spatially-dependent random fields, were discussed in [53,54]. These are phenomenological models, based on the notion of entropy (or uncertainty) [45] (see also [48]) and the maximum entropy principle for a discrete probability distribution [19][20][21], which can also be embedded in Bayesian methodologies [4] (see also [26]) for model selection and updates [35,40,44].
For these stochastic hyperelastic models, the mathematical question arises: what is the influence of probabilistic model parameters on the predicted elastic responses? This question has begun to be addressed analytically in [32,33,36], with the special cases of stochastic hyperelastic bodies with simple geometries for which the finite deformation is known. Examples include the cavitation of a sphere under uniform tensile dead load [32], the inflation of pressurised spherical and cylindrical shells [33], and the classical problem of the Rivlin cube [36]. For these fundamental problems, to which the elastic solution is known explicitly, the sensitive dependence on parameter probabilities was demonstrated by showing that, in contrast to the deterministic problem, where a single critical value strictly separates the cases where instability can or cannot occur, for the stochastic problem, there is a probabilistic interval where the two cases compete, in the sense that both have a quantifiable chance to be found. Such problems can offer significant insight into how the stochastic framework builds on the finite elasticity theory, and revisiting them from the stochastic perspective can offer also opportunities for gaining new insights into the elastic solutions. More complex, but mathematically tractable problems can then be treated in a similar manner [34]. To study the effect of probabilistic parameters in the case of more complex geometries and loading conditions, suitable numerical approaches built on rigorous analysis have started to be developed in [53,54].
In this paper, we examine the elastic responses under dead loading of a hyperelastic cylindrical tube of stochastic anisotropic material, where the model parameters are random variables, which are constant in space, and hence are particular cases of those treated in [53,54]. One can regard this stochastic tube as an ensemble of tubes with the same geometry, where each individual tube is made from a homogeneous anisotropic incompressible hyperelastic material, with the elastic parameters drawn from known probability distributions. Then, for each individual hyperelastic tube in the ensemble, the finite elasticity theory applies. The dead loading consists of simultaneous internal pressure, axial tension and torque, and the resulting universal deformation is a combined inflation, extension and torsion from the reference circular cylindrical configuration to a deformed, also circular cylindrical, state. For the anisotropic material, two non-orthogonal preferred directions are assumed, corresponding to the (mean) directions of two families of aligned fibres embedded in an isotropic matrix material. In the deterministic elastic case, the tube may undergo 'inversion' in the deformation, such that the radius first decreases and then increases, or 'perversion' whereby the torsion chirality changes from right-handed to left-handed [15,27]. The possible existence of such responses depends on the material constitutive model. For the stochastic problem, we derive the probability distribution of the deformations, and find that, due to the probabilistic nature of the material parameters, the different states always compete. In particular, at a critical load, the radius may decrease or increase with a given probability, and we refer to this phenomenon as 'likely inflation', and similarly, right-handed or left-handed torsion may occur with a given probability, and we refer to that as 'likely chirality'.
The outline of this paper is as follows: in Section 2, the kinematics of finite simple torsion superposed on uniform stretch for hyperelastic circular cylindrical tubes is briefly reviewed; Section 3 provides a summary of the stochastic elasticity modelling framework; in Section 4, the probabilistic solution for stochastic hyperelastic tubes is obtained, while some additional technical details are contained in Appendix A; concluding remarks are presented in Section 5.
Through the deformation (1), the circular plane section at Z = 0 remains fixed and each circular plane section normal to the central axis remains plane and rotates by an angle τ ζZ. Denoting the deformation gradient in terms of the current cylindrical polar coordinates (r, θ, z) is equal to and the left and right Cauchy-Green tensors are, respectively, Denoting λ a = a/A, by (2), we can write and Next, assuming that the components of the Cauchy stress, T, are independent of θ and z, the equilibrium equations are reduced to [56, p. 185], Then, by further assuming that T rθ = 0 and T rz = 0, it follows that hence, In this case, when a uniform internal pressure is applied, the outer surface of the cylinder may be rendered free of traction, i.e., T rr = 0 on r = B. Thus, for an internally pressurised tube which is free on the outer surface, where P > 0 is constant. Integration of equation (9), with respect to r, followed by substitution in the first condition of (11), then gives For the deformation (1), by (10), the resultant normal force acting upon the plane a ≤ r ≤ b, zconstant, |θ| ≤ π is calculated as follows [56, pp. 185], Equivalently, where Similarly, the resulting twisting moment is equal to [56, p. 190], After applying the change of variable (7), the expressions of P , F , and T , given by (12), (15), and (16), respectively, take the following equivalent forms, 3 Stochastic anisotropic hyperelastic material In this section, we combine finite elasticity [13,41,56] and probability theory [16,21] to analyse a tube made of a stochastic homogeneous hyperelastic material for which the model parameters are random variables characterised by probability density functions. For each model parameter, the partial information provided by the mean value and variance is commonly used [9,18,25,39]. Specifically, we assume that the tube is made of an aeolotropic material with two preferred directions with respect to the reference configuration, induced by two families of aligned extensible fibres embedded in an isotropic matrix material. The two preferred directions are taken as [14], [13, pp. 328-336] where Φ, Ψ ∈ [0, π/2]. Under the deformation (1), the stretch ratios of the fibres, λ 4 and λ 6 , respectively, are given by where C is the right Cauchy-Green tensor defined by (5).
For the corresponding homogeneous hyperelastic model, the strain-energy function requires seven independent invariants [2,3], including the principal invariants, I 1 , I 2 , I 3 , of the Cauchy-Green tensors B and C [49], the pseudo-invariants I 4 and I 6 given by (21), and two other pseudo-invariants defined as follows, For a review of a wide range of isotropic and anisotropic hyperelastic models, we refer, for example, to [6].
Here, we restrict our attention to the simple case of a stochastic incompressible anisotropic hyperelastic material with two families of extensible fibres embedded in a neo-Hookean material. Hence, we define the following strain-energy function, where µ, µ 4 and µ 6 are positive random variables, which we assume to be stochastically independent [53,54]. This assumption enables us to extend directly the analytical calculations from the deterministic case treated in [14], where µ, µ 4 and µ 6 were single-valued constants, to the stochastic problem, where these parameters are characterised by probability distributions. From the continuum mechanics perspective, in general, anisotropic materials have distinct material properties in different directions. For our model example, the shear moduli are µ in two directions and a linear combination of µ, µ 4 and µ 6 in the third direction. Therefore, assuming that µ, µ 4 , and µ 6 are independent random variables allows for the shear moduli in different directions to be derived as shown in Appendix A.
The shear moduli of other anisotropic hyperelastic materials can be treated analogously. However, the physical behaviour of stochastic anisotropic materials deserves further attention, and we hope that our theoretical analysis may serve as a motivation for future experimental work. For the stochastic model described by (23), the following mathematical constraints guarantee that the random shear modulus of the matrix material under infinitesimal deformations, µ, and its inverse, 1/µ, are second-order random variables, i.e., they have finite mean value and finite variance [32,33,35,36,[50][51][52], Then, by the maximum entropy principle, the random shear modulus, µ > 0, with mean value µ and variance Var[µ], follows a Gamma probability distribution [46,47], with shape parameter ρ 1 > 0 and scale parameter ρ 2 > 0 defined, respectively, by , The corresponding probability density function takes the form [1,22] g(µ; ρ 1 , , for µ > 0 and ρ 1 , ρ 2 > 0, (26) where Γ : R * + → R is the complete Gamma function Similarly, by setting the mathematical expectations [53] E and the random parameters µ 4 > 0 and µ 6 > 0, with mean values µ 4 and µ 6 , and variance Var[µ 4 ] and Var[µ 6 ], respectively, follow the Gamma probability distributions with shape and scale parameters defined, respectively, by and We denote by g 4 u; ρ the corresponding probability distributions. Throughout this paper, computer simulations were created by fixing the parameters given in each figure caption, and repeatedly drawing random samples from the underlying distribution. Our simulations were run in Matlab 2018a, where we made specific use of inbuilt functions for random number generation. Namely, we used "gamrnd" to generate the Gamma distributed random variables and "gamcdf" to generate the Gamma cumulative distribution function.

Inflation and torsion of a stochastic anisotropic tube
For a cylindrical tube of stochastic anisotropic hyperelastic material, with the strain-energy density function described by (23), and subject to the deformation (1), the Cauchy stress tensor takes the form where B is the left Cauchy-Green tensor described by (4), β i = 2∂W/∂I i , i = 1, 4, 6, are the material response coefficients, and p is the Lagrangian multiplier for the incompressibility constraint, det F = 1. The preferred directions given by (20) are deformed into the following directions, respectively, Thus, the non-zero components of the stress tensor given by (32) take the form, Assuming that the tube wall is thin (see Figure 2), we set A = 1 and B = 1 + , and represent P , F , and T , given by (17), (18), and (19), respectively, as the following series expansions [14], where we assume P (0) = F (0) = T (0) = 0. We then truncate the series given by (35), (36), and (37), respectively, to first order in , as follows, where Next, we define the following Jacobian matrix [14] and concentrate our attention on infinitesimal deformations near the reference configuration, with Assuming that det J| (1,1,0) = 0, we can define the matrix inverse This will be useful when exploring the critical points where an inversion in the deformation occurs. In order to obtain clear explicit results that will show the role played by the stochastic parameters, henceforth, we limit our investigation to the case where the anisotropic material has the same mechanical properties in the two preferred directions, i.e., µ 4 = µ 6 . Recalling that µ follows a Gamma probability distribution g(u; ρ 1 , ρ 2 ), defined by (26), and µ 4 follows a Gamma distribution, g 4 u; ρ 2 , the probability distribution of each entry of the Jacobian matrix (44) can be computed using the summation formula for independent Gamma-distributed random variables [24,37] (see Theorem A.1 in Appendix A).
To illustrate our theoretical results, numerical examples of the Gamma distributions considered here are represented in Figure 3, where the shape and scale parameters are ρ 1 = 405 and ρ 2 = 0.01 for µ > 0, and ρ (4) 1 = ρ 1 = 405 and ρ (4) 2 = 20ρ 2 = 0.2 for µ 4 = µ 6 = 20µ > 0. Note that, although Gamma distributions are used throughout this paper, those represented in Figure 3 appear to be approximately normal distributions. This is not a coincidence, and we discussed such a comparison in [33], in the case when ρ 1 was large compared to ρ 2 . However, despite known convergence results and the qualitative agreement between the two density functions for large values of the mean, in general, the normal distribution cannot be used to model elastic parameters. This is due to the fact that the normal distribution is defined on the entire real line, whereas elastic moduli are always positive. Approaches for the explicit derivation of probability distributions for the constitutive parameters of stochastic homogeneous isotropic hyperelastic models, calibrated to experimental data, were developed in [35,52]. As far as we are aware, similar derivations for the stochastic constitutive laws of anisotropic materials are yet to be performed. In practice, elastic model parameters can meaningfully take on different values, corresponding to possible outcomes of the experimental tests. Then, by the maximum entropy principle, one can explicitly construct probability laws using the available information.

Likely inflation of the stochastic tube
First, we assume that Φ = Ψ [14], and consider the critical point for radial stretch λ, such that The radial stretch increases if A 11 > 0 and decreases if A 11 < 0, and equation (47) is equivalent to the following quadratic equation in cos(2Φ), which has real solutions when 0 < µ ≤ µ 4 /18. As the denominator of A 11 is positive, it follows that, as the internal pressure increases, the radial stretch, λ, increases if µ > µ 4 /18, and decreases if 0 < µ < µ 4 /18.
In this case, we can express the probability distribution of stable inflation, such that the radial stretch monotonically increases when the internal pressure increases, as and that of unstable inflation, whereby the radial stretch starts to decreases under increasing pressure, as P 2 (µ 4 ) = 1 − P 1 (µ 4 ).
For example, taking ρ 1 = 405 and ρ 2 = 0.01 (see Figure 3-left), the mean value of the shear modulus is µ = ρ 1 ρ 2 = 4.05, and the probability distributions given by equations (49)-(50) are illustrated in Figure 4 (blue lines for P 1 and red lines for P 2 ). To plot those results, the interval 0, 27µ was discretised into 100 representative points, then for each value of µ 4 , 100 random values of µ were numerically generated from the specified Gamma distribution and compared with the inequalities defining the two intervals for values of µ 4 . For the deterministic elastic case, which is based on the mean value of the shear modulus, µ = ρ 1 ρ 2 = 4.05, the critical value of µ 4 = 18µ = 72.9 strictly separates the cases where radial stretch instability can, and cannot, occur. For the stochastic problem, for the same critical value, there is, by definition, exactly 50% chance that radial stretch is stable (blue lines) and 50% chance that instability occurs (red lines). To increase the probability of stable radial stretch (P 1 ≈ 1), one must consider values of µ 4 that are sufficiently smaller than the expected critical value, whereas an instability is almost certain to occur (P 2 ≈ 1) for values of µ 4 which are sufficiently greater than the expected critical value. However, the inherent variability in the probabilistic system means that there exist events where there is competition between the two cases.
In Figure 5, we represent the probability distribution of A 11 as a function of the angle Φ when µ follows a Gamma distribution with hyperparameters ρ 1 = 405 and ρ 2 = 0.01, and µ 4 follows a Gamma        Figure 3). In this case, µ 4 = 20µ < 81, and unstable radial stretch is expected. Nevertheless, the probability distribution suggests that there is also around 2% chance that the radial stretch is stable.
For the axial stretch ζ, the critical value satisfies the equation and the axial stretch increases if A 21 > 0, and decreases if A 21 < 0. Solving equation (51), we obtain where Φ m denotes the 'magic angle' [13, p. 337]. In this case, as the denominator of A 21 is always positive and (52) is independent of the random material parameters µ and µ 4 , there is no uncertainty to be resolved regarding this instability. For the deterministic cases where µ 4 = 1 and µ ∈ {0.01, 0.02, · · · , 0.1, 0.2, · · · , 1}, the values of A 11 and A 21 are presented graphically in Figure 6.
For the deterministic cases where µ 4 = 1 and µ ∈ {0.1, 0.2, · · · , 1}, the values of A 31 are illustrated in Figure 9, where Ψ 0 ∈ {0, π/30, Ψ m , π/3, π/2} (left) and Φ 0 ∈ {0, π/30, Φ m , π/3, π/2} (right). From the respective plots, we see that a pressurised tube with a right-handed family of fibres deforms by a left-handed torsion if the other fibres are kept either horizontal or vertical (top and bottom left), and a tube with a left-handed family of fibres deforms by a right-handed torsion if the other fibres are horizontal or vertical (top and bottom right). This is the expected behaviour in the case of a single family of fibres as well [14] (see Figure 10).

Conclusion
We have developed here a procedure for examining the nonlinear elastic responses of a tube of stochastic hyperelastic material with two preferred directions under three typical combined loadings. In our approach, by following the usual finite strain analysis, but with the model parameters treated as random variables defined in terms of probability distributions, instead of taking them as single-valued constants as in the deterministic case, leads to mechanical responses that are also characterised in terms of probability distributions. For a direct comparison of our stochastic results with the deterministic ones, we sampled from distributions where the parameters were set to have mean values corresponding to the deterministic system. Thus, the mean value of the distribution was guaranteed to converge to the expected value. Specifically, we found that, while for the deterministic tube, there is a strict transition between different behaviours, given by a single numerical value, for the stochastic tube, our analysis identifies probabilistic intervals where the two different states compete and both have a quantifiable chance to be presented.
In general, the phenomena of inversion and perversion in these systems are primarily driven by the geometry and fibre distribution. In order to clearly reveal the additional influence of the stochastic material parameters, we considered here the mean directions of two families of aligned fibres, without including the stochastic variations in the fibre angles. The additional influence of the fibre distribution has been studied by various authors [10,28], who homogenised the effects of stochastic dispersion, which then entered as a parameter in the deterministic setting. Nevertheless, our approach could also be extended to the problem of fibre dispersion to reveal the full effect of stochasticity.
The present study is part of an ongoing investigation where we aim to illustrate explicitly how the elastic solution of fundamental nonlinear elasticity problems can be extended to the stochastic case [32][33][34]36]. For numerical methods applied to problems that are intractable analytically, we refer to [53,54]. Important applications include soft biological tissues (e.g., plants, arterial walls) and engineered structures (e.g., soft actuators) at large strains, where mathematical models that take into account the variability in the material responses are crucial.

A Random shear moduli of stochastic anisotropic models
In this appendix, we first state an important result concerning the summation of independent Gammadistributed variables (see also Theorem 1 of [37], where a proof is provided). This result is applicable to linear combinations of independent Gamma-distributed random variables by rescaling. We then apply this to derive the random shear moduli of the stochastic anisotropic model (23) under small strains, as limiting cases of nonlinear shear moduli under large strains. The nonlinear elastic moduli presented here can be regarded as the anisotropic analogue of those discussed in [31] and [35] for isotropic materials with deterministic and stochastic parameters, respectively. In the deterministic case, where the model parameters are single-valued constants, the linear shear moduli derived here are equivalent to those obtained through a different approach in [38].
Next, we assume that the undeformed cuboid, with the preferred directions given by (63), is subjected to the simple shear deformation Following analogous calculations to those detailed above, we can define the nonlinear shear modulus and obtain µ 13 (k) = µ.
As µ is constant, the shear modulus described by (75) coincides with its linear elastic limit Finally, assuming that the original cuboid deforms by the simple shear while the preferred directions are given by (63), we define the nonlinear shear modulus In this case also, we obtain µ 23 (k) = µ, with its linear elastic limit µ 23 = lim k→0 µ 23 (k) = µ.
For the deterministic problem, the formulae of the shear moduli (72), (76), (80), under infinitesimal deformation, are equivalent to those described by equation (6.1) of [38]. In addition, for the stochastic problem, the material parameters are random variables defined by probability distribution. In particular, if µ, µ 4 , and µ 6 are statistically independent, and each of them is drawn from a Gamma distribution, then, by Theorem A.1, one can obtain the probability distributions of the shear moduli described by (72), (76) and (80), respectively.