- Email: [email protected]
Computational methods in constrained multibody dynamics: Matrix formalisms
Computers & Srrucrures Vol. 29. No. 2, pp. 331-338. Printed in Great Britain.
1988
0045-7949/88 $3.00 + 0.00 Pergamon Press plc
COMPUTATIONAL METHODS IN CONSTRAINED MULTIBODY DYNAMICS: MATRIX FORMALISMS J. T. WANG? and R. L. HUSTON~ tEngineering Mechanics Department, General Motors Research Laboratory, Warren, MI 48090, U.S.A. $,Department of Mechanical and Industrial Engineering, University of Cincinnati, Cincinnati, OH 45221, U.S.A. (Receiued 6 July 1987)
matrix-based analysis of multibody dynamics is presented. The analysis is applicable to large systems subjected to both geometric and kinematic constraints. The matrix formulation leads directly to algorithms for computational procedures. It also presents the analysis in a format where the dynamical
Abstract-A
behavior of the system is evident. An example of a towing chain is presented as an illustration.
1. INTRODUCTION
Recently there has been a marked increase in dynamical analyses of multibody systems. This interest stems from the ability to model a broad class of physical systems by multibody systems. The use of multibody system methods in dynamical analyses is analogous to the use of finite element methods in structural analyses. A multibody system may be defined as a collection of connected bodies such that there is at least one degree of freedom at the connecting joint. Figure 1 depicts a typical multibody system. The system may be internally constrained with closed loops and motion restricting joints, leading to geometric (or holonomic) constraint equations. The system also may be externally constrained by specifying the motion of some of the bodies, leading to kinematic (or non-holonomic) constraints. Multibody systems have been used to model manipulators, robots, cables, beams, and biosystems. References [l-lo] document a few of the applications of multibody systems methods. The successful application of these methods is based upon two factors: (1) advances in computational procedures and (2) advances in system modelling and analysis. References [1l-301 document some of these advances. In this paper we present an analysis which we believe provides a basis for an efficient computational procedure. The analysis employs a matrix formalism which readily leads to algorithms for computation. The paper is divided into five parts with the following part discussing some concepts useful in the sequel. The next two parts present the dynamic analysis and the governing equations. The fourth part presents an illustrative example and the final part contains a discussion and concluding remarks.
center of Bk. Let R be an inertial reference frame in which S moves. Let S have N bodies and n degrees of freedom characterized n generalized coordinates q, (r = 1, . . . , n). Then, Huston et al. [15, 161 have shown that the velocity of Gk and the angular velocity of Bk in R may be expressed in the forms: “k =
vmlk !h %a
and
Ok
=
%dk
Ch’btn
9
(1)
where h,,, (m = 1,2,3) are unit vectors fixed in R. The coefficients v,,,,~and w,,,[k(m = 1,2,3; 1 = 1, . . . , n; N) are components of the partial velocity k=l,..., and partial angular velocity vectors as defined by Kane et al. [15, 16,20,3 11. They may be viewed as ‘base vectors’ in the n-dimensional space of the 4,, in which S moves. (Regarding notation, in eqn (1) there is a sum on repeated indices over the range of the index.) In matrix notation, let the 3 x n array ok be defined as:
vk =
01Ik
v12k
U2lk
v22k
L?lk
%k . . . V2”k
hk
%,k
=
bmnkl.
(2)
I
Then v, may be expressed in matrix form as: vkl “k =
ok2
=
(3)
okq,
.
[Ivk3
where 4 is the array whose transpose dT is 4r= [4,,42,.
. .,61.
(4)
Next, let v be the set of N velocity arrays each having the form of eqn (3). Then v may be written as
2. PRELIMINARY CONSIDERATIONS
v = V&
(5)
2.1 Partial velocities and partial angular velocities
Consider again the multibody system S of Fig. 1. Let Bk be a typical body of S and let Gk be the mass
where V is the N x 3 x n block array of partial velocities.
331
332
J. T. WANGand R. L.
HUSTON
Fig. 1. A typical multibody system.
Similarly, the angular velocities may be expressed as:
and SJK SKL = SJL
o,=R,g
or
o=Qq,
where I is the 3 x 3 identity matrix and L refers to a third body B,. By repeated use of the last of eqn (1 1), we can obtain a transfo~ation matrix between Bk and the inertia frame R (designated by 0), i.e.
where the partial angular velocity arrays !& and Q are defined in the same manner as V, and V are defined. 2.2 Transformation matrices The relative orientations of the bodies may be defined in terms of orientation angles (such as Euler angles). The generalized coordinates q, may be identified with these angles. Transformation matrices relating unit vectors fixed in adjoining bodies then may also be defined in terms of these angles. Specifically, let Bj and Bk be typical adjoining bodies of the system and let n,i and nki (i = 1,2,3) be mutually perpendicular unit vector sets fixed in Bj and B,. Then let the 3 x 3 transformation matrix SJK be de&red as: SJK, = njp-nkq.
SOK=SOl
SOK, = qp.nkq.
-
[
sk2
ck3
sklsk2
sklck2ck3- sk3cki
sk3 skl
-cklckZsk3
-
Ck3skl
-sklck2sk3
-
ck3ckl
nk = SJKTnj.
(10)
It is readily seen that the transformation matrices obey the following identity and transitive relations: SJK SKJ=SJK
(8)
I I
SO&%,
= e~~~~~~SOK~~, and
SJKr=SJK
SJK-‘=I
(14)
sk2sk3
ckl ck2 ck3 -
where ski and g, represent the sine and cosine of @ki. SJK may be used to relate the unit vector sets by Hence, the derivative of SOK, is the expressions: dSOK,,/dt = l$,‘Wk X nkq njP= SJK,n, and nk#= SJK,ul,, (9) = ~p’@h,,kcir%m X
nj = SJKn,
(13)
dn,,/df = ok x nkg (no sum on k).
cklsk2
or in matrix form as:
(12)
Equation (13) is useful for developing algorithms for computing the time derivative of SOK: since the n,, are fixed in 4 their derivatives in inertia frame R are:
(7)
ck2
S12,...,SJK,
where the indicted product involves the bodies of the branches of S containing Bk, From eqn (7) SOK is seen to have the elements
For example, if ‘1,3, 1’ Euler angles 6,, (i = 1,2,3) (see [32]) are used to orient B, relative to B,, then SJK is
SJK =
(11)
(6)
(IS)
where ep,,,,is the permutation symbol [33]. In matrix form eqn (15) may be written as SdK = WOK SOK,
(16)
where WOK is the matrix whose ‘dual vector’ 1333is
Computational methods in constrained multibody dynamics ok,
i.e.
0 =
-
wk3
0
wk3
[-
OkI
Ok2
“k2 -wkl
0
7
(17)
I
where oki (i = 1,2,3) are the II,,~components of Ok. Equations (14) and (16) are especially significant in the numerical analyses of multibody system dynamics. They provide multiplication algorithms for computing derivatives-operations ideally suited for digital computation.
333
popular option is the use of Euler parameters for this purpose: as before, let Bj and Bk be typical adjoining bodies of the system. Then at any instant the change in orientation of Bk relative to Bj may be represented as a rotation about a line fixed in each body. If 12, is a unit vector parallel to this line and if ok is the rotation angle, four Euler parameters cki (i=l,... ,4) may be defined as: 6ki= Lkisin(8,/2),
i = 1,2,3
and ck4 =
wek/2h
(20)
2.3 Generalized speeds and Euler parameters
In eqns (1) (5) and (6) the generalized coordinate derivations q, are the coefficients of the partial velocity and partial angular velocity vectors. If the partial velocity and partial angular velocity vectors are viewed as base vectors in the n-dimensional space in which S moves, then the B (r = 1, . . . , n) may be viewed as magnitudes or measures of the movement. On many occasions it is convenient to use independent linear combinations of the 4, as measures of the motion of S. Such linear combinations are called ‘generalized speeds’. If U, (r = 1,. . . , n) is a set of generalized speeds, the transformation between 4, and the U, may be expressed in the matrix form 4 = Au
or
u = A T4,
and
fi=QA.
2 -1 .
~:,+&+&+t,-
(21)
An advantage of using Euler parameters is that there is a linear relation between the Euler parameters and the relative angular velocity components. This avoids singularities which arise when orientation angles are used [14]. Specifically, if e, and djk are the arrays of the Euler parameters and the relative angular velocity components, defined as
(18)
where A is an orthogonal matrix. By substituting from eqn (18) we can define partial velocity and partial angular velocity arrays B and fi associated with the generalized speeds. Specifically, P=VA
where the Akiare the components of & relative to unit vectors nji fixed in Bj. The Euler parameters are not unique. They are related by the equation
(19)
A convenient set of generalized speeds for the system of Fig. 1 is the set of local components of the relative angular velocity vectors [15]. There is, however, a nonlinear relationship between the orientation angles and the relative angular velocity components. This means that the orientation angles cannot simply be differentiated to produce the relative angular velocity components. Indeed, in general there are no functions which can be differentiated to produce the relative angular velocity components. Hence, if we
use relative angular velocity components as generalized speeds, specific generalized coordinates do not exist. To overcome this problem it is convenient to use variables other than orientation angles to define the relative orientations of the bodies. An increasingly
and ~i)kT=[~kl,~k2r~kj,o],
G-9
then 4 and Ok are related by the expressions 4 = Eok
and
where the orthogonal
Ek3
-
ck4
-;;; -6kl
ck2
ekl
ckl
ck2
6k4
ck3
Ek3
&k4
. -Ekl
I
(23)
matrix E is Ek4
E = (l/2)
uik = ETik,
-
ck2
-
(24)
I
Finally, when Euler parameters are used the transformation matrix SJK has the form
2.4 Derivatives of the partial angular velocity arrays
velocity
and partial
An advantage of using relative angular velocity components as generalized speeds is that the resulting partial angular velocity components may be directly identified with the elements of the transformation
J. T. WANG and R. L. HUSTON
334
matrices. Hence, since we have an algorithm for the differentiation of the transformation matrices [eqn (1611 we can then readily differentiate the partial angular velocity array. To see this, consider again the angular velocity wk of a typical body Bk of the system. From the addition theorem for angular velocities we can express o, as o,=~,+fi*+“‘+&,+&
(26)
where, as before, the carats over the symbols designate the relative angular velocity of the subscript body with respect to the adjacent lower numbered body. The sum in eqn (26) is taken over the bodies in the branch of S containing Bk. Consider Ok: it may be expressed as:
In matrix notation ak and a, may be expressed as [see eqns (3) and (6)] ak = hktf + ti,ri
and
ak = &ikuf vkti,
(30)
where we have dropped the carats for notational simplicity. The four arrays @,,&,ti,,,[k,v,,,,~and d, , or alternatively, R, h, V and v, play a central role in the analysis of multibody dynamics. 3. DYNAMICS
To explore the systems acting on be separated into forces; constraint
dynamics of S consider the force the bodies. Let these force systems three categories: externally applied forces; and inertia forces.
3.1 Applied forces where as before the n,, are unit vectors fixed in adjoining body B,. Let the e& be identified with the generalized speeds up where p is 3(k - 1) + i. Also, from eqn (9) let the nji be expressed as SO~~i~~. Then by substituting into eqns (27) and (26), and by comparing with eqns (1) and (18), we see that the wmIk may be identified with the elements of SOJ and the transformation matrices of the other bodies in the branch of Bk. Hence, the derivatives of the o,,,,~ are obtained directly from eqn (16). Additional details are given in [15]. Finally, it is seen that the elements of the partial velocity array can also be expressed in terms of the elements of the transformation matrices. To see this, consider again a typical body Bk with mass center Gk. Let pk locate Gk relative to a fixed point 0 in inertia frame R. Then pR may be expressed in terms of vectors fixed in the bodies in the branch of Bk. Let rj be one of these vectors. Let 4, be fixed in B,. Then in computing the velocity of Bk it is necessary to differentiate cj in R. This derivative fj is simply +r+x
pj.
(28)
Let the applied (or ‘active’) forces on each body Bk force system consisting of a Single force Fk passing through the mass center Gk together with a couple having torque Mk. Then the contribution by Fk and Mk to the generalized active force F,, associated with the generalized speed u,, is [15,31] be replaced by an equivalent
= Fkn%k +
Mkm %nik 3
(31)
where Fkm and Mkm are the n, components of Fk and mk. Observe that in eqn (31) the terms represent the contribution to F, from Fk and M,. Hence, if there is a sum on k we have the ~nt~bution from all the bodies leading to the total F,. Observe further that if the partial velocity and partial angular velocity vectors are viewed as base vectors defining the directions in which S can move, then F, represents for u, the projections of F, and Mk along these directions. Equation (31) may be written in the matrix form:
(32) Then by eqns (1) tj and similar derivatives may be expressed in terms of the transfo~ation matrix where (Pkand & are the arrays of rbm components of elements. Hence, the velocity of G, and the elements Fk and Mk. of the partial velocity arrays may also be expressed in terms of the transfo~ation matrix elements. Then 3.2 Constraint forces the derivatives of the partial velocity array elements Constraints and thus constraint forces arise when may be computed using the algorithm of eqn (16). the movement of S is restricted. This may occur in Additional details are given in 1151. By knowing the derivatives of the o,,,~~and =?&,k several ways: (1) there are fewer than six degrees of arrays we can compute the angular acceleration ak freedom at a joint; (2) S may have one or more closed and mass center acceleration ak of the bodies in R, i.e. loops; or (3) the movement of some of the bodies of S relative to each other or relative to inertia space R may be specified. As with the applied forces, let the constraint forces on each body be replaced by an equivalent force and system consisting of a single force FL passing through (29) mass center Gk (or alternatively through joint center
Computational methods in constrained multibody dynamics 0,) together with a couple having torque ML (or Ma). Then the contribution by FL and Mi to the generalized constraint force F; associated with the generalized speed u, is of the form [see eqn (32)]: F; = V:d;
+ i-l,‘&.
(33)
where 4; and II; are the arrays of no,,,components of Fk and Mk. The restriction of the movement of S may be described by equations of the form Bijuj=gi,
i=l,...,
m; j=l,...,
n; m
or more generally the forms F+F’+F*
=0
or
F+B’l+F*=O,
(40)
where F, F’ and F* are the arrays of generalized applied, constraint and inertia forces. Finally, observe that by substituting from eqns (29) and (37), Fy may be expressed as
or as
(34)
or in matrix form as
335
F* = -Mti
+J
(41)
where the mb and f; are Bu =g
(35)
where B is the M x n ‘constraint assay’ with elements b,.
Each of the individual equations of eqns (34) may be associated with a single constraint force or moment component. Indeed, it has been shown that if 1 is the combined array of constraint force and moment components then Fj and B are related by the expression [35] B’A = F;.
(36)
(42) and fi = -
crnk %lk
+
The inertia forces may be treated in the same manner as the applied ‘and constraint forces. Specifically, let the inertia force system on each body Bk be replaced by a single force F! passing through Gk together with a couple having torque Mz . Then F: and Mt may be expressed as
Ih
%,lk
ersmlksn%~k
&“,?,k up “rqkonpk
uqup)T
(43)
where the Ik,,,,,are the no,,,and nr,” components of I, and e,, is the permutation symbol [33]. The array IV, called the ‘generalized mass array’, is seen to be symmetric and positive definite. 4. GOVERNING
3.3 Inertia forces
i?npk up +
EQUATIONS
AND SOLUTION
METHODS
From eqns (34), (40) and (41) we seen that the governing equations for constrained multibody systems can be expressed in the matrix forms: Mti=F+f+BTL
(44)
Bu =g.
(45)
and
and Mf = -Ik.ak -ok
x (Ikeok),
(37)
where mk and Ik are the mass and central inertia dyadic [31] of Bk. Then in view of eqns (32) and (33) the generalized inertia force F: associated with generalized speed a, may be expressed as F:=
V:4:+Q,Tpc,+,
(38)
where as before 4: and p: are the arrays of b,,, components of Ff and Mt.
Since M is an n x n array and since B is an m x n array, these equations are equivalent to a set of m + n scalar equations for the n generalized speeds u and the m constraint force and moment components 1. The equations thus form a coupled differential-algebraic system of equations. A convenient solution method is as follows. Suppose an n x (n - m) matrix C can be found such that BC is zero. Then C is an ‘orthogonal complement’ [12, 18, 191 of B. Hence, we can eliminate 1 be premultiplying eqn (44) by CT, i.e. C’Mti = C’F+ C’f:
3.4 Dynamical equations Using Kane’s method take the simple form
Then by differentiating the dynamical
F,+F;+F:=O,
(46)
eqn (45) we have
equations Bti=g-h. (3%
(47)
Equations (46) and (47) constitute a set of n equations
J. T. WANGand R. L.
336
-
4 hlsec’
(1.22 mlsec’)
i;
Fig. 2. Towing chain-initial
configuration.
for the n unknowns u. By solving for u we can then use eqn (44) to find 1. 5. EXAMPLE:
A TOWING
HUSTON
given an acceleration of 4.0 ft/sec2 (1.22 m/set*) to the right. This system is thus a constrained multibody system with specified motion of two points of the system. The response of the system may be obtained by formulating the dynamics in the form of eqns (44) and the constraints in the form of eqns (45). The solution of these equations using the procedures of eqns (46) and (47) produces the time interval response as in Fig. 3. In developing this solution we used the computer code DYNOCOMBS developed by Kamman and Huston 1341.
Equations (46) are equivalent to the standard form of Kane’s equations. They represent a projection of eqns (38) onto a direction orthogonal to the constraint force array. This can be depicted geometrically by the triangle shown in Fig. 4. In Fig. 4 K and K* are the generalized force arrays employed by Kane et al. ]27,31,32]. Specifically, K and K* are
K=CTF
and
K*=CTF*.
(48)
CHAIN
As a simple example illustrating the use of some of these ideas consider the towing chain depicted in Fig. 2. The chain consists of six identical pinconnected rods. They each have length 1.0 ft (0.3048 m) and mass 1.0 slug (14.6 kg). One end of the chain is attached to a fixed point (representing an object to be towed) and the other end is given an acceleration (simulating the towing vehicle). Initially, the two ends are placed together and the chain links hang in a vertical gravity field. The towing end is then
Fig. 4. Geometric interpretation of generalized force arrays.
,t=1.m
Fig. 3. Towing chain ~nfiguration
at successive time intervals.
Computational
methods in constrained multibody dynamics
putational assistance of Mr Timothy P. King and Mr Yung
governing equations are then
The reduced
337
Sheng Liu. K+K+=O.
(49)
Observe that the constraint forces do not contribute to these generalized forces or to the reduced governing equations. Hence, the constraint forces may be ignored in analyses using Kane’s equations. This is a principal advantage of Kane’s equations. Next, the constraint force and moment array 1 in eqns (36), (40) and (44) may be considered as Lagrange multipliers. Hence, eqns (40) may be considered as a generalization of Kane’s equations to include constraint forces-a procedure long used with Lagrange’s equations [36,371. Finally, eqns (44) and (45) are in a form which is convenient for multibody dynamics analyses. For example, we can formally eliminate 1. If we solve eqn (44) for zi we obtain ti = M_‘(F +f+
BTA).
(50)
Then by substituting into eqn (47) [the differentiated form of eqn (45)] we obtain BM-‘(F+f+B=A)=g
-Lb.
(51)
Solving for I we have 1 = (BM-‘P--‘[g
- lh - BM-‘(F
Finally, by substituting
zi= M-‘(F
+f)
+
+f)].
(52)
into eqn (50) we obtain
M-‘Br{(BM-‘BT)-’
[g - Lh - BW’(F
+f)]j.
The
term M-‘BT(BM-‘Br)-’ is the pseudoinverse’ [38,39] of B. That is, B[M-‘Br(BM-‘Br)-‘I=
Z,
(53)
‘weighted
(54)
where I is the m x m identity
matrix. The pseudoinverse has been employed by a number of investigators in studying inverse kinematic and inverse dynamic problems in robotics [38-41]. Finally, it is believed that the range of applications of the procedures exposited herein and the governing equations have yet to be developed. The advantage of the matrix formalism is its simplicity and its ready conversion to computer algorithms. This, coupled with the simplicity and eloquence of Kane’s equations, makes the governing equations very attractive for application to constrained, redundant multibody systems. Specific applications in robotics, biosystems and structural analysis are anticipated. Acknowledgements-The research for this paper was partially supported by the National Science Foundation under Grant MSM-8612970. This support is gratefully acknowledged. The authors also gratefully acknowledge the com-
REFERENCES
1. G. E. Fliescher, Multibody systems applications and simulations at the jet propulsion laboratory. In Dynamics of Multibody Systems (Edited by K. Magnus), pp. 3647. Springer, New York (1978). 2. J. G. Gallenstein and R. L. Huston, Analysis of swimming motions. Human Factors 15, 91-98 (1973). 3. J. C. Huston, M. W. Harlow and R. L. Huston, A comprehensive, three-dimensional head-neck model for impact and high acceleration studies. Aviation, Space Envir. Medicine 49, 205-210 (1978). 4. R. L. Huston, A summary of three-dimensional, grossmotion, crash-victim simulators. In Structural Mechanics Software Series, Vol. I, pp. 611622. University Press of Virginia, Charlotttesville, VA (1977). 5. R. L. Huston, R. E. Hessel and J. M. Winget, Dynamics of a crash victim-a finite segment model. AIAA Jnll4, 173-178 (1976). 6. R. L. Huston and C. E. Passerello, On the dynamics of a human body model. J. Biomech. 4, 369378 (1971). 7. J. W. Kamman and R. L. Huston. Modellina of submerged cable dynamics. Cotnput. Srruct. 20;623-629 (1985). 8. T. R. Kane and D. A. Levinson, The use of Kane’s dynamical equations in robotics. Int. J. Robotics Res. 2, 3-21 (1983). 9. F. A. Kelly and R. L. Huston, The development of equations of motion of single arm robots. IEEE Trans. Systems, Man, Cybernetics, SMC-12, 259-266 (1982). 10. J. M. Winget and R. L. Huston, Cable dynamics-a finite segment approach. Comput. Struct. 6, 475480 (1976). 11. Y. 0. Bayazitoglu and M. A. Chace, Methods of automated dynamic analysis of discrete mechanical systems. J. appl. Mech. 40, 809-819 (1973). 12. H. Hemami and F. C. Weimer, Modelling of nonholonomic dynamic systems with applications. J. appl. Mech. 48, 177-182 (1981). 13. W. W. Hooker, A set of r dynamical attitude equations for an arbitrary n-body satellite having r rotational degrees of freedom. AIAA Jnl8, 1205-1207 (1970). 14. R. L. Huston and C. E. Passerello, Eliminating singularities in governing equations of mechanical systems. Mech. Res. Commun. 3, 361-365 (1976). 15. R. L. Huston, C. E. Passerello and M. W. Harlow, Dynamics of multirigid-body systems. J. appl. Mech. 45, 889-894 (1978). 16. R. L. Huston and C. E. Passerello, On multi-rigid-body system dynamics. Comput. Struct. 19, 439446 (1979). 17. R. L. Huston and C. E. Passerello, Multibody structural dynamics including translation between the bodies. Comput. Srrucr. 11, 715-720 (1980). 18. J. W. Kamman and R. L. Huston, Dynamics of constrained multibody systems. J. uppl. Mech. 51, 899-904 (1984). 19. J. W. Kamman and R. L. Huston, Constrained multibody system dynamics. Comput. Srruct. 18, 999-1003 (1984). 20. T. R. Kane and D. A. Levinson, Formulation of equations of motion for complex spacecraft. J. Guidance Control 3, 99112 (1980). 21. P. W. Likins, Guasi-coordinate eauations for flexible spacecraft. AIAA Jnl 13, 524-526 (1975). 22. P. W. Likins, Dynamic analysis of a system of hingeconnected rigid bodies with nonrigid appendages. In?. J. Solids Struct. 9. 1473-1487 (1973). 23. P. W. Likins, Point-connected rigid bodies in a topological tree. Celestial Mech. 11, 301-317 (1975).
J. T. WANG and R. L. HUSTON
338
24. L. K. Lilov and V. A. Chirikov, On the dynamics equations of systems of interconnected bodies. PMM (J. anal. Mum. Mech.) 45. 383-390 (1982). 25. L. L& and M. Lorer; Dynamic analysis of multirigidbody system based on the Gaus principle. ZAMM 62, 539-545 (1982). 26. P. E. Nikravesh and I. S. Chung, Application of Euler parameters to the dynamic analysis of threedimensional constrained mechanical systems. J. Mech. Des. 104, 785-791 (1982). 21. C. E. Passerello and R. L. Huston, Another look at nonholonomic systems. J. uppl. Mech. 40, 101-104 (1973). 28. R. Schwertassek and R. E. Roberson, A state-space dynamical representation for multibody mechanical systems, Part I: systems with tree configurations. Actu Mechanica 50, 141-161 (1984). 29. R. Schwertassek and R. E. Roberson, A state-space dynamical representation for multibody mechanical systems, Part II: systems with closed loops. Actu Mechun icu 51, 15-29 (1984).
30. M. Vukobratovic, Computer method for dynamic model construction of active articulated mechanisms using kinetostatic approach. Me& Muchine Theory 13, 19-39 (1978).
31. T. R. Kane and D. A. Levinson, Dynamics: Theory and Applications, pp. 15, 45. McGraw-Hill (1985).
32. T. R. Kane, P. W. Likins and D. A. Levinson, Spucecraft Dynamics, pp. 422-43 1. McGraw-Hill (1983). 33. L. Brand, Vector und Tensor Analysis, pp. 329, 354.
John Wiley, New York (1947). 34. J. W. Kamman and R. L. Huston, User’s manual for
35. 36. 37. 38. 39.
UCIN DYNOCOMBS. National Technical Information Service Report No. PB85-240075/AS (1985). J. T. Wang and R. L. Huston, Kane’s equations with undetermined multipliers--application to constrained multibody systems. J. appl. Mech. 54, 424429 (1987). L. Meirovitch, Methodsif Analytical Dynamics, p. II McGraw-Hill, New York (1970). L. A. Pars, A Treatise on A&y&l Dynamics, p. 75. Ox Bow Press, Woodbridge, CT (1979). D. E. Whitney, Resolved motion rate control of manipulators and human prostheses. IEEE Trans. Man-Machine Syst. MM!+10, 47-53 (1969). D. E. Whitney, The mathematics of coordinated control of prosthetic arms and manipulators. J. Dynum. Sysr.
Measurement, Control 303-309 (1972). 40. M. Vukobratovic and M. Kircanski,’ A dynamic
approach to nominal trajectory synthesis for redundant manipulators. IEEE Trans. Systems. Man, Cybernetics SMC-14, 580-586 (1984). 41. J. M. Hollerback and K. C. Suh, Redundancy resolu-
tion of manipulators through torque optimization. Massachusetts Institute of Technology-Artificial Intelligence Laboratory, AI. Memo 882 (1986).
1988
0045-7949/88 $3.00 + 0.00 Pergamon Press plc
COMPUTATIONAL METHODS IN CONSTRAINED MULTIBODY DYNAMICS: MATRIX FORMALISMS J. T. WANG? and R. L. HUSTON~ tEngineering Mechanics Department, General Motors Research Laboratory, Warren, MI 48090, U.S.A. $,Department of Mechanical and Industrial Engineering, University of Cincinnati, Cincinnati, OH 45221, U.S.A. (Receiued 6 July 1987)
matrix-based analysis of multibody dynamics is presented. The analysis is applicable to large systems subjected to both geometric and kinematic constraints. The matrix formulation leads directly to algorithms for computational procedures. It also presents the analysis in a format where the dynamical
Abstract-A
behavior of the system is evident. An example of a towing chain is presented as an illustration.
1. INTRODUCTION
Recently there has been a marked increase in dynamical analyses of multibody systems. This interest stems from the ability to model a broad class of physical systems by multibody systems. The use of multibody system methods in dynamical analyses is analogous to the use of finite element methods in structural analyses. A multibody system may be defined as a collection of connected bodies such that there is at least one degree of freedom at the connecting joint. Figure 1 depicts a typical multibody system. The system may be internally constrained with closed loops and motion restricting joints, leading to geometric (or holonomic) constraint equations. The system also may be externally constrained by specifying the motion of some of the bodies, leading to kinematic (or non-holonomic) constraints. Multibody systems have been used to model manipulators, robots, cables, beams, and biosystems. References [l-lo] document a few of the applications of multibody systems methods. The successful application of these methods is based upon two factors: (1) advances in computational procedures and (2) advances in system modelling and analysis. References [1l-301 document some of these advances. In this paper we present an analysis which we believe provides a basis for an efficient computational procedure. The analysis employs a matrix formalism which readily leads to algorithms for computation. The paper is divided into five parts with the following part discussing some concepts useful in the sequel. The next two parts present the dynamic analysis and the governing equations. The fourth part presents an illustrative example and the final part contains a discussion and concluding remarks.
center of Bk. Let R be an inertial reference frame in which S moves. Let S have N bodies and n degrees of freedom characterized n generalized coordinates q, (r = 1, . . . , n). Then, Huston et al. [15, 161 have shown that the velocity of Gk and the angular velocity of Bk in R may be expressed in the forms: “k =
vmlk !h %a
and
Ok
=
%dk
Ch’btn
9
(1)
where h,,, (m = 1,2,3) are unit vectors fixed in R. The coefficients v,,,,~and w,,,[k(m = 1,2,3; 1 = 1, . . . , n; N) are components of the partial velocity k=l,..., and partial angular velocity vectors as defined by Kane et al. [15, 16,20,3 11. They may be viewed as ‘base vectors’ in the n-dimensional space of the 4,, in which S moves. (Regarding notation, in eqn (1) there is a sum on repeated indices over the range of the index.) In matrix notation, let the 3 x n array ok be defined as:
vk =
01Ik
v12k
U2lk
v22k
L?lk
%k . . . V2”k
hk
%,k
=
bmnkl.
(2)
I
Then v, may be expressed in matrix form as: vkl “k =
ok2
=
(3)
okq,
.
[Ivk3
where 4 is the array whose transpose dT is 4r= [4,,42,.
. .,61.
(4)
Next, let v be the set of N velocity arrays each having the form of eqn (3). Then v may be written as
2. PRELIMINARY CONSIDERATIONS
v = V&
(5)
2.1 Partial velocities and partial angular velocities
Consider again the multibody system S of Fig. 1. Let Bk be a typical body of S and let Gk be the mass
where V is the N x 3 x n block array of partial velocities.
331
332
J. T. WANGand R. L.
HUSTON
Fig. 1. A typical multibody system.
Similarly, the angular velocities may be expressed as:
and SJK SKL = SJL
o,=R,g
or
o=Qq,
where I is the 3 x 3 identity matrix and L refers to a third body B,. By repeated use of the last of eqn (1 1), we can obtain a transfo~ation matrix between Bk and the inertia frame R (designated by 0), i.e.
where the partial angular velocity arrays !& and Q are defined in the same manner as V, and V are defined. 2.2 Transformation matrices The relative orientations of the bodies may be defined in terms of orientation angles (such as Euler angles). The generalized coordinates q, may be identified with these angles. Transformation matrices relating unit vectors fixed in adjoining bodies then may also be defined in terms of these angles. Specifically, let Bj and Bk be typical adjoining bodies of the system and let n,i and nki (i = 1,2,3) be mutually perpendicular unit vector sets fixed in Bj and B,. Then let the 3 x 3 transformation matrix SJK be de&red as: SJK, = njp-nkq.
SOK=SOl
SOK, = qp.nkq.
-
[
sk2
ck3
sklsk2
sklck2ck3- sk3cki
sk3 skl
-cklckZsk3
-
Ck3skl
-sklck2sk3
-
ck3ckl
nk = SJKTnj.
(10)
It is readily seen that the transformation matrices obey the following identity and transitive relations: SJK SKJ=SJK
(8)
I I
SO&%,
= e~~~~~~SOK~~, and
SJKr=SJK
SJK-‘=I
(14)
sk2sk3
ckl ck2 ck3 -
where ski and g, represent the sine and cosine of @ki. SJK may be used to relate the unit vector sets by Hence, the derivative of SOK, is the expressions: dSOK,,/dt = l$,‘Wk X nkq njP= SJK,n, and nk#= SJK,ul,, (9) = ~p’@h,,kcir%m X
nj = SJKn,
(13)
dn,,/df = ok x nkg (no sum on k).
cklsk2
or in matrix form as:
(12)
Equation (13) is useful for developing algorithms for computing the time derivative of SOK: since the n,, are fixed in 4 their derivatives in inertia frame R are:
(7)
ck2
S12,...,SJK,
where the indicted product involves the bodies of the branches of S containing Bk, From eqn (7) SOK is seen to have the elements
For example, if ‘1,3, 1’ Euler angles 6,, (i = 1,2,3) (see [32]) are used to orient B, relative to B,, then SJK is
SJK =
(11)
(6)
(IS)
where ep,,,,is the permutation symbol [33]. In matrix form eqn (15) may be written as SdK = WOK SOK,
(16)
where WOK is the matrix whose ‘dual vector’ 1333is
Computational methods in constrained multibody dynamics ok,
i.e.
0 =
-
wk3
0
wk3
[-
OkI
Ok2
“k2 -wkl
0
7
(17)
I
where oki (i = 1,2,3) are the II,,~components of Ok. Equations (14) and (16) are especially significant in the numerical analyses of multibody system dynamics. They provide multiplication algorithms for computing derivatives-operations ideally suited for digital computation.
333
popular option is the use of Euler parameters for this purpose: as before, let Bj and Bk be typical adjoining bodies of the system. Then at any instant the change in orientation of Bk relative to Bj may be represented as a rotation about a line fixed in each body. If 12, is a unit vector parallel to this line and if ok is the rotation angle, four Euler parameters cki (i=l,... ,4) may be defined as: 6ki= Lkisin(8,/2),
i = 1,2,3
and ck4 =
wek/2h
(20)
2.3 Generalized speeds and Euler parameters
In eqns (1) (5) and (6) the generalized coordinate derivations q, are the coefficients of the partial velocity and partial angular velocity vectors. If the partial velocity and partial angular velocity vectors are viewed as base vectors in the n-dimensional space in which S moves, then the B (r = 1, . . . , n) may be viewed as magnitudes or measures of the movement. On many occasions it is convenient to use independent linear combinations of the 4, as measures of the motion of S. Such linear combinations are called ‘generalized speeds’. If U, (r = 1,. . . , n) is a set of generalized speeds, the transformation between 4, and the U, may be expressed in the matrix form 4 = Au
or
u = A T4,
and
fi=QA.
2 -1 .
~:,+&+&+t,-
(21)
An advantage of using Euler parameters is that there is a linear relation between the Euler parameters and the relative angular velocity components. This avoids singularities which arise when orientation angles are used [14]. Specifically, if e, and djk are the arrays of the Euler parameters and the relative angular velocity components, defined as
(18)
where A is an orthogonal matrix. By substituting from eqn (18) we can define partial velocity and partial angular velocity arrays B and fi associated with the generalized speeds. Specifically, P=VA
where the Akiare the components of & relative to unit vectors nji fixed in Bj. The Euler parameters are not unique. They are related by the equation
(19)
A convenient set of generalized speeds for the system of Fig. 1 is the set of local components of the relative angular velocity vectors [15]. There is, however, a nonlinear relationship between the orientation angles and the relative angular velocity components. This means that the orientation angles cannot simply be differentiated to produce the relative angular velocity components. Indeed, in general there are no functions which can be differentiated to produce the relative angular velocity components. Hence, if we
use relative angular velocity components as generalized speeds, specific generalized coordinates do not exist. To overcome this problem it is convenient to use variables other than orientation angles to define the relative orientations of the bodies. An increasingly
and ~i)kT=[~kl,~k2r~kj,o],
G-9
then 4 and Ok are related by the expressions 4 = Eok
and
where the orthogonal
Ek3
-
ck4
-;;; -6kl
ck2
ekl
ckl
ck2
6k4
ck3
Ek3
&k4
. -Ekl
I
(23)
matrix E is Ek4
E = (l/2)
uik = ETik,
-
ck2
-
(24)
I
Finally, when Euler parameters are used the transformation matrix SJK has the form
2.4 Derivatives of the partial angular velocity arrays
velocity
and partial
An advantage of using relative angular velocity components as generalized speeds is that the resulting partial angular velocity components may be directly identified with the elements of the transformation
J. T. WANG and R. L. HUSTON
334
matrices. Hence, since we have an algorithm for the differentiation of the transformation matrices [eqn (1611 we can then readily differentiate the partial angular velocity array. To see this, consider again the angular velocity wk of a typical body Bk of the system. From the addition theorem for angular velocities we can express o, as o,=~,+fi*+“‘+&,+&
(26)
where, as before, the carats over the symbols designate the relative angular velocity of the subscript body with respect to the adjacent lower numbered body. The sum in eqn (26) is taken over the bodies in the branch of S containing Bk. Consider Ok: it may be expressed as:
In matrix notation ak and a, may be expressed as [see eqns (3) and (6)] ak = hktf + ti,ri
and
ak = &ikuf vkti,
(30)
where we have dropped the carats for notational simplicity. The four arrays @,,&,ti,,,[k,v,,,,~and d, , or alternatively, R, h, V and v, play a central role in the analysis of multibody dynamics. 3. DYNAMICS
To explore the systems acting on be separated into forces; constraint
dynamics of S consider the force the bodies. Let these force systems three categories: externally applied forces; and inertia forces.
3.1 Applied forces where as before the n,, are unit vectors fixed in adjoining body B,. Let the e& be identified with the generalized speeds up where p is 3(k - 1) + i. Also, from eqn (9) let the nji be expressed as SO~~i~~. Then by substituting into eqns (27) and (26), and by comparing with eqns (1) and (18), we see that the wmIk may be identified with the elements of SOJ and the transformation matrices of the other bodies in the branch of Bk. Hence, the derivatives of the o,,,,~ are obtained directly from eqn (16). Additional details are given in [15]. Finally, it is seen that the elements of the partial velocity array can also be expressed in terms of the elements of the transformation matrices. To see this, consider again a typical body Bk with mass center Gk. Let pk locate Gk relative to a fixed point 0 in inertia frame R. Then pR may be expressed in terms of vectors fixed in the bodies in the branch of Bk. Let rj be one of these vectors. Let 4, be fixed in B,. Then in computing the velocity of Bk it is necessary to differentiate cj in R. This derivative fj is simply +r+x
pj.
(28)
Let the applied (or ‘active’) forces on each body Bk force system consisting of a Single force Fk passing through the mass center Gk together with a couple having torque Mk. Then the contribution by Fk and Mk to the generalized active force F,, associated with the generalized speed u,, is [15,31] be replaced by an equivalent
= Fkn%k +
Mkm %nik 3
(31)
where Fkm and Mkm are the n, components of Fk and mk. Observe that in eqn (31) the terms represent the contribution to F, from Fk and M,. Hence, if there is a sum on k we have the ~nt~bution from all the bodies leading to the total F,. Observe further that if the partial velocity and partial angular velocity vectors are viewed as base vectors defining the directions in which S can move, then F, represents for u, the projections of F, and Mk along these directions. Equation (31) may be written in the matrix form:
(32) Then by eqns (1) tj and similar derivatives may be expressed in terms of the transfo~ation matrix where (Pkand & are the arrays of rbm components of elements. Hence, the velocity of G, and the elements Fk and Mk. of the partial velocity arrays may also be expressed in terms of the transfo~ation matrix elements. Then 3.2 Constraint forces the derivatives of the partial velocity array elements Constraints and thus constraint forces arise when may be computed using the algorithm of eqn (16). the movement of S is restricted. This may occur in Additional details are given in 1151. By knowing the derivatives of the o,,,~~and =?&,k several ways: (1) there are fewer than six degrees of arrays we can compute the angular acceleration ak freedom at a joint; (2) S may have one or more closed and mass center acceleration ak of the bodies in R, i.e. loops; or (3) the movement of some of the bodies of S relative to each other or relative to inertia space R may be specified. As with the applied forces, let the constraint forces on each body be replaced by an equivalent force and system consisting of a single force FL passing through (29) mass center Gk (or alternatively through joint center
Computational methods in constrained multibody dynamics 0,) together with a couple having torque ML (or Ma). Then the contribution by FL and Mi to the generalized constraint force F; associated with the generalized speed u, is of the form [see eqn (32)]: F; = V:d;
+ i-l,‘&.
(33)
where 4; and II; are the arrays of no,,,components of Fk and Mk. The restriction of the movement of S may be described by equations of the form Bijuj=gi,
i=l,...,
m; j=l,...,
n; m
or more generally the forms F+F’+F*
=0
or
F+B’l+F*=O,
(40)
where F, F’ and F* are the arrays of generalized applied, constraint and inertia forces. Finally, observe that by substituting from eqns (29) and (37), Fy may be expressed as
or as
(34)
or in matrix form as
335
F* = -Mti
+J
(41)
where the mb and f; are Bu =g
(35)
where B is the M x n ‘constraint assay’ with elements b,.
Each of the individual equations of eqns (34) may be associated with a single constraint force or moment component. Indeed, it has been shown that if 1 is the combined array of constraint force and moment components then Fj and B are related by the expression [35] B’A = F;.
(36)
(42) and fi = -
crnk %lk
+
The inertia forces may be treated in the same manner as the applied ‘and constraint forces. Specifically, let the inertia force system on each body Bk be replaced by a single force F! passing through Gk together with a couple having torque Mz . Then F: and Mt may be expressed as
Ih
%,lk
ersmlksn%~k
&“,?,k up “rqkonpk
uqup)T
(43)
where the Ik,,,,,are the no,,,and nr,” components of I, and e,, is the permutation symbol [33]. The array IV, called the ‘generalized mass array’, is seen to be symmetric and positive definite. 4. GOVERNING
3.3 Inertia forces
i?npk up +
EQUATIONS
AND SOLUTION
METHODS
From eqns (34), (40) and (41) we seen that the governing equations for constrained multibody systems can be expressed in the matrix forms: Mti=F+f+BTL
(44)
Bu =g.
(45)
and
and Mf = -Ik.ak -ok
x (Ikeok),
(37)
where mk and Ik are the mass and central inertia dyadic [31] of Bk. Then in view of eqns (32) and (33) the generalized inertia force F: associated with generalized speed a, may be expressed as F:=
V:4:+Q,Tpc,+,
(38)
where as before 4: and p: are the arrays of b,,, components of Ff and Mt.
Since M is an n x n array and since B is an m x n array, these equations are equivalent to a set of m + n scalar equations for the n generalized speeds u and the m constraint force and moment components 1. The equations thus form a coupled differential-algebraic system of equations. A convenient solution method is as follows. Suppose an n x (n - m) matrix C can be found such that BC is zero. Then C is an ‘orthogonal complement’ [12, 18, 191 of B. Hence, we can eliminate 1 be premultiplying eqn (44) by CT, i.e. C’Mti = C’F+ C’f:
3.4 Dynamical equations Using Kane’s method take the simple form
Then by differentiating the dynamical
F,+F;+F:=O,
(46)
eqn (45) we have
equations Bti=g-h. (3%
(47)
Equations (46) and (47) constitute a set of n equations
J. T. WANGand R. L.
336
-
4 hlsec’
(1.22 mlsec’)
i;
Fig. 2. Towing chain-initial
configuration.
for the n unknowns u. By solving for u we can then use eqn (44) to find 1. 5. EXAMPLE:
A TOWING
HUSTON
given an acceleration of 4.0 ft/sec2 (1.22 m/set*) to the right. This system is thus a constrained multibody system with specified motion of two points of the system. The response of the system may be obtained by formulating the dynamics in the form of eqns (44) and the constraints in the form of eqns (45). The solution of these equations using the procedures of eqns (46) and (47) produces the time interval response as in Fig. 3. In developing this solution we used the computer code DYNOCOMBS developed by Kamman and Huston 1341.
Equations (46) are equivalent to the standard form of Kane’s equations. They represent a projection of eqns (38) onto a direction orthogonal to the constraint force array. This can be depicted geometrically by the triangle shown in Fig. 4. In Fig. 4 K and K* are the generalized force arrays employed by Kane et al. ]27,31,32]. Specifically, K and K* are
K=CTF
and
K*=CTF*.
(48)
CHAIN
As a simple example illustrating the use of some of these ideas consider the towing chain depicted in Fig. 2. The chain consists of six identical pinconnected rods. They each have length 1.0 ft (0.3048 m) and mass 1.0 slug (14.6 kg). One end of the chain is attached to a fixed point (representing an object to be towed) and the other end is given an acceleration (simulating the towing vehicle). Initially, the two ends are placed together and the chain links hang in a vertical gravity field. The towing end is then
Fig. 4. Geometric interpretation of generalized force arrays.
,t=1.m
Fig. 3. Towing chain ~nfiguration
at successive time intervals.
Computational
methods in constrained multibody dynamics
putational assistance of Mr Timothy P. King and Mr Yung
governing equations are then
The reduced
337
Sheng Liu. K+K+=O.
(49)
Observe that the constraint forces do not contribute to these generalized forces or to the reduced governing equations. Hence, the constraint forces may be ignored in analyses using Kane’s equations. This is a principal advantage of Kane’s equations. Next, the constraint force and moment array 1 in eqns (36), (40) and (44) may be considered as Lagrange multipliers. Hence, eqns (40) may be considered as a generalization of Kane’s equations to include constraint forces-a procedure long used with Lagrange’s equations [36,371. Finally, eqns (44) and (45) are in a form which is convenient for multibody dynamics analyses. For example, we can formally eliminate 1. If we solve eqn (44) for zi we obtain ti = M_‘(F +f+
BTA).
(50)
Then by substituting into eqn (47) [the differentiated form of eqn (45)] we obtain BM-‘(F+f+B=A)=g
-Lb.
(51)
Solving for I we have 1 = (BM-‘P--‘[g
- lh - BM-‘(F
Finally, by substituting
zi= M-‘(F
+f)
+
+f)].
(52)
into eqn (50) we obtain
M-‘Br{(BM-‘BT)-’
[g - Lh - BW’(F
+f)]j.
The
term M-‘BT(BM-‘Br)-’ is the pseudoinverse’ [38,39] of B. That is, B[M-‘Br(BM-‘Br)-‘I=
Z,
(53)
‘weighted
(54)
where I is the m x m identity
matrix. The pseudoinverse has been employed by a number of investigators in studying inverse kinematic and inverse dynamic problems in robotics [38-41]. Finally, it is believed that the range of applications of the procedures exposited herein and the governing equations have yet to be developed. The advantage of the matrix formalism is its simplicity and its ready conversion to computer algorithms. This, coupled with the simplicity and eloquence of Kane’s equations, makes the governing equations very attractive for application to constrained, redundant multibody systems. Specific applications in robotics, biosystems and structural analysis are anticipated. Acknowledgements-The research for this paper was partially supported by the National Science Foundation under Grant MSM-8612970. This support is gratefully acknowledged. The authors also gratefully acknowledge the com-
REFERENCES
1. G. E. Fliescher, Multibody systems applications and simulations at the jet propulsion laboratory. In Dynamics of Multibody Systems (Edited by K. Magnus), pp. 3647. Springer, New York (1978). 2. J. G. Gallenstein and R. L. Huston, Analysis of swimming motions. Human Factors 15, 91-98 (1973). 3. J. C. Huston, M. W. Harlow and R. L. Huston, A comprehensive, three-dimensional head-neck model for impact and high acceleration studies. Aviation, Space Envir. Medicine 49, 205-210 (1978). 4. R. L. Huston, A summary of three-dimensional, grossmotion, crash-victim simulators. In Structural Mechanics Software Series, Vol. I, pp. 611622. University Press of Virginia, Charlotttesville, VA (1977). 5. R. L. Huston, R. E. Hessel and J. M. Winget, Dynamics of a crash victim-a finite segment model. AIAA Jnll4, 173-178 (1976). 6. R. L. Huston and C. E. Passerello, On the dynamics of a human body model. J. Biomech. 4, 369378 (1971). 7. J. W. Kamman and R. L. Huston. Modellina of submerged cable dynamics. Cotnput. Srruct. 20;623-629 (1985). 8. T. R. Kane and D. A. Levinson, The use of Kane’s dynamical equations in robotics. Int. J. Robotics Res. 2, 3-21 (1983). 9. F. A. Kelly and R. L. Huston, The development of equations of motion of single arm robots. IEEE Trans. Systems, Man, Cybernetics, SMC-12, 259-266 (1982). 10. J. M. Winget and R. L. Huston, Cable dynamics-a finite segment approach. Comput. Struct. 6, 475480 (1976). 11. Y. 0. Bayazitoglu and M. A. Chace, Methods of automated dynamic analysis of discrete mechanical systems. J. appl. Mech. 40, 809-819 (1973). 12. H. Hemami and F. C. Weimer, Modelling of nonholonomic dynamic systems with applications. J. appl. Mech. 48, 177-182 (1981). 13. W. W. Hooker, A set of r dynamical attitude equations for an arbitrary n-body satellite having r rotational degrees of freedom. AIAA Jnl8, 1205-1207 (1970). 14. R. L. Huston and C. E. Passerello, Eliminating singularities in governing equations of mechanical systems. Mech. Res. Commun. 3, 361-365 (1976). 15. R. L. Huston, C. E. Passerello and M. W. Harlow, Dynamics of multirigid-body systems. J. appl. Mech. 45, 889-894 (1978). 16. R. L. Huston and C. E. Passerello, On multi-rigid-body system dynamics. Comput. Struct. 19, 439446 (1979). 17. R. L. Huston and C. E. Passerello, Multibody structural dynamics including translation between the bodies. Comput. Srrucr. 11, 715-720 (1980). 18. J. W. Kamman and R. L. Huston, Dynamics of constrained multibody systems. J. uppl. Mech. 51, 899-904 (1984). 19. J. W. Kamman and R. L. Huston, Constrained multibody system dynamics. Comput. Srruct. 18, 999-1003 (1984). 20. T. R. Kane and D. A. Levinson, Formulation of equations of motion for complex spacecraft. J. Guidance Control 3, 99112 (1980). 21. P. W. Likins, Guasi-coordinate eauations for flexible spacecraft. AIAA Jnl 13, 524-526 (1975). 22. P. W. Likins, Dynamic analysis of a system of hingeconnected rigid bodies with nonrigid appendages. In?. J. Solids Struct. 9. 1473-1487 (1973). 23. P. W. Likins, Point-connected rigid bodies in a topological tree. Celestial Mech. 11, 301-317 (1975).
J. T. WANG and R. L. HUSTON
338
24. L. K. Lilov and V. A. Chirikov, On the dynamics equations of systems of interconnected bodies. PMM (J. anal. Mum. Mech.) 45. 383-390 (1982). 25. L. L& and M. Lorer; Dynamic analysis of multirigidbody system based on the Gaus principle. ZAMM 62, 539-545 (1982). 26. P. E. Nikravesh and I. S. Chung, Application of Euler parameters to the dynamic analysis of threedimensional constrained mechanical systems. J. Mech. Des. 104, 785-791 (1982). 21. C. E. Passerello and R. L. Huston, Another look at nonholonomic systems. J. uppl. Mech. 40, 101-104 (1973). 28. R. Schwertassek and R. E. Roberson, A state-space dynamical representation for multibody mechanical systems, Part I: systems with tree configurations. Actu Mechanica 50, 141-161 (1984). 29. R. Schwertassek and R. E. Roberson, A state-space dynamical representation for multibody mechanical systems, Part II: systems with closed loops. Actu Mechun icu 51, 15-29 (1984).
30. M. Vukobratovic, Computer method for dynamic model construction of active articulated mechanisms using kinetostatic approach. Me& Muchine Theory 13, 19-39 (1978).
31. T. R. Kane and D. A. Levinson, Dynamics: Theory and Applications, pp. 15, 45. McGraw-Hill (1985).
32. T. R. Kane, P. W. Likins and D. A. Levinson, Spucecraft Dynamics, pp. 422-43 1. McGraw-Hill (1983). 33. L. Brand, Vector und Tensor Analysis, pp. 329, 354.
John Wiley, New York (1947). 34. J. W. Kamman and R. L. Huston, User’s manual for
35. 36. 37. 38. 39.
UCIN DYNOCOMBS. National Technical Information Service Report No. PB85-240075/AS (1985). J. T. Wang and R. L. Huston, Kane’s equations with undetermined multipliers--application to constrained multibody systems. J. appl. Mech. 54, 424429 (1987). L. Meirovitch, Methodsif Analytical Dynamics, p. II McGraw-Hill, New York (1970). L. A. Pars, A Treatise on A&y&l Dynamics, p. 75. Ox Bow Press, Woodbridge, CT (1979). D. E. Whitney, Resolved motion rate control of manipulators and human prostheses. IEEE Trans. Man-Machine Syst. MM!+10, 47-53 (1969). D. E. Whitney, The mathematics of coordinated control of prosthetic arms and manipulators. J. Dynum. Sysr.
Measurement, Control 303-309 (1972). 40. M. Vukobratovic and M. Kircanski,’ A dynamic
approach to nominal trajectory synthesis for redundant manipulators. IEEE Trans. Systems. Man, Cybernetics SMC-14, 580-586 (1984). 41. J. M. Hollerback and K. C. Suh, Redundancy resolu-
tion of manipulators through torque optimization. Massachusetts Institute of Technology-Artificial Intelligence Laboratory, AI. Memo 882 (1986).