VOL. 10, NO. 2, MARCHAPRIL 1987
J. GUIDANCE
139
Dynamics of a Cantilever Beam Attached to a Moving Base
T. R. Kane* and R. R. Ryanf
Stanford University, Stanford, California
and
A. K. Banerjeerj: Lockheed Missiles and Space Company, Sunnyvale, California
The behavior of a cantilever beam built into a rigid body that is performing a specified motion of rotation and translation is studied with two objectives in mind. First, because the subject is of interest in connection with spacecraft antennae, helicopter rotor blades, robot arms, and other systems that perform complex motions, we create an algorithm that can be used to predict the behavior of the beam when the base undergoes general threedimensional motions. Effects such as centrifugal stiffening and vibrations induced by Coriolis forces are accommodated automatically, rather than with the aid of ad hoc provisions. The second objective is to draw attention to fundamental flaws in certain multibody computer programs currently under development or already in use. To this end, we construct a second simulation algorithm, one that embodies the procedure apparently employed in the programs in question, and then compare simulation results produced by computer programs based on the two algorithms. Conflicts between the two approaches that thus come to light are discussed in detail.
S
I. Introduction TUDY of the behavior of flexible bodies attached to moving supports has been vigorously pursued for over fifty years in connection with a number of diverse disciplines, such as machine design, robotics, aircraft dynamics, and spacecraft dynamics. In particular, beams attached to moving bases have received attention in hundreds of technical papers dealing variously with elastic linkages, rotating machinery, robotic manipulator arms, aircraft propellors, helicopter rotor blades, and flexible satellites. Indeed, the existing literature is so voluminous as to preclude comprehensive review within the confines of a technical paper. (A 1974 review, by Modi,1 of the literature on rigid bodies with flexible appendages contained more than two hundred references.) Accordingly, we shall mention only representative publications that permit us to direct attention to issues of primary concern. In the area of machine design, the published analyses2"5 deal mainly with planar bending of rotating links or with pure torsion of shafted systems. For the most part, the goal of these analyses is to determine forces or stresses of interest in connection with machine component sizing. In robotics,6"9 displacements, such as those giving rise to bending of a manipulator arm, are at the center of attention, and various beam effects are frequently neglected in an attempt to increase online computational speed. Analyses performed in connection with aircraft dynamics10"14 address questions concerning tapered, twisted, and rotating beams; and in the field of spacecraft dynamics,15"21 the effect of vehicle elasticity on attitude motions is of particular interest. What is missing from this vast literature is a comprehensive theory for dealing with small vibrations of a general beam attached to a base that is undergoing an arbitrary, prescribed motion (translation and rotation). By a general beam we mean one that has an arbitrary cross section and material properties that vary from point to point along the beam; and a compre
Received Aug. 30, 1985; revision received April 8, 1986. Copyright ? American Institute of Aeronautics and Astronautics, Inc., 1987. All rights reserved. * Professor of Applied Mechanics. fResearch Assistant. Member AIAA. :j:Staff Engineer.
hensive theory is one that accommodates all of the following: extension, bending in two principal planes, torsion (uniform or nonuniform), shear displacements, and warping. In the existing literature, one often finds that a particular effect is introduced in an ad hoc fashion. For example, there are many analyses that make special provisions for centrifugal stiffening, Coriolis forces, etc. Most of these analyses are based on one of two methodologies. The first22 involves starting with equations applicable to a beam with an inertially fixed support, and then adding terms intended to reflect an effect of interest, such as the Coriolis acceleration associated with a rotational motion of the base that supports the beam. In the second approach,21'23 which is based on a RayleighRitz technique, one works with trial functions or assumed modes, as these functions are often called, some of which are directly associated with the gyroscopic effects or centrifugal effects present in a beam performing a simple spinning motion at a constant spin rate. These two approaches share a basic shortcoming: they do not account, in a general way, for arbitrary base motion. While most previous analyses deal either with spinning beams or beams attached to a base undergoing rectilinear translation, some attempts to construct more comprehensive theories have been undertaken. For example, Benedetti,24 formulating equations of motion for a uniform EulerBernoulli beam with a symmetric cross section, let the base supporting the beam move with arbitrary angular and translational accelerations; however, terms associated with stiffening due to rotation were neglected because the analysis was intended for the study of the behavior of flexible projectiles, in which connection it is appropriate to drop these terms since the dominant angular velocity component is parallel to the axis of the projectile. The fact that this work is restricted to beams with symmetric cross sections precludes its use for the purpose of analyzing motions of open section beams, such as STEM (selfstoring tubular extensible module) booms, which are commonly used on satellites. Laskin et al.25 also attempted td develop a general theory for a freefree, symmetric cross section Rayleigh beam, but stopped short of the mark and developed final equations only for the BernoulliEuler beam with a quasisteadystate stretch. Modi ended his literature review with the suggestion that in future work "more attention should be directed towards evolution of general analyses rather than exploration of specific
140
KANE, RYAN, AND BANERJEE
J. GUIDANCE
configurations." It is the principal objective of this paper to present such an analysis. In addition, we wish to draw attention to the fact that rather subtle errors can lead to theories with serious defects. Specifically, there exist multibody computer programs, such as DISCOS,26 NBOD2,27ALLFLEX,28'29 and TREETOPS,30 for which a claim is made in the associated documentation that these programs can produce accurate simulation results of "large" motions of systems containing flexible bodies that are themselves concurrently undergoing "small" deformations. In fact, unless ad hoc corrective measures are taken, these programs can produce totally incorrect results in certain situations. The remainder of this paper is organized as follows. The system to be studied is described in the next section. In Sec. III, the new theory is set forth. An associated algorithm that can be implemented on a digital computer is described in Sec. IV. An illustrative example involving a robot manipulator arm with general base motion is presented in Sec. V, and the present approach is compared with the "conventional" one in Sec. VI. Finally, in Sec. VII, conclusions are drawn, and plans to incorporate the new theory into a multibody computer program are discussed.
section for the b2 and b3 directions. The parameters e2(x) and e3(x) are measure numbers of the b2 and b3 components of the eccentricity vector, which extends from the elastic center of the generic cross section to the centroid. The elastic center, the flexural center, and the center of twist all are assumed to coincide in the present analysis.32
III. Equations of Motion
The goal of the following analysis is to produce equations governing the extension, transverse bending, shear, and torsion of the beam introduced in the previous section. The reader who is interested primarily in the resulting computational algorithm rather than details of the equation formulation should proceed directly to the following section. In the sequel, an underlined symbol denotes a vector, a doubly underlined symbol signifies a dyadic, and a symbol in boldface type represents a matrix. Numbers listed over equal signs refer to previous equations, and a triangle over an equal sign indicates a definition. The motion of A in N is characterized by the six quantities vl9 v2, v3,'(?)l9 ? 2 > anc^ W3 defined as 0' = 1,2,3) (1) (2)
II. System Description
The system to be analyzed, shown in Fig. 1, consists of a cantilever beam B built into a rigid body A whose motion in a Newtonian reference frame N is prescribed as a function of time. The beam is characterized by a natural length L, material properties E0(x)9 G(x)9 p(x)9 and crosssectional properties
where v° is the inertia! velocity of the point O (see Fig. 2) on the elastic axis§ of B at the root and NwA is the angular velocity of A in N. Equations (1) and (2) imply that
N
N
e3(x)9 defined as follows. Let x be the distance from a point 0, located at the root of B, to the plane of a generic cross section of B, when B is undeformed. Then E0(x)9 G(x), and p ( x ) are the modulus of elasticity, the shear modulus, and the mass per unit length of the beam, respectively, as functions of x. The area of the cross section located at a distance x from point O is denoted as A0(x)9 and the Saint Venant torsion factor and the warping factor are represented by the symbols K(X) and F(JC), respectively. In order to define the symbols I2(x)9 I3(x)9 <*2(x), a3(x)9 e2(x)9 and e3(x)9 introduce a dextral set of mutually perpendicular unit vectors Jbl9b29b39 fixed in the plane of the cross section located at a distance x from point O9 and oriented such that b± is parallel to the centroidal axis of B, while b2 and b3 are parallel to central principal axes of the cross section. Additionally, introduce a similar set of unit vectors ql9 a2, a3, fixed in A and parallel, respectively, to bl9b29b3 when B is undeformed. The symbols I2(x) and I3(x) denote the central principal second moments of area of the generic cross section for unit vectors b2 and b3, respectively, and ct2(x) and a3(x) are the shear area ratios (also called shear correction factors, see Ref. 31, p. 351) of the
v° = vlal + v2a2 + V3q3
(3) (4)
and
To describe the motion of B in N, introduce a rigid differential element dB of B9 with centroid C and elastic center E located at a distance x from point O when B is undeformed. The differential element dB can be brought into a general orientation in A by aligning blfb2)b3 with *h,*?2>^3> respectively, and subjecting dB to successive dextral rotations in A of amounts 0l9 029 and 03 about lines parallel to a, q29 and q39 respectively. Similarly, point E of dB can be brought into a general position in A by subjecting E to successive translations MX^X, U2q29 and w 3 0 3 . If the vectors r and w are defined as
r^xal (5)
u=
(6)
then the position vector p°E from O to E can be expressed as
N
generic cross section
(7)
The position vector pEC from the elastic center E to the centroid C is the eccentricity vector; that is,
(8)
Therefore, the position vector poc from O to C is given by
pOC=pOE+pEC (9)
Fig. 1 Beam attached to moving rigid base.
§The elastic axis is the line along which transverse loads must be applied in order to produce bending unaccompanied by torsion of the beam at any section.32 This axis passes through the elastic center E of every section.
MARCHAPRIL 1987
CANTILEVER BEAM ATTACHED TO A MOVING BASE
141
The previous definition of 0l9 02, and 03 leads to the following expression for the angular velocity of dB in A:
(16) The angular velocity of dB in N, found by using the angular velocity addition theorem,34 is thus given by
(4,12,16)
. 1  W352 + Bl  03S2\ b±
Fig. 2 Deformed beam with differential element.
nOC(58)x x
,
,
,
,
(10)
The orientation of dB in A can be described by reference to a direction cosine matrix ACdB such that the element C/7 in the ith row and y'th column is defined as
] b3
(17)
C* =Qi'b
(/, y = l,2,3)
(11)
The matrix ACdB9 referred to as a spacethree, 123 rotation matrix, is given by33
A
The quantities iil9 u2, u3, 6l9 02, and 03 appearing in Eqs. (15) and (17) ultimately must be expressed in terms of system generalized coordinates. With this in mind, we introduce f as
$±p°Eg
an(
(18)
CdB = c2s3 \ —s s2
sls2s3 + c3q Sjc
+ J 3 *il CiS2s3  c3sl cc
(12)
with c, and 5, defined as (/ = 1,2,3) (13)
i ^et s(x> 0 denote the stretch in the beam along the elastic axis; tnen tne distance (see Fig. 2) along the deformed elastic axis; then the distance (see Fig. 2) along the axis from O to point E, which point is at a distance x from O when B is undeformed, is
The velocity of C in N is found by using the relationship
(19)
Now it is convenient to define J(a,t) as
where Avc9 the velocity of C in A9 can be formed by differentiating the position vector poc with respect to time in reference frame A. Substituting the resulting expression for Ayc, along with the expressions for Ny°9 NuA, and poc from Eqs. (3), (4), and (10), respectively, into Eq. (14) yields
e3clC2]
<o3[ w 2
S1S2S3
(20)
whereupon differentiation of Eq. (19) with respect to time leads to
, 01* (21)
and solution of Eq. (21) for ? produces
rf { y2 + co3 [ x + u± + e2 ( s^^
(22)
+ e3(c1s2c3
From Eqs. (57) and (18), it is apparent that
? = X + ?!
(23)
+ e3 [ (sls2s3 + c3q) ^
+
Hence,
(23) .(22)
V + 0)
W f ^
5
^
e2
a3
(15)
(24)
142
KANE, RYAN, AND BANERJEE
J. GUIDANCE
If the quantities s(x, t\ u2(x, /), u3(x,t), O^x, /), 02(x, /), and 03(x, t) are expressed! as
(25) 72,3) (26)
1=1
il
= 1,2,3) (27)
il
where $tj (i = l,...,6; j = l,...,v) are as yet totally unrestricted spatial functions, qi (i = 1,..., P) are generalized coordinates, and v is an arbitrary integer, then the term dJ(a, t)/dt which appears in Eq. (24) can be written as
2 E
il
iiyi
(28)
where a prime denotes differentiation with respect to the dummy variable a. Substituting this expression along with the one for s(x, t) from Eq. (25) into Eq. (24) yields
4 J3
+
\
v3 +col[u2\e2(sls2s3 +c^) + e3(cls2s3c3sl)]
?lV
?f
(29)
/?I
1=1
il
while substituting (31) jnto ?q (17)
* Us
'(32)
Expressions for the remaining time derivatives appearing in Eqs. (15) and (17) follow directly from Eqs. (2627). Specifically, (.7 = 2,3) (30)
il
the expressions for 0'1? 02, and 03 from Eq.
7 = 1,2,3) (31)
il
il
and substitution of the expressions for iilt w 2 , i/3, ^ j 02, and ^3 from Eqs. (2931) into Eq. (15) results in
il
/I
I, [il
(33)
Terms referred to as partial velocities (Ref. 34, p. 45) of C in TV, and partial angular velocities of dB in N will be denoted by the symbols Nyf and NutfB (/ = 1,..., P), respectively, and will play a major role in the subsequent formulation of equations of motion. The partial velocity Nyf is defined as the coefficient of qt (/ = !,...,?>) in the velocity expression shown in Eq. (32), and the partial angular velocity N B o)f is similarly defined as the coefficient of qt (z = 1,..., P) in the angular velocity expression given in Eq. (33). Consequently, these partial velocities and partial angular velocities
il
1JThe assumed mode method is employed here with a concomitant assumption that the displacement functions can be separated into products of spatial and temporal functions.
MARCHAPRIL 1987
CANTILEVER BEAM ATTACHED TO A MOVING BASE
143
can be written
An expression for u^ is formed as follows. From Eqs. (19) and (20), write
(38) Linearize both sides of Eq. (38), which yields
x + s(x,t)  f = X + MX
' (37) (23)
(39)
so that
_ (39)
Ul
= s(X)t) ^(25) =
(40)
Define two new symbols, )8/y and y/7, as
3i( a )$3/( a ) ^a 0*>/ = ! > ? ? ? > *0 (41)
(34)
and express the linearized partial velocity of C in TV as
. (34,3641) /
and
(35)
Since the goal of this analysis is to produce equations governing motions of B in N during which s(x, t\ u2(x9 1), B2(x,t\ and 03(.x, 0 remain small, it is 3(x9t),appropriate at this point in the analysis to linearize all kinematic quantities in the variables ' ql9 . . . , ^, ql9 . . . , qv. When linearizations are performed subsequent to the formation of expressions for partial velocities and partial angular velocities, as is being done here, the resulting dynamical equations are identical to those that would be obtained if linearizations were performed after completely nonlinear equations have been formulated (Ref. 33, p. 277). The great advantage of linearizing at the present stage of the analysis is that this relieves one of the considerable burden of developing complicated nonlinear expressions for translational and angular accelerations. It is tempting to linearize still earlier, that is, prior to the formation of expressions for partial velocities and partial angular velocities; but such linearization would be premature and would lead to dynamical equations lacking certain "dynamic stiffness" terms. In the following analysis, we indicate that a quantity has been linearized in the generalized coordinates and their time derivatives by placing a tilde over the associated symbol.
71
y1
71
E [
(42)
The linearized partial angular velocity of dB in N is
B(3541)
When linearizing velocity and partial velocity expressions, it is beneficial to note that
/(*,/)> i+E E
(36)
(20,26)
"
v
from which it follows that
(37)
The trigonometric quantities st and cf'(i = 1,2,3) defined in Eq. (13) can be linearized by expanding in Maclaurin's series and dropping terms of degree 2 and higher in the variables
/=!
(43)
144
KANE, RYAN, AND BANERJEE
J. GUIDANCE
while the linearized velocity of C in N and linearized angular velocity of dB in N can be expressed as
E
where 7n, 722, 733, expressed in terms of the central principal second moments of area 72 and 73, are given by
T> sln
J
33T, si.,.
III ^IT1I*+ 1*3 A
(49)
Since only internal forces contribute to the generalized active force F,, an expression for F/ can be derived from the
+ E {4ik(x)e2hk(x)+.e&sk(x)}qk U
strain energy function U by taking advantage of the fact that
F,where U is given by
d\J
(50)
A0
dx
03
V3
(44)
and
and the symbols Px, P^, P^, Tl9 M2i and M3 are defined as follows. Replace the set S of all forces acting on the particular cross section of B which is located at a distance * from O when B is undeformed with a couple of torque T together with a force R whose line of action passes through the elastic center E of the cross section. Then R can be resolved into an axial component, P^a^ and a shear component, V2Q2 + V3a3. Similarly, the torque T can be resolved into an axial torque, /i^!, and a bending moment, M2 a2 + M3a3. All other symbols in Eq. (51), with the exception of /c have been defined previously. The symbol k denotes an effective torsional constant3536 which would permit the sanie relative twist between beam ends for a given torque if one used the firstorder relationship (52)
+ E'
co 2 +
k=l
as would result from the use of the exact thirdorder equation for torsion,37'38 namely,
(53)
(45)
To produce linearized equations of motion, the expressions for "jyf, *<$B, Njf9 and NudB from Eqs. (4245) will be used in conjunction with the relationship (Ref. 34, p. 158)
where JF.* and ^ are the ith linearized generaHzed inertia force and ith linearized generalized active force for the system, respectively. The linearized generalized inertia force Ft* can be written
F*   f V§f 'NSC dx  {LN$B  ( V ?/ v y0 A) (47)
where K is the Saint Venant torsional factor for pure torsion arid T is the warping factor introduced earlier. The effective torsional constant K depends on the degree of fixity and warping at each end of the beam and can be calculated with the aid of Eqs. (52) and (53) in conjunction with tables in Ref. 32 (pp. 301312), or by utilizing publicly available computer programs (SDRC FRAME,39 SPOTS340). When effects of shear are included in the analysis, 02 and 03 are related to u3 and w 2 , respectively, as follows:
dx
(54)
and
du2(x9't) dx
(55)
where pdx is the mass of dB,Idx represents the central inertia dyadic of dB, Nac denotes the linearized acceleration of C in N9 Na.dB is the linearized angular acceleration of dB in N9 and all other symbols have been defined previously. The central inertia dyadic can be represented as
where the symbols i//2 and i//3 denote angles of shear of a cross section, measured at the elastic axis and associated with a2 and 0 3 , respectively. With these relationships, the shearing force measure numbers V2 and V3 appearing in Eq. (51) can be expressed as41
'
(56)
MARCHAPRIL 1987
CANTILEVER BEAM ATTACHED TO A MOVING BASE
145
A0G ,
A0G\du3(x,i)
?+?,
(57)
The other force and moment measure numbers present in Eq. (5i) are given by
PEA *1  t<0A dx
(59)
together with the ones for Nfg, NZ?B9 "&dB9 and /, from Eqs. (42), (43), (45), and (48) into Eq. (47), and equating to zero the sum [see Eq. (46)] of the resulting generalized inertia force and the generalized active force given in Eq. (62). These operations lead to a set of linear dynamical equations that can be represented in matrix form as
Mq +
(63)
and substitution of these expressions for Pl9 V29 V3, Tl9 M29 and M3 into Eq. (51) yields
where M, G9 and K are a mass matrix, skewsymmetric gyroscopic matrix, and stiffness matrix, respectively, all of dimension v X v\ ij, q, and q are v X 1 column matrices having elements qi9 qi9 and qh respectively, in the zth row; and F is a v X 1 column matrix. (If viscous "modal damping" is desired, it can be accounted for in thei usual manner by adding a term Cq to the above equation, with C representing a symmetric matrix of suitable damping coefficients.) The four matrices Af, G9 K9 and F are given explicitly in the next section.
IV. Simulation Algorithm
What follows is a description of a stepbystep {procedure for simulating motions of the beam system described earlier, and thus deterriiining values of s(x91)9 u2(x91)9 u3(x91)9 and B l ( x 9 1 ) for t > 0, once all system parameters, base motion temporal functions, and the initial beam displacement and velocity spatial functions have been specified. 1) Specify beam length L, material properties EQ(x)9 G(x)9 p(x) and crosssectional properties A0(x)9 I2(x), Ii(x)> ?2C^X *3(x), k(x)9 e2(x)> and e3(x). 2) Specify base motion temporal functions %(*), w 2 (0>
(60)
M
dx
3) Specify initial displacement function g(x) and initial velocity function h(x) defined as
? 3 (r) } ^(0, v2(t), and v3(t).
If the expressions for 5, t/ 2 , 3> ^i> ^2> ^^ ^3 froni Eqs. (2527) are now substituted into Eq. (60), and the resulting expression for U, in which the generalized coordinates qi9..,,'q,, appear explicitly, is differentiated in the manner shown in Eq. (50), one finds
?2(x,6)
?3(*>°)
03(x,0)
.,..,; *3(*,o)
(64)
?  E
dx
(61)
The term in braces in the above equation will be denoted by H,,. Thus
~
4) Specify the number v of assumed modal functions to be employed in the simulation. 5) Solve the eigenproblem associated with the cantilever beam with the root fixed in a Newtonian reference frame. Determine eigenfunctions <j>tj(x) (i = 1,..., 6; j = 1,..., v) and eigenvalues Xy (j — 1,...., v). (Note: For uniform beams of constant symmetrical cross section, one may wish to neglect shear effects and solve the associated eigenproblem analytically to obtain the necessary eigenfunctions and eigenvalues. For more complex beams, one should resort to a numerical technique, such as the finite element method, which yields discrete eigenvectors. These eigenvectors can be used in their discrete form in the subsequent calculations, or they can be converted into continuous functions through the use of shape functions.) 6) Calculate modal initial values as
(65)
;  I ^
A
1,...;,)
(62)
and
At this point in the derivation, equations of motion could be written in literal form by differentiating Nvc and N&dB in TV to obtain Nqc and Nac9 substituting the resulting expressions,
(66>
146
KANE, RYAN, AND BANERJEE
J. GUIDANCE
C, where E,? 15 U.1C ?;iCllltULl 111 lilt tia matrix ?" defined as
Table 1 Simplifications of general equations
J HI 1\J\\ CU1U A. U
To Neglect
Set a ^ = 0, ./ = 1,..., y
"l 0 0 0 0 0
0 1 0
0 0 1
e2
0
0 — e3
2
0 0 0
I
0 /22 + ?32 e2??3
Q 0 0
,
~i 0
e
r^ 1 4. ^2 + e
0
*2*3
e
0
Extension Bending in a2dir. Bending in 03dir. Torsion Rotatory inertia Shear Warping restraint Nonuniform torsion
^3i = 0, / = !,...,?' ^>4i = 0, / = 1,...,^
722 = 7 3 3=0
?2 = ?3 = 0
r =o
K =K
0
0
0
h^ \_ (67)
.?'Scittiii^'^OXiil, . . . ,6) implies setting WkHj = XkHj = Yk ijzkii. = Wki. = Xki = y = Z^, = 0 (k or / = n; j, j = 1,.. ,',V), as well as setting"2 //; associated primarily with direction n equal to zero. [See Eqs. (6875) and (78).]
and gk and hk are the elements in the &th row of the
respective column matrices g and h introduced in step 3. 7) Evaluate modal integrals:
(68)
/ (*)*//*)<**
l,..:, ?;*,/!,.. .,6) (69) (I'l,.'..,*; *!,... ,6) (70)
(Note; Prior to evaluation of these integrals by numerical quadrature techniques, analytical integration by parts should be performed on the integrals appearing in E<qs. (7677) in order to bring them into a simpler form. Additionally, the quantity Htj in Eq. (78) usually can be evaluated by using modal orthogonality properties in such a way that integration becomes unnecessary. For instance, if the finite element method is employed in the eigensolution, and the resulting orthogonal eigenvectors are normalized with respect to the finite element mass matrix, then it can be shown that HfJ = XjSij, where 8^ is the Kronecker delta.) 8) Assemble system mass, gyroscopic, stiffness, and force matrices:**
dx
= l,...,v, JM1.....6) (71)
(il,...,i;*4,5,6) (72)
)d*
^44/7 + ^66/yJ
w34ij + W43lj + e2(W44iJ+ W66ij)]
7l,...,F;fc,/4,5,6)
(73)
+ e3[(W15IJ WWf Wwj+Wnj)
+ e3 ( WMij + Waij)] ~ *2*3 [ W56ij + W65iJ] }
(ii,...,i;*4,5,6) (74)
(79)
(75)
+ "2
13/7^31,7 + Z^ij ~ Z64//]
(76)
(77)
^j^W^j + Y54ij  Y45ij]
UIJ
+ W36iJ  Wni}  W63iJ)
 W26ij} + e2U2(W^j WMij)\
?x ( W43iJ  WWJ) + ?2 ( W53IJ  W3Sij)
GAM
( x) #j
MJ
 Wmj)\  e2e3 [ a>2 ( W45ij  W54iJ) (i,j = l,...,v) (80)
dx
(78)
**The purpose of underlining certain terms is to facilitate the discussion in Sec. V concerning discrepancies between the present algorithm and various "multibody dynamics" algorithms.
MARCHAPRIL 1987
CANTILEVER BEAM ATTACHED TO A MOVING BASE
j
147
+ W64iJ)
lt 3
o ( wl3IJ + w3lij)  ( <4 + <*l}w22ij
.7 + W32ij)  ( <*l + <ol) WMJ
W56iJ
( <b2  ?1
+ "2 "i ( Z4SiJ + Z54ij  745,^  Y54,7)
+ 2 Wl<o2 ( W^y + WMiJ) ~ Z46,7 
sIJ + W54ij)
+ <o2<o3 ( Y^j + Y65ij)  ( ?i  <o?) Z66,7
 ?i \W23iJW32ij + YMJ + y65,7]
*i(WUIJ ~ Wmj)  *2(W45ij W54ij)
+ &3(W?iJ  WMij) + o>2(W45iJ + WS4iJ)
+ "2 \Wl3iJW3liJ + Y46ij + YMiJ + Zmj  Z64iJ]
 Y54ij + z45,7 + z54,7]
O'.yi,. ..,?) (si)
(>! + U2V3 ~ ^V^Wu + (f>2 + tQ3t)! ~ (0^3) W2i
Z4,.  74i)
?l( y4, + Z 4( ) + ( W2 + W 3Wl ) ys/ + ( W3  W!W2) Z6/
+ <o2co3 ( WWi
J
+ W62ij)
( CO + M) IL + ( CJ3 + (OltQ2) 12, ~ ( (02 ~ 6>!CJ3) I3i
^2 [  ( "3  ?1?2) ?U' ( "1 + "3 ) ^2,
u>ie2) W4i
(w 2 w 1 w 3 )i 4 ,
~(v2 + u3vl  Wly3  e
 u3e2 ) W6i  W63ij)
+ *2(WuljWntj)
 *3 ( W2kij  W62ij) +
 ^2 ( W45ij + WMJ)
 w2w3 ? 2 ,  w
— wxe3) W4i

?3
w1w3e3 + o)2e3) WSi
(82)
lco2 ( W25iJ
+ W52ij}  ?2co3 ( WMIJ
9) Solve the following system of equations for q.
(83)
63
10) Calculate physical displacements s, u2, u3, and Ol by using Eqs. (2527). All the equations in the above algorithm are written in a
J + W65ij) j + WMij)
 W43ij)  <b2 ( WxtJ  W53ij)
way that facilitates reduction of the general equations to deal with restricted cases. For instance, if the beam to be analyzed has a symmetric cross section, all terms in M, G, K, and F that are preceded by the scalar eccentricities e2 and e3 may be neglected. Additionally, if rotatory inertia is to be neglected, all modal integrals involving 722 and /33 will be zero, and system matrix terms including these integrals then need not be evaluated. Hence, by setting appropriate terms to zero in the complete equations provided here and subsequently
148
KANE, RYAN, AND BANERJEE
J. GUIDANCE
performing numerical simulations, one can recover the vast majority of previously published simulation results for beams. Table 1 shows how one proceeds to eliminate various beam effects. Note that, for purposes of clarity, the above equations and modal integrals are developed for a beam with constant eccentricity values. When the eccentricities are variable, they are simply included in the integrands of the modal integrals appearing in Eqs. (6876).
10 modes of the eigensolution were used to obtain these results.)
VI. Conventional Approach
The purpose of this section is to bring to light certain deficiencies of various existing multibody dynamic simulation computer programs, such as DISCOS,26 NBOD2,27 ALLFLEX,29 and TREETOPS,30 and to assess the significance of these deficiencies in connection with practical problems. Toward these ends, we first examine a fundamental issue
V. Illustrative Example
The generality of the algorithm presented in Sec. IV permits one to obtain solutions to a wide range of problems that span diverse disciplines. One representative problem involving a beam subjected to general base motion is that of simulating the behavior of a spacebased robotic manipulator, such as the one shown in Fig. 3, which consists of three links Ll9 L2, and L3 , connected by revolute joints. The outboard link L3 of the manipulator consists of a base A and two distinct segments Bl and B2. Segment Bl is 2f meters long and has a symmetric box cross section, while B2 is a 5 1/3 meterlong channel. Both segments are made of a material for which E0 = 6.895 X 1010 N/m2, G = 2.6519 X1010 N/m2, and p/A0 = 2766.67 kg/m3. The section properties for Bl are AQ = 3.84 X 10~4 m% 72 = 1.50 X 1(T^ m4, 73 = 1.50 X 10~7 m4, ?2 = 2.09, ?3 = 2.09, K = 2.2 X lO'7 m4, F = 0, e2 = e3 = 0, while the corresponding properties for B2 are A0 = 7.3 X 10~5 m2, 72 = 4.8746 X 10~9 m4, 73 = 8.2181 X10 ~ 9 m = 3.174, a3 = 1.52, /c = 2.4330 X 10~n m4, T = 5.0156 X 10~ 13 m6, e2 = 0, and e3 = 0.01875 m. To demonstrate the use of the simulation algorithm, we examine the behavior of L3 during deployment of the manipulator from a stowed configuration to a fully operational configuration. The deployment process is presumed to last for 15 seconds (T = 15), during which time the values of the angles <frl9 ^2, and ^3, shown in Fig. 3, change from 180 to 90 deg, 180 to 45 deg, and 180 to 90 deg, respectively. If the two inboard links L± and L2, each of length *?, are treated as rigid, this maneuver results in a prescribed motion of the base A of the outboard link, a motion characterized by the following six temporal functions:
Vl = <
+ C2)
Fig. 3 Spacebased robotic manipulator.
V2 =
V = 
(84)
0.2
(85)
where cz. = We take /= 8 m, and let (86)
0.1 10
(deg)
 
IT (

T . 2irt\ , T — rad
n
0< / < T
10
15 Time (sec)
 rad
T
IT
t>T
. 2vt\
(87)
Fig. 4 Manipulator deployment simulation results.
(88)
rad
4
t>T
2irt
(89)
k
Orad
t>T
Figure 4 shows, for the tip of B2, the transverse displacements u2 and w 3 , and the twist Ol during a time interval of 30 s, which corresponds to twice the deployment time. (The first
Fig. 5 Generic body of multibody system.
MARCHAPRIL 1987
CANTILEVER BEAM ATTACHED TO A MOVING BASE
149
that arises when one seeks to formulate equations to be used for the simulation of large overall motions of systems containing deformable bodies that experience only "small" deformations. Next, we comment briefly on matters of concern solely in connection with beams. Finally, we present and discuss numerical results for both a repositioning maneuver and a spinup maneuver of a flexible beam attached to a rigid base. The system consisting of a flexible beam attached to a rigid base constitutes one of the simplest systems within the application range of multibody programs and thus serves here as a " model" problem with which to check multibody dynamics formalisms. While the derivations of the equations of motion underlying the aforementioned multibody programs differ from each other in many details, they all rely on a common approach to characterize "small" deformations of individual flexible bodies. This conventional approach, as we shall call it, involves developing an expression for the velocity of a generic particle, say C, of a deformable body B (Fig. 5) by writing
V7 =Nv° + V X (r + u) +Avc
(90)
where AT is a Newtonian reference frame; O is the point at which B is connected to another body of the system; A is a reference frame in which a differential element of B containing O is fixed; Nv° is the velocity of O in AT, V4 is the angular velocity of A in AT, r is the position vector from O to C, that point of A with which C coincides when B is undeformed; u is the displacement vector extending from C to C; and Avc is the velocity of C in A. Once a set of unit vectors ql9 q2, q3 fixed in A has been introduced, the vectors r and u are written
r — x^ + x2q2 + x3q3
u = ?!?! + U2q2 + U3q3
(91)
(92)
and ul9 u2, and u3 are expressed as [see Eq. (1124) in Ref. 26, Eq. (43) in Ref. 27, Eqs. (6) and (8) in Ref. 29, and Eq. (2) in Ref. 30] 01,2,3) 71 where v is the number of terms in the series expansion, $ij(xl, x2, x3) (i = 1,2,3; 7 = 1, ...,*>) are independent spatial (modal) functions, and q}(t) (y = 1, ...,*>) are generalized coordinates. It follows immediately that
V
involving fjitj and TJO. We shall discuss the practical significance of this difference between the present theory and that corresponding to the conventional approach presently. First, however, we wish to point out additional differences between the two theories. In the formalisms underlying the aforementioned multibody programs [see Eqs. (II18a) and (1121) in Ref. 26, Eqs. (20) and (21) in Ref. 27, Eq. (22) in Ref. 29, and Eq. (7) in Ref. 30], as well as in the derivation underlying the present theory [see Eq. (47)], the development of expressions for inertia forces for a flexible body necessitates either an implied or an explicit integration over the space occupied by the body. When this space is threedimensional, it is appropriate for the integrand to involve solely quantities associated with translation of a generic point of the body. But when the space is onedimensional or twodimensional, as in connection with beams, plates, shells, etc., provisions must be made not only for translation of a generic point of a differential element, but also for rotation of the differential element. It is precisely for this reason that we introduce the angular velocity NudB in Eq. (17). The equation derivations in Refs. 26, 29, and 30 fail to incorporate such a term. If this term is absent also from a computer program based on one of these references, then rotational inertia effects associated with deformable structural elements are neglected in this program. For beams, this means that torsional effects are neglected entirely, which is a serious matter when one is dealing with beams connected to each other. Also, it precludes the accommodation of rotatory inertia effects, but this is relatively unimportant since often these effects are quite small. Additional differences between the present theory and others reflect refinements rather than conflicts. Specifically, provisions are made in the present theory to deal with shear and warping by employing the concepts of a shear area ratio and an effective torsional factor calculated with the aid of a warping factor; coupling of bending and torsion is made possible through the introduction of an eccentricity vector extending from the shear center to the centroid of a cross
(93)
0.5
*r E <M*i>*2,* 3 HO)
71
01,2,3)
(94)
so that, for / = 1, one has (95)
0.05
Now, this relationship, which is intended to be valid for a "general" structure, is in direct conflict with Eq. (29), developed specifically for a beam. The source of the conflict is the absence from the conventional approach of a counterpart to Eq. (19), which reflects an important physical fact, namely, that every transverse displacement of a point on the neutral axis of a beam (that is, a nonzero value of either u2 or u3) gives rise to an axial displacement (that is, a nonzero value of M X ). In other words, the conventional approach treats u± exactly like u2 and t/ 3 , which is precisely what is wrong with it. The failure of the conventional approach to deal adequately with the interdependence of displacements ultimately leads to a theory in which the counterpart to Eq. (81) lacks the terms
0.10
50 25
l(deg) (
10
15
Time (sec)
20
30
Fig. 6 Spinup maneuver.
150
KANE, RYAN, AND BANERJEE
J. GUIDANCE
section. These effects are not considered in the conventional approach. In summary, to see the difference in final equation form between the present theory and that corresponding to the conventional approach applied to a beam, one can refer to Eqs. (7982), where the doubly underlined terms represent the conventional theory. The remaining terms in these equations account for dynamic stiffness effects (p^ and rj /y ), torsion (Wuij, XMiJ, YMiJ, Zntj, and W4i, X4i, 74/, Z4/), rotatory inertia (Y55iJ and Z66/y), coupling of bending and torsion (e2 and e3 terms, as well as terms with indices 4,5 or 4,6), and coupling of differential element rotations B2 ^d ^3 [^ terms in Eqs. (81) and (82) with indices k, /= 5,6]. While many of the existing multibody programs may, in fact, incorporate the torsional terms mentioned above, none of these programs accommodate the dynamic stiffness effects that are crucial for certain simulations. While the comparison of equations is enlightening, it furnishes no indication regarding the quantitative significance of the terms one loses when taking the conventional approach. Therefore, we consider a numerical example which, it turns out, permits one to reach some conclusions regarding this matter. Returning to the manipulator system introduced earlier, we now examine the response of the tip of B2 when ^Er1, ^ ^3 (see Fig. 3) are given by (96) ^ ( 0 =  l 5 [ 2 + U  J l\ C ^ l ^ VJ^' l_ \ /
k (6*45)rad,
from those obtained by using an algorithm based on the conventional approach underlying various multibody dynamics programs. The present theory has been extended to incorporate "free" base motion as well; this work will be discussed in a forthcoming paper dealing with the attitude behavior of an elastic satellite subjected to gravitational and control forces. In an effort to develop a capability to simulate multibody systems correctly with beamlike components, the theory presented here is currently being incorporated into a new multibody computer program. We believe that recent advances in the use of simple continuum models in conjunction with modal parameter identification methods warrant such an extensive and comprehensive development activity. Certain features of the present derivation, such as the inclusion of differential element rotation, facilitate the process of connecting flexible bodies, and thus render this theory ideally suited for inclusion in such programs.
References
^odi, V.J., "Attitude Dynamics of Satellites with Flexible Appendages—A Brief Review," Journal of Spacecraft and Rockets, Vol. 11, Nov. 1974, pp. 743751. 2 Dubowsky, S. and Maatuk, J., "The Dynamic Analysis of Elastic Spatial Mechanisms," Proceedings of the Fourth World Congress on the Theory of Machines and Mechanisms, University of Newcastle upon Tyne, England, Sept. 812,1975, pp. 927932. 3 Turcic, D.A. and Midha, A., "Generalized Equations of Motion for the Dynamic Analysis of Elastic Mechanism Systems," Journal of Dynamic Systems, Measurement, and Control, Vol. 106, Dec. 1984, pp. 243248. 4 Turcic, D.A. and Midha, A., "Dynamic Analysis of Elastic Mechanism Systems. Part I: Applications," Journal of Dynamic Systems, Measurement, and Control, Vol. 106, Dec. 1984, pp. 249254. 5 Turcic, D.A., Midha, A., and Bosnik, J.R., "Dynamic Analysis of Elastic Mechanism Systems. Part II: Experimental Results," Journal of Dynamic Systems, Measurement, and Control, Vol. 106, Dec. 1984, pp. 255260. 6 Luh, J.Y.S., Walker, M.W., and Paul, R.P.C, "OnLine Computational Scheme for Mechanical Manipulators," Journal of Dynamic Systems, Measurements, and Control, Vol. 102, June 1980, pp. 6976. 7 Huston, R.L., "MultiBody Dynamics Including the Effects of Flexibility and Compliance," Computers and Structures, Vol. 14, No. 56, 1981, pp. 443451. 8 Geradin, M., Robert, G., and Bernardin, C, "Dynamic Modelling of Manipulators with Flexible Members," Advanced Software in Robotics Seminar, Liege, Belgium, Proceedings published by NorthHolland, 1983. 9 Singh, R.P. and Likins, P.W., "Manipulator Interactive Design with Interconnected Flexible Elements," Automatic Control Conference, San Francisco, CA, June 1983. 10 Gupta, R.S. and Rao, S.S., "Finite Element Eigenvalue Analysis of Tapered and Twisted Timoshenko Beams," Journal of Sound and Vibration, Vol. 56, No. 2, 1978, pp. 187200. 11 Putter, S. and Manor, H., "Natural Frequencies of Radial Rotating Beams," Journal of Sound and Vibration, Vol. 56, No. 2,1978, pp. 175185. 12 Kaza, K.R.V. and Kielb, R.E., "Effects of Warping and Pretwist on Torsional Vibration of Rotating Beams," Journal of Applied Mechanics, Vol. 51, Dec. 1984, pp. 913920. 13 Kaza, K.R.V. and Kvaternik, R.G., "Nonlinear FlapLagAxial Equations of a Rotating Beam," AIAA Journal, Vol. 15, June 1977, pp. 871874.
14 Hodges, D.H., "Proper Definition of Curvature in Nonlinear Beam Kinematics," AIAA Journal, Vol. 22, Dec. 1984, pp. 18251827. 15 Likins, P.W., "Dynamics and Control of Flexible Space Vehicles," JPL TR 321329, 1970. 16 Kulla, P., "Dynamics of Spinning Bodies Containing Elastic Rods, Journal of Spacecraft and Rockets, Vol. 9, April, 1972, pp. 246253. 17 Bodley, C.S. and Park, A.C., "The Influence of Structural Flexibility on the Dynamic Response of Spinning Spacecraft," AIAA Paper 72348, 1972. 18 Likins, P.W., Barbera, F.J., and Baddeley, V., "Mathematical Modeling of Spinning Elastic Bodies for Modal Analysis," AIAA Journal, Vol. 14, Sept. 1973, pp. 12511258. 19 Likins, P.W., "Geometric Stiffness Characteristics of a Rotating
(9?)
0<?<15s t> 15 s
(98)
During the corresponding manipulator motion, which we call a spinup maneuver, link Lt remains at rest, L2 and A remain at rest relative to each other, and the angular speed of A, namely, 4^, increases from zero to 6 rad/s in a time interval of 15 s, and remains constant thereafter. (The manner in which L^ and L2 are shown to be attached to each other precludes this sort of motion, but this is irrelevant for present purposes.) Such a motion is of practical interest in connection with the operation of helicopter rotor blades, aircraft propellors, satellites with flexible appendages, etc. The behavior of link L3 during this maneuver was predicted using two Fortran computer simulation programs. One program is based on the present theory and includes the full equations given earlier as Eqs. (7983). The second program is based on the conventional approach discussed above; in other words, it includes only those terms in Eqs. (7982) which are doubly underlined. For the particular spinup maneuver mentioned, the conventional approach leads to the dashed curves for u2, u3, and Ol in Fig. 6, which indicate a cataclysmic instability of the motion under consideration, whereas use of our theory produces the solid curves in the same figure, which show just the opposite, namely, an entirely wellbehaved response. On physical grounds, the latter is precisely what is to be expected, for we are dealing here with a simple case of "centrifugal stiffening." Thus one is forced to conclude that the conventional approach can lead to totally misleading predictions of the behavior of elastic bodies in certain situations.
VII. Conclusion
A comprehensive theory for dealing with small vibrations of a general beam attached to a base that is undergoing an arbitrary, prescribed motion has been presented. Numerical results show that predictions of this theory can differ markedly
MARCHAPRIL 1987
CANTILEVER BEAM ATTACHED TO A MOVING BASE
151
Elastic Appendage," International Journal of Solids and Structures, Vol. 10, 1974, pp. 161167. 20 Vigneron, F.R., "Comment on Mathematical Modeling of Spinning Elastic Bodies for Modal Analysis," AIAA Journal, Vol. 13, Jan. 1975, p p . 126127. . ' ' . . ' ? ?1 Hablani, H.B., "Constrained and Unconstrained Modes; Some Modeling Aspects of Flexible Spacecraft," Journal of Guidance, Control, and Dynamics, Vol. 5, MarchApril 1982, pp. 164173, 22 Meirovitch, L., Analytical Methods in Vibrations, Macmillan, New York, 1967, pp. 443445'.' 23 Hughes, P.C., "Dynamics of Flexible Spinning Spacecraft," XXIV Congress of the International Astronautical Federation, Baku, USSR, Oct. 1973. 24 Benedetti, G.A., "Derivation of the Coupled Equations of.Motion for a Beam Subjected to Three Translational Accelerations and Three Rotational Accelerations," Sandia National Laboratories, Rept. SAND798202, Feb. 1979. 25 Laskin, R.A., Likins, P.W., and Longman, R.W., "Dynamical Equations of a FreeFree Beam Subject to Large Overall Motions," Journal of the Astronautical Society, Vol. XXXI, Oct.Dec. 1983, pp. 507528. 26 Bodley, C.S., Devers, A.D., Park, A.C., and Frisch, H,P, "A Digital Computer Program for the Dynamic Interaction Simulation of Controls and Structure (DISCOS)," Vols. 1&2, NASA TP1219, May 1978. 27 Frisch, H.P., "A VectorDyadic Development of the Equations of Motion for ATCoupled Flexible Bodies and Point Masses," NASA TN D8047, Aug. 1975. 28 Ho, J.Y.L., "Direct Path Method for Flexible Multibody Spacecraft Dynamics," Journal of Spacecraft and Rockets, Vol. 14,
Jan.Feb. 1977, pp. 102110. 29 Ho, J.Y.L. and Herber, D.R., "Development of Dynamics and Control Simulation of Large Flexible Space Systems," Journal of Guidance, Control, and Dynamics, Vol. 8, MayJune 1985, pp. 374383. 30 Singh, R.P., VanderVoort, R.J., and Likins, P.W., "Dynamics of Flexible Bodies in Tree Topology—A Computer Oriented Approach," AIAA Paper 841024,1984. 31 Oden, J.T. and Ripperger, E.A., Mechanics of Elastic Structures, 2nd ed., McGrawHill, New York, 1981. 32 Roark, R.J. and Young* W.C., Formulas for Stress and Strain, 5th ed., McGrawHill, New York, 1975, pp. 512. 33 Kane, T.R., Likins, P.W., and Levinson, D.A., Spacecraft Dynamics, McGrawHill, New York, 1983, p. 425. 34 Kane, T.R and Levinson, D.A., Dynamics, Theory and Applications, McGrawHill, New York, 1985, p. 24. 35 Hay, J.K. and Blew, J.M., "Dynamic Testing and Computer Analysis of Automotive Frames," SAE Paper 720046, Jan. 1972. 36 McClelland, W.A, Hay, J.K., and Klosterman, A.L., "Frame Design Analysis Under Complete Vehicle Boundary Conditions," SAE Paper 741142, Nov. 1974. 37 Timoshenko, S., Strength of Materials, 3rd ed., Vol. 2, Van Nostrand, Princeton, NJ, 1955, p. 265. 38 Heins, C.P., Bending and Torsional Design in Structural Members, Lexington Books, Lexington, MA, 1976, pp. 824, 65113. 39 SDRC Frame User's Manual, Vol. 1, Structural Dynamics Research Corp., 1980. 40 SDRC MechanicalDesign Library'User's Manual, Vol. 1, Structural Dynamics Research Corp., 1978. 41 Timoshenko, S.P., Young, D.H., and Weaver, W., Vibration Problems in Engineering, 4th ed., Wiley, New York, 1974, pp. 432435.
From the AIAA Progress in Astronautics and Aeronautics Series... AEROTHERMODYNAMICS AND PLANETARY ENTRY—v. 77 HEAT TRANSFER AND THERMAL CONTROL—v. 78
Edited by A. L. Crosbie, University of MissouriRolla The success of a flight into space rests on the success of the vehicle designer in maintaining a proper degree of thermal balance within the vehicle or thermal protection of the outer structure of the vehicle, as it encounters various remote and hostile environments. This thermal requirement applies to Earthsatellites, planetary spacecraft, entry vehicles, rocket nose cones, and in a very spectacular way, to the U.S. Space Shuttle, with its thermal protection system of tens of thousands of tiles fastened to its vulnerable external surfaces. Although the relevant technology might simply be called heattransfer engineering, the advanced (and still advancing) character of the problems that have to be solved and the consequent need to resort to basic physics and basic fluid mechanics have prompted the practitioners of the field to call it thermophysics. It is the expectation of the editors arid the authors of these volumes that the various sections therefore will be pf interest to physicists, materials specialists, fluid dynamicists, and spacecraft engineers, as well as to heattransfer engineers. Volume 77 is devoted to three main topics, Aerothermodynamics, Thermal Protection, and Planetary Entry. Volume 78 is devoted to Radiation Heat Transfer, Conduction Heat Transfer, Heat Pipes, and Thermal Control. In a broad sense, the former volume deals with the external situation between the spacecraft and its environment, whereas the latter volume deals mainly with the thermal processes occurring within the spacecraft that affect its temperature distribution. Both volumes bring forth new information and new theoretical treatments not previously published in book or journal literature.
Published in 1981, Volume 77—444 pp., 6x9, illus., $35.00 Mem., $55.00 List Volume 78—538 pp., 6x9, illus., $35.0(9 Mem., $55.00 List
TO ORDER WRITE: Publications De.pt.. AIAA. 1633 Broadway. New York. N.Y. 10019