Finite deformation of a rotating orthotropic cylinder with linear elasticity

conrpurcrs & structrrres, vol. 4. pp. 58 l-591.

Pergamon Press 1974. Printed in Great Britain

FINITE DEFORMATION OF A ROTATING ORTHOTROPIC CYLJ NDER WITH LINEAR ELASTICITY B. E. SANDMAN Analysis Staff, Naval Underwater

Systems Center, Newport, Rhode Island 02840, U.S.A.

Abstract-The effects of finite deformation upon a rotating, orthotropic cylinder with linear elasticity is investigated. The governing equation and boundary conditions form a non-linear two-point boundaryvalue problem. A numerical integration technique is used in conjunction with the related initial-value problem and Newton’s method to obtain approximate solutions. The problem and method of solution provide an example for treatment of non-linear boundary conditions. The stresses induced by finite deformation and the influence of the degree of orthotropy are studied.

NOMENCLATURE

r, 0 polar coordinates a, b inner and outer radii of cylinder, respectively 0

IT,

fJtlee radial and circumferential 0

u %8

ClL, P a9 P

n

F, G, Y, Z 0

SK S' R

u x

stresses, respectively angular velocity radial displacement radial and circumferential strains, respectively elastic constants mass density of cylinder ratios of elastic constants relative inertia force (Z x 1) vector functions (Z x 1) null vector sequence of approximation vectors sequence of roots corresponding to the 0’ value of R radii ratio (=b/a) dimensionless radial displacement dimensionless radial position. INTRODUCTION

flexibility is emphasized and applied to an advantage in modern, light-weight designs. The analysis of flexible structures often becomes exceedlingly complex due to the difficulties involved in the consideration of finite deformation. Recently, considerable effort was put forth in an attempt to develop a high-performance flywheel [l, 21. Research and development has revealed that highly orthotropic cylinders can be utilized to achieve a desirable distribution of stress during rotation. A reduced radial modulus diminishes the radial stress and shifts the maximum circumferential stress toward the outer radius of the cylinder. However, a highly elastic wheel results in larger magnitudes of strain. As indicated in Ref. [2], a composite filament fiber and elastomer matrix cylinder will exhibit STRUCTURAL

581

582

B. E. SANDMAN

strains of the order of 10 per cent when performing at design speeds of rotation. Strains of this order of magnitude may influence the stresses considerably due to the higher order effects of finite strains. An accurate prediction of the stresses and displacements developed in such a cylinder requires an analysis which accounts for finite deformations. Green and Zerna [3] have solved the title problem for an incompressible and isotropic material. The consideration of orthotropy and compressibility results in considerable complexity and the present method offers a direct means of solution. In this paper the equation of equilibrium is written in a form which accounts for finite radial displacements. Linear elasticity is assumed and the Lagrangian non-linear strain tensor is introduced to obtain the equilibrium equation in terms of the radial displacement. Thus, the governing equation of elastostatics is a nonlinear form of Navier’s equation. In particular the equation is one-dimensional and constitutes a non-linear two-point boundaryvalue problem when the appropriate boundary conditions are introduced. A numerical study of these equation is conducted by introducing the related initial-value problem and implementing an initial-value method [4, 51. The influence of finite deformations is studied and related to the degree of orthotropy.

GOVERNING

EQUATION

OF EQUILIBRIUM

Consider the homogeneously orthotropic annular cylinder of constant thickness shown in Fig. 1.

FIG. 1.

The equation of equilibrium of symmetrical forces acting on the differential element shown in Fig. 1 is [6] f (J-a,,) - (Tee+ prF, = 0

(1)

where F, is the body force per unit mass. Since the element located at the undeformed radius, r, is displaced to the position, r+zt, during rotation, the body force is written in the form F,=wZ(r+u)

Finite Deformation

of a Rotating Orthotropic

Cylinder with Linear Elasticity

which includes the second-order increase in inertia force induced by deformation. the equilibrium equation becomes

583

Thus,

as found in Ref. [7]. The stress-strain relations for a cylindrically orthotropic cylinder are given by (W g&-=C, 1% + Cl 2%? (3b)

where .srr and see are the normal strains in the radial and circumferential directions, respectively a,, and a, are the corresponding components of stress, and the C, are the elastic constants. It is noted that equations (3) imply perfectly elastic behavior. The straindisplacement relations for finite, symmetrical deformations are [6j

with u=u(r) describing the radial displacement as a function of position. equations (3 and 4) into equation (2)yields

which is the equation of equilibrium in terms of displacements. By introducing the dimensionless quantities

XC! Q

equation (5) can be written in the form

“==.~22!Gl

Substitution of

584

B. E.

SANDMAN

Equation (6) governs the finite deformation of a rotating orthotropic cylinder with linear elasticity. For R=O, the problem reduces to the classical linear case. The problem becomes completely defined when equation (6) is accompanied by a set of boundary conditions. BOUNDARY CONDITIONS The boundary conditions for equation (6) are prescribed in terms of the displacement at the dimensionless radii, x= R and x= 1. For a cylinder with unconstrained edges the radial stress must vanish at the boundaries or equivalently (7) where R=b/a is the ratio of inner/outer ratii. Equations (6 and 7) constitute a non-linear, two-point boundary-value problem. The fact that the boundary conditions form a nonlinear functional adds to the complications involved in the solution of the problem. METHOD OF SOLUTION The solution of the boundary-value problem posed by equations (6 and 7) would be attempted by the use of a perturbation method based on the parameter II. However, the algebra becomes exceedingly difficult and the practical result is only valid for small values of Sk Hence, the ‘shooting’ or initial-value method is implemented for a numerical study of approximate solutions to the non-linear two-point boundary-value problem. The procedure is outlined as follows. The governing equation (6) can be written as the system of two first-order equations

$- =F(x,

Y; a),

R
(8a)

where

and F is the appropriately defined (2 x 1) vector function. are written in the generalized form GPXR), Y(l)l=O

The boundary conditions (7)

@b)

where G is a (2 x 1) vector function. A detailed description of the vector functions F and G is given in the Appendix. The following discussion is a theoretical format for the numerical solution of equations (8). The non-linear boundary-value problem described by equations (8) may be studied conveniently by introducing the related initial-value problem.

dT;=F(x,

Z; i2)

Z(R)=S

(94

Finite Deformation

of a Rotating

Orthotropic

Cylinder with Linear Elasticity

585

where S is a (Z x 1) vector of initial values. Let a solution to equations (9) be denoted by

s, 0).

Z=Z(x;

(10)

The secondary arguments indicate that the solution depends upon the initial values S and the parameter a. Substitution of equations (10) into equations (8b) yields G[Z(R; S, Q), Z(l;

S, a)]=0

(11)

which may be considered as two equations in terms of the two unknowns S with parametric dependence in R. With the value of R fixed, say R=Q”, equations (11) reduce to two equations G[Z(R; S,

Sz’),Z(l; S, S-lo)]=0

(12)

in terms of two unknowns S. Newton’s method [8] provides a direct means of obtaining a root So of equation (12). Starting from an initial guess So, the iterative sequence Skfl=Sk+ASk; is generated.

k=O,

1,2,.

. .

UW

The correction vector A& is obtained by solving the linear equation

Equation (13b) is obtained by rejecting the higher-order terms of a Taylor series expansion about Sk. The (z x 2) Jacobian matrix

gives the rate of change of the error with respect to a change of initial values S at the kth step of iteration. If the starting guess So is located in a sufficiently small neighborhood of So, then convergence of the sequence SC to the root So can be expected. The nonlinear character of the vector function F is equations (9) prohibits an explicit solution. Hence, the Jacobian matrix (14) cannot be evaluated directly, and a method of constructing the Jacobian matrix must be devised. Differentiation of the initial-value problem (9) with respect to S yields the variational problem

%(E)=(E)(g az

(as>

=(I)

(15b)

x=R

which is a system of four first-order equations, and (Z) is a (2 x 2) identity matrix. For a given vector S and Q=R’, this derived problem and the initial-value problem (9) may be

586

B. E.

SANDMAN

integrated simultaneously, at least numerically, from x=R to x= 1. Evaluation of the resulting solution to the variational problem at x= 1 and substitution into

provides the Jacobian matrix corresponding to the given values of S and Q=R”. Thus, by setting S=S, and integrating equations (9) and (15) from x= R to x= 1, the first corrected vector Sr can be calculated using equations (13 and 16). Returning with S=Sr and repeating the same operations, one obtains the second corrected vector Sz. Successive repetition of this procedure yields the desired sequence, S,. In summary, the two initialvalue problems (9 and 15) in conjunction with Newton’s method provide a convenient mode of investigating the two-point boundary-value problem (8) when the value of Cl is held fixed. The analysis of the non-linear system (20) is complete when the functional relation S=S(n) and accompanying solutions

Y(x;R)=Z(x; are established. manner, i.e.

s(n),

n>

The initial-value method may be used to accomplish this in a discrete S’=S’(@);

i=O, 1,2, . . . , m.

After having obtained a root, So, which corresponds perturbed O=ClO+Ano=Q1.

to R=R’,

(17) the value of Q can be

For this value of R iteration is reinstated starting from S=S”. If AR0 is not exceedingly large, then So is contained in the convergence domain of Newton’s method and iteration yields the root S’ corresponding to O=R r. Successive repetition of this analytic continuation approach leads to the relation given in (17) with fii+l=Qi+A@;

i=o, 1, . ,

. ) HZ- 1.

It is noted that the present analysis avoids the singular point in equation (6) by monitoring the magnitude of R. NUMERICAL

INTEGRATION

The theoretical analysis presented above suggests the employment of a numerical integration technique. Thus, by integrating the initial-value problems (9 and 15) simultaneously with a fourth-order Runge-Kutta-Gill method and performing the successive corrections of S according to Newton’s method, approximate solutions to the governing equation (8) can be obtained. First, the problem is reduced to the linear case by setting Q equal to zero, R’=O. In this case, a numerical solution is obtained in one corrective step due to the truly linear character of equations (8) when R=O. The corresponding approximate vector So is stored and the solution is recorded. From this linear solution the effects of finite deformations can be studied with the continual use of the concept of neighboring solutions. By gradually

Finite Deformation of a Rotating Orthotropic Cylinder with Linear Elasticity

587

incrementing the value of n and restarting the correction and integration procedure from the solution vector S obtained in the solution corresponding to the previous value of R discrete representations of the Q-dependent family of solutions are found. Numerical solutions to equations (9 and 13)were obtained with the Runge-Kutta-Gill method using a step-size, Ah= l/40. For each element Sz’ of a sequence of discrete values of R successive corrections of S were performed until the initial values S and the final values Z(1) satisfied the norm max , I I If712 (gp) so.1 x lo-” (

I

where (gp)=G. The resulting solutions to equations (8) were recorded in sequential order. The R-dependence was established by letting ai+ ’ =n’+ AQ’, 0 5 IAQil-< 1 .O.

10.0 -

“i *.$z

ISOTROPIC

8.0 1 t ;i

H

a= I.D.

a f 0.3

R -0.5 6.0 -

E Y)

P

._5; s B ‘El

4.0 -

2.0 -

I 0.0

I 0.2

0.4

, 0.6

Radial Coordinate

0.8

I.1

x

Fit;. 2. Radial stress across dimensionlessradius s. Figures 2-7 illustrate the results that were obtained for a rotating cylinder with R=O.S and various degrees of orthotrophy. The stresses were computed according to the following equations

CT& =

2+y[l+y+,~[l+~]. a

PW

These equations were obtained by substituting equations (4) into equations (3) and converting the results to dimensionless forms. Figs. 2-7 illustrate that the stresses and displacements are not truly linear functions of the inertial force po2u2. A linear increase in pw2a2 results in an increase of stresses and displacements which is greater than that predicted by the classical linear theory. By comparing Figs. 4 and 7, it is seen that an

588

B. E.

SANDMAN

increase in the orthotropocity factor a decreases the magnitude and non-linear Rdependence in the displacement u(6) when one considers equivalent values of Cr,. For large values of a the influence of fi is small. High orthotropocity in the form of large a (Fig. 6) causes the m~imum values of circumferential stress and greatest Q-dependence to occur near the outer radius. This is in contrast with the result shown in Fig. 3 for an isotropic cylinder. Likewise, high orthotropocity shifts the maximum radial stress and its greatest n-dependence toward the outer radius of the cylinder (compare Figs. 2 and 5). Thus, for a cylinder with R=0.5 and a= 100, the likelihood of catastrophic failure in the form of crack propagation from the interior of the cylinder is reduced.

0.0

0.4

a2

L

0.6

0.8

.O

Radia!Coordisatex

FIG. 3. Circumferential

stress across dimensionless

1.3 _*Cr 9 _& 1.2 c* >,a -G i

0.9

__--

_e--

ISOTROPIC

radius x.

_..---

a- 1.0.p.0.s

R =0.5

~

Fro. 4. Relative displacement

induced by inertial force.

Finite Deformation

of a Rotating Orthotropic

Cylinder with Linear Elasticity

“5 L

uL

4.0

b*

ORTHOTRO+'IC

I

I

a0

0.4

0.2

I

0.6

0.8

I

RadialCoordinatex

FIG. 5. Radial stress across dimensionless

radius x. -I

ORTHOTROPIC '3 *"Q b"

a*Kl.o. 8.0 -

,I“

p=o.,

R-O.5

2 E z 3

6.0-

g J

..l

0.2

0.0

0.4

0.6

0.8

1.0

Radial Coordinatex

FIG. 6. Ci~umferential

stress across dimensionless

radius x.

0’ IW.0. p. 0.5 E F2 0.004 ____________-______

Y :: P 'Z 2 6

Q~IW.0,

p.o.3

O.M3

j ORTHOTROPIC

0.0

R.O.5

10.0 Relative tnertial

Fro. 7. Relative displacement

20.0

Force _Q

induced by inertial force.

589

590

IL E.

SANDMAN

CONCLUSION The initial-value or ‘shooting’ method provides a convenient means of computing numerical solutions to the complex two-point boundary-value problem given by equations (6 and 7). A variable thickness cylinder could be considered without difficulty. The results show that an elastic cylinder which can withstand considerable deformation is subject to a more rapid increase in stresses than that predicted by linear theory as the angular velocity is increased. This factor must be considered in proper design. High values of the orthotropicity factor a are desirable for reducing the rapid growth of radial displacement u as w increases. Increasing the value of a results in a shift of the maximum stresses and maximum R-dependence of stresses toward the outer radius of the rotating cylinder. Although these results are presented for the specific radii ratio R=OS, they exhibit the general character of results obtained for 0.1 I RsO.5. Previous calculations of maximum :tress levels in filament-elastomer flywheels utilizing linear deformation theory are 15-20 per cent below the levels predicted by the present non-linear theory which accounts for significant radial expansion. Consequently, indicated capabilities for energy storage must be revised to include the effect of increased stresses due to flywheel growth. It should be mentioned that the above analysis is based upon the assumption of linear elasticity. Many extremely elastic materials have moduli which are strain sensitive. A modification of this investigation to include non-linear elasticity would improve the results and extend its applicability. However, the treatment of the non-linear boundary conditions without difficulty offers the possibility of generalized applications outside the field of elasticity. REFERENCES [I] A. E. WETHERBEE,Feasibility of flywheel high energy storage.

First Progress Report, PWA-2717, Contract Nonr-4880(00), Pratt & Whitney Aircraft, November (1965). [2] A SVEN~~~Nand A. E. WIXHERBEE.Feasibility of flywheel high energy storage. Final Report, PWA-2676, Contract Nonr-4880(00), Pratt & Whitney Aircraft, April (1969). [3] A. E. GREEN and W. ZERNA, Theoreticul Elasticity. Oxford University Press, Oxford (1954). [4] H. B. KELLER, Numerical Methods for tie-Point Boundary-Value Problems. Blaisdell, New York (1968). [5] A. N. SHERBOURNE,Numerical methods in bending and buckling of plates. J. Str. Div. ASCE 89, (1963). [6] H. L. LANGHAAR,Energy Methods in AppIied Mechanics. Wiley, New York (1962). [7] E. J. BRUNELLE,Stress redistribution and instability of rotating beams and disks. AIAA J., 9, (1971). [S] L. COLLATZ,Functional Analysis and Numerical Mathematics. Academic Press, New York (1966). (Received 16 March 1973)

APPENDIX The explicit form of the system of first-order equations which define the related initialvalue problem used in the Runge-Kutta-Gill integration is dZ1 -dX

-z

2

Finite Deformation

of a Rotating Orthotropic

Cylinder with Linear Elasticity

591

where .=_“;‘[1+;Z2]+-2L[1+;

~]+R~[Z,-+A2Z,

D=l+QZ, and Z,(x)= U(x). The initial conditions are

Z2(R)=S2.

The variational equations as given symbolically in equations (15) are obtained by a formal differentiation of the above equations with respect to S, and S,. The vector function G which defines the boundary conditions is given as follows in component form

g,=z*[ 1+nlz,]

+,:[

1+C3]b;I

=o

Differentiation of (gI, g,) with respect to (S,, S,) yields the Jacobian matrix ofequation (16. I) which is exemplified below

(iq=(Z ;) 1

2

where

ag2_ a!?$) [l +nz,(l>] --

as, a2 -=

as2

+/I$$)

1

a?

[i +z~z,(i)] 2

[l +nz,(1)] 1

+j3%[1

+Rz,(i 2

j].