University of Southern California Dissertations and Theses

On the synthesis of dynamics and control for complex multibody systems

(USC Thesis Other)

On the synthesis of dynamics and control for complex multibody systems

PDF

Download a page range

Download transcript

Copy asset link

Request this asset

Transcript (if available)

Content ON THE SYNTHESIS OF DYNAMICS AND CONTROL FOR COMPLEX MULTIBODY SYSTEMS by Aaron Dane Schutte A Dissertation Presented to the FACULTY OF THE GRADUATE SCHOOL UNIVERSITY OF SOUTHERN CALIFORNIA In Partial Fulllment of the Requirements of the Degree DOCTOR OF PHILOSOPHY (AEROSPACE ENGINEERING) May 2009 Copyright 2009 Aaron Dane Schutte Dedication To my wife Jen and my son Cole ii Acknowledgements I would like to convey my gratitude for my advisor Professor Firdaus Udwadia. This dissertation was possible through his guidance and perpetual encouragement. He has a unique ability to understand the eld of mechanics in a simple way, and it has motivated me to pursue my goals in the theoretical areas of my interest. I have always highly valued his advice and I can never thank him enough. I would also like to give a special thanks to Dr. Jay Bernard, Tom McLaughlin, and Timothy Paulitz of The Aerospace Corporation for supporting my endeavors. They pro- vided me with a source of intellectual capacity and aided in my professional development. I could not have carried out my studies without their support. I am also greatly indebted to Dr. Jordi Puig-Suari for introducing me to the eld of spacecraft dynamics and control. Finally, most importantly, I want to thank my family. My wife Jennifer has been an inspiration for me from the beginning. Her love and support has guided me along through the most important parts of my life, and she has been there for me every step of the way. My parents, Dane and Lynda, and my sister, Amber, have always been a source of unconditional love and support, and they always did everything they could to help me achieve my goals. I also owe a great deal of gratitude to the Tomlinson family for their continual and enthusiastic support. iii Table of Contents Dedication ii Acknowledgements iii List of Tables vi List of Figures ix Abstract x Chapter 1: Introduction 1 1.1 Complex Multibody Dynamics and Control . . . . . . . . . . . . . . . . . 3 1.2 Outline . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6 Chapter 2: Fundamental Equation for Multibody Systems 8 2.1 Kinetic Energy of N Non-Interacting Systems of Particles . . . . . . . . . 8 2.2 Kinetic Energy of N Non-Interacting Rigid Bodies . . . . . . . . . . . . . 11 2.3 Multibody Coordinates . . . . . . . . . . . . . . . . . . . . . . . . . . . . 14 2.4 Equations of Motion for N Rigid Bodies . . . . . . . . . . . . . . . . . . . 17 2.5 Modeling a Complex N Rigid Body System . . . . . . . . . . . . . . . . . 20 2.6 Explicit Equations of Motion for Complex N Rigid Body Systems . . . . 21 2.7 Example: Two Spacecraft Linked at an Arbitrary Point . . . . . . . . . . 24 2.8 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 28 Chapter 3: Explicit Control for Multibody Systems 30 3.1 Multibody Control Formulation . . . . . . . . . . . . . . . . . . . . . . . . 30 3.2 Trajectory Control . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 37 3.3 Stabilization . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 40 3.4 Underactuated Control . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 42 3.5 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 52 Chapter 4: Rotational Dynamics and Control Using Quaternions 54 4.1 Rotational Kinetic Energy in Terms of Quaternions . . . . . . . . . . . . . 55 4.2 Lagrange's Rotational Equations of Motion . . . . . . . . . . . . . . . . . 57 iv 4.3 Control of Rotational Motion in Terms of Quaternions . . . . . . . . . . . 62 4.4 Rotational Controller 1 . . . . . . . . . . . . . . . . . . . . . . . . . . . . 65 4.5 Rotational Controller 2 . . . . . . . . . . . . . . . . . . . . . . . . . . . . 69 4.6 Numerical Examples . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 75 4.6.1 Rotational Controller 1 . . . . . . . . . . . . . . . . . . . . . . . . 75 4.6.2 Rotational Controller 2 . . . . . . . . . . . . . . . . . . . . . . . . 79 4.7 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 85 Chapter 5: Spacecraft Precision Pointing 89 5.1 Precision Pointing Description . . . . . . . . . . . . . . . . . . . . . . . . . 89 5.2 Model Development . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 91 5.3 Multibody Spacecraft Pointing Control . . . . . . . . . . . . . . . . . . . . 97 5.4 Numerical Examples . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 99 5.5 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 106 Chapter 6: Multibody Spacecraft Tumbling Control 107 6.1 Multibody Spacecraft Tumbling Description . . . . . . . . . . . . . . . . . 108 6.2 Model Development . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 111 6.3 Tumbling Control . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 114 6.4 Numerical Examples . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 116 6.5 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 123 Chapter 7: Conclusions 125 Bibliography 130 Appendix A 133 Appendix B 136 Appendix C 138 v List of Tables 1 Eigenvalues of J F for Rotational Controller 2 with identical i . . . . . . . 137 2 Eigenvalues of J F for Rotational Controller 2 with distinct i . . . . . . . . 137 3 Mars physical properties. . . . . . . . . . . . . . . . . . . . . . . . . . . . 138 4 Zonal harmonic coecients up to degree four (k = 4) for the planet Mars. 138 5 Tesseral and sectoral harmonic coecients up to degree four (k = 4) for the planet Mars. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 139 vi List of Figures 2.1 Schematic of N independent and non-interacting systems of particles. The i th system of particles contains a collection of ` i particles dened relative to an inertial coordinate frame with axes R i 1 , R i 2 , and R i 3 . . . . . 9 2.2 Schematic of N independent and non-interacting rigid bodies. The i th rigid body V i rotates with an absolute angular velocity ! i and contains a body-xed coordinate frame with axes i , i , and i xed at its center of mass. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11 2.3 A two spacecraft system|denoted by rigid bodiesV 1 andV 2 |elastically linked together at the points P 1 and P 2 by a spring whose constants k l andk nl denote the linear and cubically nonlinear restoring forces exerted by the spring. The two spacecraft are free to rotate independently about the line P 1 P 2 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 25 3.1 A prescribed trajectory (q c ; _ q c ) and the actual trajectory (q; _ q) over a time duration t2 (t 0 ;t f ) along the modeling constraint manifold dened by ' m \ m \ m . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 39 3.2 An actual trajectory (q; _ q) converging along the modeling constraint man- ifold' m \ m \ m to the desired trajectory (q c ; _ q c ), which travels along the control constraint manifold dened by c \ c . . . . . . . . . . . . . 41 3.3 A surface vessel vehicle with mass m and inertia J. . . . . . . . . . . . 43 3.4 Surface vessel vehicle trajectories (dotted lines) near the desired ellipse path with eccentricity e = 0:745. . . . . . . . . . . . . . . . . . . . . . . 50 3.5 Control action required to approach and move along the desired ellipse path. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 51 3.6 Error in control constraints 1;c and 2;c while moving along the ellipse. 51 4.1 Rotational Controller 1 { Rotational trajectories (dotted lines) on a unit 2-sphere. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 76 vii 4.2 Rotational Controller 1 { The rotation (components of the 4-vector,(t)) generated by the control acceleration that is explicitly given by equation (4.43) as a function of time. . . . . . . . . . . . . . . . . . . . . . . . . . 78 4.3 Rotational Controller 1 { Components of the control torque about the body-xed axes given by equation (4.41). . . . . . . . . . . . . . . . . . 78 4.4 Rotational Controller 1 {The modeling constraints' m and _ ' m as a func- tion of time. Note the vertical scale. . . . . . . . . . . . . . . . . . . . . 79 4.5 Rotational Controller 2 i = 1 8 ; i = 0; 1; 2; 3 { Rotational trajectories (dotted lines) on the unit 2-sphere. . . . . . . . . . . . . . . . . . . . . . 80 4.6 Rotational Controller 2 i = 1 8 ; i = 0; 1; 2; 3 { The rotation (compo- nents of the 4-vector, (t)) generated by the control acceleration that is explicitly given by equation (4.55) as a function of time. . . . . . . . . . 81 4.7 Rotational Controller 2 i = 1 8 ; i = 0; 1; 2; 3 { Components of the con- trol torque about the body-xed axes given by equation (4.41). . . . . . 82 4.8 Rotational Controller 2 i = 1 8 ; i = 0; 1; 2; 3 { The modeling constraints ' m and _ ' m as a function of time. Note the vertical scale. . . . . . . . . 82 4.9 Rotational Controller 2 0 = 1 8 ; 1 = 1 3 ; 2 = 2 7 ; 3 = 1 2 { Rotational tra- jectories (dotted lines) on the unit 2-sphere. . . . . . . . . . . . . . . . . 84 4.10 Rotational Controller 2 0 = 1 8 ; 1 = 1 3 ; 2 = 2 7 ; 3 = 1 2 { The rotation (components of the 4-vector, (t)) generated by the control acceleration that is explicitly given by equation (4.55) as a function of time. . . . . . 85 4.11 Rotational Controller 2 0 = 1 8 ; 1 = 1 3 ; 2 = 2 7 ; 3 = 1 2 { Components of the control torque about the body-xed axes given by equation (4.41). 86 4.12 Rotational Controller 2 0 = 1 8 ; 1 = 1 3 ; 2 = 2 7 ; 3 = 1 2 { The modeling constraints ' m and _ ' m as a function of time. Note the vertical scale. . . 86 5.1 Schematic of a two spacecraft system denoted by rigid bodies V 1 andV 2 in orbit around a rotating central body with a non-uniform gravity eld. A payload is xed to each spacecraft located by the position vectors r i P , i = 1; 2: The desired pointing axis ^ n i for the i th spacecraft is required to point along the relative distance vector . . . . . . . . . . . . . . . . . . 90 5.2 The body-xed coordinate axes (x b ;y b ;z b ) xed to the central body and the coordinate frame (^ x i ; ^ y i ; ^ z i ) relative to the inertial coordinate frame (x;y;z). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 93 viii 5.3 Orbital trajectories of spacecraft V 1 and V 2 . Spacecraft V 1 is in an elliptic equatorial orbit with an eccentricity e = 0:5 and spacecraft V 2 is in a polar circular orbit. . . . . . . . . . . . . . . . . . . . . . . . . . . . 101 5.4 Pointing error over the time period t2 (0; 150)s. . . . . . . . . . . . . . 102 5.5 Control torques, i c , generated by equation (5.30) required to approach the payload pointing requirement over the time interval t2 (500; 38600). 103 5.6 Control forces, i c , generated by equation (5.30) required to approach the payload pointing requirement. . . . . . . . . . . . . . . . . . . . . . . . . 104 5.7 Pointing error over the time period t2 (500; 38600)s. . . . . . . . . . . . 104 5.8 Error in T 1 and T _ over the time interval t2 (0; 38600)s. . . . . . 105 5.9 Control torques, i c , generated by equation (5.30) required to maintain the payload pointing requirement over the time interval t2 (500; 38600)s. 105 6.1 Schematic of the multibody spacecraft system in orbit relative to the rotating central body. . . . . . . . . . . . . . . . . . . . . . . . . . . . . 108 6.2 Schematic of the orbital and tumbling trajectory of the multibody space- craft system wheny cm = 0. The dark gray circle and the light gray circle denote spacecraft V 1 and V 2 , respectively. The inertial coordinate sys- tem is oriented so that the y-axis is directed into the orbital plane. The spacecraft system is orbiting in the counter-clockwise direction and the principal rotational states and _ are listed at intervals A, B, C, and D, which are located at 90 degree increments along the orbit. . . . . . . 110 6.3 Schematic of spacecraftV 1 andV 2 showing spring connections atP 1 and P 2 . The tumbling axis is depicted along the axes i , i = 1; 2: . . . . . . 117 6.4 Error in the 8 modeling constraints ' m and m over the time duration t2 (0; 13645)s. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 119 6.5 Error in the 12 control constraints c showing their satisfaction through- out the integration. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 120 6.6 Control forces, i c , required to converge to the desired tumbling trajectory.121 6.7 Control torques, i c , required to converge to the desired tumbling trajectory.121 6.8 Control forces, i c , required to maintain the desired tumbling trajectory over 2 orbits. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 122 ix 6.9 Control torques, i c , required to maintain the desired tumbling trajectory over 2 orbits. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 122 6.10 Vibrational motion between spacecraft V 1 and V 2 throughout the tum- bling trajectory. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 123 6.11 Angular velocity, ! i , of the i th spacecraft throughout the tumbling tra- jectory. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 124 x Abstract This dissertation develops in a unied manner a new and simple approach for the modeling and control of complex multibody systems. Complex multibody systems are those nonlinear mechanical systems consisting of individual subsystems that are describ- able by nonlinear dierential equations. The interaction of the various subsystems is, in general, governed by nonlinear 'elements.' The nonlinear analysis of general multibody systems presents theoretically challenging problems in both its dynamics and controls as- pects. The theoretical developments obtained in this dissertation permit the modeling of multibody systems so that no restrictions are imposed in its formulation, except that its physical model is continuous. The idea of permissible multibody control is developed, and an eective method is provided for generating nonlinear controllers that satisfy a general set of control objectives. The generalized acceleration describing the rotational motion of a rigid body in terms of quaternions is directly obtained thus providing a rotational description for the general motion of multibody systems. Utilizing this formulation, two new control strategies are obtained that explicitly yield the nonlinear control torque required to re-orient a rigid body from an arbitrary rest orientation to another arbitrary rest orientation. The ap- plication of this new approach to problems in complex multibody spacecraft systems is carried out for spacecraft precision pointing and to the precise tumbling control of an elas- tically connected multibody spacecraft system. Numerical examples are provided showing the ease of implementation of the methodology and the accuracy with which the control objectives are satised. These results illustrate the capability to easily generate physical models and provide exact control of highly nonlinear, complex multibody systems. xi Chapter 1 Introduction T HE central problem of determining the motion and the forces that generate the motion of discrete dynamical systems is a long studied and well developed subject. The principal importance of this problem has been evident in numerous elds of study not only for the sake of knowledge, but in its application to solve real life problems. In mod- ern times, the dynamic and control analysis of systems with many degrees of freedom is an important aspect in engineering and science applications. Systems with many degrees of freedom include multibody systems, which are systems containing multiple interacting bodies. Practical examples of multibody systems include aerospace vehicles and robotic systems. In fact, these applications alone have been highly in uential to the study of multibody systems since they are built to solve dicult problems such as exploring the surface of other planets. The dynamics of multibody systems are generally highly non- linear, introducing complex problems that are often only understood through numerical simulations. Indeed, the study of multibody dynamics is not a specialized eld, for in reality, the motion of a multibody system is an extension of the motion of a single particle. The fundamental principles laid out by Newton have dictated the evolution of the eld of dynamics from its inception. In its most basic form, the study of dynamics is concerned 1 with the motion of a single free particle with mass m under the action of given forces by assuming the existence of an inertial frame of reference. Relative to this inertial frame of reference the particle is displaced by the position 3-vector q = (x;y;z) and has a velocity _ q = ( _ x; _ y; _ z). The 'given' forces F (q; _ q;t) that act on the particle are known functions of the particles position and velocity, and also of time. An important problem in particle dynamics is the orbital motion of a satellite around a massive central body. The given force in this problem is the force of gravity. Although the exact origins of gravity are not fully understood, we nevertheless know that in its simplest form, its force acts at a distance and follows from the inverse-square law and, accordingly, it is a given force. In nature, such systems typically exhibit certain symmetry properties, which are identied as the invariants of motion. The invariants of motion for an orbiting satellite are its orbital angular momentum and its total energy, which is given by the sum of its kinetic and gravitational potential energy. In classical dynamics, these invariants are extractions of the laws of conservation and they appear in the form '(q; _ q) = 0. It follows that the orbital motion of the satellite evolves over time such that these known invariant functions of position,q, and velocity, _ q, remain invariant. Similarly, we may think of the invariants of motion as constraints on the motion wherein the presence of the constraints dictate the evolution of the particles position and velocity. We can then think of the particle as being acted upon by the given forces and controlled by constraints [17]. Taking a step beyond particle dynamics, we can extrapolate the idea of constraints that control the motion of a particle to constraints that control the motion of systems of particles. This is a direct result of contributions by D'Alembert and Lagrange [12]. A system of particles forms a basis for the study of multibody dynamics and the idea of constraints are intrinsic to its description. However, the point of view that constraints 2 control the motion of multibody systems is often excluded in favor of the view that constraints restrict the motion of multibody systems. In this dissertation, we shall exploit both of these viewpoints with the aim of developing a general and simple methodology for the modeling and control of multibody systems. Obtaining equations of motion that govern the motion of constrained discrete dynam- ical systems has been pursued by many scientists including Appell [3], Dirac [6], Gauss [8], Gibbs [9], and Lagrange [12]. Recently, a fundamental result obtained by Udwadia and Kalaba [27, 28] has opened up the traditional boundaries of mechanics. They ob- tain explicitly, general equations of motion for constrained discrete dynamical systems in terms of the generalized coordinates that describe the systems conguration. With the aid of this general formulation, we shall be concerned throughout this dissertation with describing the general motion of multibody systems and determining explicitly the general control action that is required to control the multibody system along a desired constraint, where the desired constraints will be seen as a set of control objectives. 1.1 Complex Multibody Dynamics and Control As already stated, the motion of a multibody system is generally highly nonlinear and also complex since their equations of motion can easily become dicult to evaluate. For- mulating equations of motion using existing methods is a nontrivial task, and there is no single technique that yields a direct and simple way of modeling general multibody systems, especially for systems with many degrees of freedom. In the area of spacecraft dynamics, this aspect was indicated in [11] and a number of formulations were provided to evaluate their advantages and disadvantages. In more recent monographs [2, 16], we nd 3 that most multibody systems are described using equations of motion beginning from rst principles. Upon adding or removing structural connections at the so-called nodes|often termed joint kinematics|that exist between various subsystems, a complete re-modeling of the multibody system is nearly always required. Moreover, many methods to date have dicultly gaining insight into the multibody problem at hand because they are not general in implementation, but rather assume that the system takes a special form with possibly many coordinate transformations. The control of multibody systems is most often reduced to some form of linear feedback control in order to achieve stabilization and tracking objectives. Furthermore, the linear controller utilized almost always relies on PID-type controllers as pointed out in [4]. There are many open questions in the control of general multibody systems since they are generally nonlinear, and a linear feedback controller will generally not be eective. In addition to nonlinearities, multibody systems can have multiple time scales since both translational and rotational coordinates are required to describe the systems conguration. This further complicates the control design because coupled translational and rotational motions may occur. A theme in this dissertation is to provide a general way of dealing with these types of problems in multibody systems, and provide an ecient and eective way of constructing controllers with no simplifying assumptions in its formulation. Another theme in this dissertation is the application of the general methodology devel- oped herein to problems in multibody spacecraft dynamics and control. These problems epitomize a complex multibody system since they operate at vastly dierent time scales (orbital, rotational, and elastic motions), and they are generally described by highly non- linear and non-autonomous equations of motion. There has been some eort recently in studying the coupling between rotational and orbital motions in spacecraft systems. 4 Quadrelli et al. [18] demonstrate the complex attitude and orbital dynamics interac- tions for a conceptual Jovian moon tour spacecraft. The interaction between attitude and orbital motions is a result of their requirement for a low-thrust propulsion system to follow complex orbital trajectories such as spiraling maneuvers. Using a linear con- troller, they provide simulation results with reasonably accurate assessments of pointing performance under certain criteria. Sanyal et al. [19, 20] consider the dynamics and control of an elastic dumbbell-shaped spacecraft in a central gravitational eld. Using a Lagrangian approach, they formulate equations of motion for the dumbbell spacecraft in two dimensions. Through linearization they present controllability results of the orbital and rotational motions of the dumbbell spacecraft showing that it is possible to obtain eective control by using coupling between orbital, rotational, and elastic degrees of free- dom. Although it wasn't demonstrated, they point the way towards the possibility of obtaining certain orbital maneuvers using only rotational and elastic control inputs. The formation ight of spacecraft is another growing area of research that requires advances in the dynamics and control of multibody systems. Early spacecraft formation ying researchers were focused on relative motion, interception, and docking control, and concern was mostly focused on the orbital motion and not the rotational motion. Many attempts in developing a control theory for these problems were made by linearizing the dynamics about a reference orbit. In the case of the Clohessy-Wiltshire equations [5], an ideal non-perturbed circular orbit is assumed as a reference orbit, and the relative motion of bodies (particles) was analyzed about this reference orbit. In more recent studies, Schaub et al. [22] present two nonlinear feedback control laws for a spacecraft formation problem that maintains aJ 2 invariant relative orbit of a "deputy" satellite with respect to a "chief" satellite using mean orbital elements and Cartesian coordinates. This particular 5 approach utilizes Lyapunov stability theory by dening a Lyapunov function in terms of the tracking error, which is parameterized by mean orbital elements and Cartesian coordinates. The theoretical developments in this dissertation are motivated by complex multibody problems of this nature. The general approach developed herein fundamentally diers from current methods in multibody dynamics in that both the modeling and control aspects of complex multibody systems are carried out in a direct and unied manner with no simplifying assumptions or approximations. 1.2 Outline This dissertation is organized as follows. Chapter 1 provides an introduction to complex multibody systems and discusses related work. In chapter 2, the fundamental equation for multibody systems is developed and we dene the modeling constraints required to de- scribe general multibody systems. The permissible control action of a multibody system is derived in chapter 3, and the development of control constraints, which are classied as the control objectives, is carried out. In chapter 4, the generalized accelerations govern- ing the rotational motion of a rigid body in terms of unit quaternions is derived and two new methods for explicitly obtaining the nonlinear control torque required to re-orient a rigid body are obtained. The modeling and control of a multibody spacecraft precision pointing problem is carried out in chapter 5, and numerical examples are provided to show the capability of the development. In chapter 6, the precise tumbling control is de- rived for a complex multibody system consisting of two rigid bodies elastically connected 6 at arbitrary locations relative to their body-xed coordinate frames. Numerical exam- ples demonstrating the ecacy of the approach are again provided. Finally, concluding remarks are provided in chapter 7. 7 Chapter 2 Fundamental Equation for Multibody Systems In this chapter, a general formulation is developed to describe the dynamics of multiple rigid body systems. In multibody dynamics, the complexity of the equations of motion grows rapidly as the number of bodies in the system increases. In most methods available to date, the equations of motion are generated using recursive methods that treat the multibody system as a tree-like structure wherein each body is connected at hinge or joint locations xed at arbitrary points on each body. Furthermore, the coordinates in these equations of motion are taken relative to the so-called joint coordinates so that a minimum set of equations of motion are obtained [23]. Modeling general multibody systems using these methods is dicult if not impossible since applying general constraint equations or forcing functions is not straightforward. In the following, a new and simple approach of systematically modeling a complex system of N rigid bodies is developed. 2.1 Kinetic Energy ofN Non-Interacting Systems of Parti- cles Consider theN independent and non-interacting systems of particles in Figure (2.1). The 8 1 R N R 1 1 m N m 1 1 1 r N r 1 1 m i i m l i j m i m 1 i r 1 i j r i i r l 1 j m 1 i m l 1 i r l 1 j r i j p i R N i m l N j m N j r N i r l i m N m i R 3 i R 2 i R 1 Figure 2.1: Schematic of N independent and non-interacting systems of particles. The i th system of particles contains a collection of ` i particles dened relative to an inertial coordinate frame with axes R i 1 , R i 2 , and R i 3 . i th system of particles contains a collection of` i particles dened relative to a xed inertial coordinate frame with axes R i 1 , R i 2 , and R i 3 . For each of the N systems of particles, let the position vector R i = R i 1 R i 2 R i 3 T denote the inertial position of the center of mass for the i th system. The position vector r i j = [ i j i j i j ] T then denotes the position of the j th particle with mass m i j measured from the center of mass of the i th system. The individual particles are located relative to the inertial coordinate frame by the position vector p i j =R i +r i j : (2.1) 9 The total kinetic energy of the i th system is then found by the sum of the individual kinetic energies of each of the ` i particles as T i = 1 2 ` i X j=1 m i j ( _ p i j _ p i j ) (2.2) = 1 2 ` i X j=1 m i j ( _ R i + _ r i j ) ( _ R i + _ r i j ) (2.3) = 1 2 _ R i _ R i ` i X j=1 m i j + _ R i ` i X j=1 m i j _ r i j + 1 2 ` i X j=1 m i j ( _ r i j _ r i j ): (2.4) The rst term in equation (2.4) is simplied by noting that the summation of m i j over the ` i particles is just the total mass of the i th system given by m i = ` i X j=1 m i j ; (2.5) and the second term is zero sincer i j is measured from the center of mass of thei th system, which yields ` i X j=1 m i j r i j = 0: (2.6) The total kinetic energy of the i th system of particles is thus T i = 1 2 m i ( _ R i ) 2 + 1 2 ` i X j=1 m i j ( _ r i j ) 2 ; (2.7) where ( _ R i ) 2 = _ R i _ R i and ( _ r i j ) 2 = _ r i j _ r i j . 10 2.2 Kinetic Energy ofN Non-Interacting Rigid Bodies The total kinetic energy in equation (2.7) was derived for N independent and non- interacting systems of particles. Let us now assume the N systems of particles make-up N independent and non-interacting rigid bodies denoted by V 1 ;V 2 ;:::;V N . The term rigid body indicates a collection of a large number of particles where the distance between any two particles is constant. In this N rigid body system, we shall assume the i th rigid body has a total massm i and also a coordinate frame attached to it, which is xed at its center of mass given by the axes i , i , and i as shown in Figure (2.2). Furthermore, 1 R N R 1 m N m 1 η 1 μ 1 γ N γ N μ N η 1 V N V i η i μ i γ i j r i V i R i m i R 1 i R 2 i R 3 i ω 1 ω N ω Figure 2.2: Schematic of N independent and non-interacting rigid bodies. The i th rigid bodyV i rotates with an absolute angular velocity! i and contains a body-xed coordinate frame with axes i , i , and i xed at its center of mass. each rigid body rotates with an absolute angular velocity ! i whose components in the 11 i th body-xed coordinate frame are given by ! i = [! i ! i ! i ] T . Consequently, the j th particle of the i th rigid body has a velocity _ r i j =! i r i j (2.8) relative to the origin of the body-xed coordinate frame. Since we have assumed the i th body is rigid, we can modify equation (2.7) so that the second term becomes the limit given by T i = 1 2 m i ( _ R i ) 2 + 1 2 lim ` i !1 ` i X j=1 m i j ( _ r i j ) 2 (2.9) = 1 2 m i ( _ R i ) 2 + 1 2 Z V i ( _ r i ) 2 dm; (2.10) where the integral is taken over the mass of body V i . Using equation (2.8), we obtain T i = 1 2 m i ( _ R i ) 2 + 1 2 Z V i (! i r i ) (! i r i )dm: (2.11) The second term in equation (2.11) denotes an integral over the i th rigid body where the j th particle is equivalent to a dierential element. The term inside the integral is written in matrix notation as (! i r i ) (! i r i ) = (! i ) T 2 6 6 6 6 6 4 ( i ) 2 + ( i ) 2 i i i i i i ( i ) 2 + ( i ) 2 i i i i i i ( i ) 2 + ( i ) 2 3 7 7 7 7 7 5 ! i ; (2.12) 12 where ! i = [! i ! i ! i ] T : It then follows that Z V i (! i r i ) (! i r i )dm = (! i ) T J i ! i ; (2.13) where J i is the symmetric inertia matrix for the i th rigid body given by J i = 2 6 6 6 6 6 4 J i J i J i J i J i J i J i J i J i 3 7 7 7 7 7 5 : (2.14) The moments of inertia are J i = Z ( i ) 2 + ( i ) 2 dm; J i = Z ( i ) 2 + ( i ) 2 dm; J i = Z ( i ) 2 + ( i ) 2 dm; (2.15) and the products of inertia are J i = Z i i dm; J i = Z i i dm; J i = Z i i dm: (2.16) The total kinetic energy for the i th rigid body is thus T i = 1 2 m i ( _ R i ) 2 + 1 2 (! i ) T J i ! i : (2.17) In summary, we see the i th rigid body contains a coordinate frame attached to it, which is xed at its center of mass, and thei th rigid body is free to rotate relative to an inertial coordinate frame. The total kinetic energy of the i th rigid body is then the sum of two parts: (1) the kinetic energy of a single particle that has a mass equal to the total mass of 13 the i th rigid body and moves with the velocity of the center of mass given by 1 2 m i ( _ R i ) 2 , and (2) the kinetic energy due to the motion of the rigid body relative to its center of mass given by 1 2 (! i ) T J i ! i . When the body axes of the i th rigid body are chosen as the principal axes the products of inertia become zero (J = J = J = 0). In this case, the symmetric inertia matrix J i of thei th rigid body given by equation (2.14) becomes a diagonal matrix containing the principal moments of inertia along the diagonal. 2.3 Multibody Coordinates The total kinetic energy of the i th independent and non-interacting rigid body given in equation (2.17) reveals that its motion is described by coordinates involving the velocities _ R i and the angular velocities ! i . The conguration of the i th rigid body is described by the location of its center of mass R i = R i 1 R i 2 R i 3 T and also by its orientation i = i 1 i 2 i 3 T relative to the inertial coordinate frame. Thus, the i th rigid body has a conguration space that is completely described by the 3-vectorR i and the 3-vector i |a total of 6 degrees of freedom|so that each rigid body contains 6 independent coordinates. For theN rigid body system, let us assume thatq i 2R n i is the generalized coordinate vector that describes the conguration of thei th rigid body. At this juncture, we make no assumptions on the nature of then i generalized coordinates for thei th rigid body; i.e., the generalized coordinates may or may not represent the minimum number of coordinates. Note that ! i is not a time derivative of any particular angle coordinate. 14 In general, we can represent the components ofR i in terms of the generalized coordinates by R i 1 = ^ f i 1 ( 1 ; 2 ;:::; l i;t) := ^ f i 1 (;t) (2.18) R i 2 = ^ f i 2 ( 1 ; 2 ;:::; l i;t) := ^ f i 2 (;t) (2.19) R i 3 = ^ f i 3 ( 1 ; 2 ;:::; l i;t) := ^ f i 3 (;t); (2.20) and for the components of i by i 1 = ^ g i 1 ( 1 ; 2 ;:::; s i;t) := ^ g i 1 (;t) (2.21) i 2 = ^ g i 2 ( 1 ; 2 ;:::; s i;t) := ^ g i 2 (;t) (2.22) i 3 = ^ g i 3 ( 1 ; 2 ;:::; s i;t) := ^ g i 3 (;t); (2.23) where l i 3 and s i 3. The generalized coordinate vector q i is thus composed of the coordinate vectors i 2R l i and i 2R s i so that the n i -vector q i = ( i ) T ; ( i ) T T , where n i = l i +s i . In general, given the chosen n i 6 generalized coordinates we know there must existn i 6 constraints due to the choice of coordinates since the i th rigid body has only 6 degrees of freedom. In the case that n i = 6, the generalized coordinates yield the minimum number of independent coordinates necessary to describe the conguration of the i th rigid body, and constraints are non-existent in its conguration prescription. Upon dierentiating R i with respect to time, we can obtain the components of the velocity _ R i by _ R i k = l i X j=1 @ ^ f i k @ i j _ i j + @ ^ f i k @t :=f i k ( i ; _ i ;t); k = 1; 2; 3: (2.24) 15 Similarly, the components of the velocity _ i are _ i k = s i X j=1 @^ g i k @ i j _ i j + @^ g i k @t := g i k ( i ; _ i ;t); k = 1; 2; 3: (2.25) The components of the absolute angular velocities ! i are then ! i k = h i ! ( i k ; _ i k ); k = 1; 2; 3; (2.26) = h i ! ^ g i k ( i ;t); g i k ( i ; _ i ;t) ; k = 1; 2; 3; (2.27) or ! i k =g i k ( i ; _ i ;t); k = 1; 2; 3: (2.28) Typically, the generalized coordinates i are given by the Cartesian coordinates x, y, andz, or by spherical coordinates, which are often represented by the componentsr,, and . There are also a number of commonly used generalized coordinates i , which include Euler angles, Rodrigues parameters, and unit quaternions. In the case of Euler angles and Rodrigues parameters, the orientation of the body-xed coordinate frame is described by three parameters, which is the minimal number of coordinates necessary to describe orientation. However, the choice of using the minimal number of coordinates comes at some cost since any minimal set of three rotational coordinates will become singular for certain orientations [21]. On the the other hand, the four parameter unit quaternion does not suer from this problem. But in this case, we have one more generalized coordinate than (independent) rotational degrees of freedom. As a result, an additional equation of constraint exists, and it must be satised in order to describe the orientation of the rigid 16 body. An extensive review of these and other possible rotational parameterizations, i , is given in [24]. 2.4 Equations of Motion forN Rigid Bodies To obtain the equations of motion for the i th rigid body let us rst consider its kinetic energy. Given the generalized coordinates in equations (2.24) and (2.28) the kinetic energy of the i th rigid body becomes T i = 1 2 m i (f i ) T f i + 1 2 (g i ) T J i g i (2.29) for f i = f i 1 f i 2 f i 3 T and g i = g i 1 g i 2 g i 3 T . The governing equations of motion are subsequently found by application of Lagrange's equation. Under the assumption that all the generalized coordinates are independent, Lagrange's equation becomes d dt @T i @ _ q i k @T i @q i k =Q i k k = 1; 2;:::;n i ; (2.30) where Q i 2 R n i is the generalized force vector (part of which may be derived from a potential) applied to the i th rigid body. The n i -vector Q i is composed of the generalized forces i 2 R l i and generalized torques i 2 R s i so that Q i = [( i ) T ; ( i ) T ] T . Thus, Lagrange's equation yields m i " @f i @ _ i k T d dt f i + d dt @f i @ _ i k @f i @ i k T f i # = i k ; k = 1; 2;:::;l i ; (2.31) 17 and @g i @ _ i k T J i d dt g i + d dt @g i @ _ i k @g i @ i k T J i g i = i k ; k = 1; 2;:::;s i : (2.32) In general, equations (2.31) and (2.32) yield a set ofn i second-order nonlinear dierential equations for thei th rigid body. Note that these equations, which are specied for the i th rigid body, may describe interactions due to conservative and/or non-conservative forces between other bodies in the N rigid body system. Since equations (2.31) and (2.32) are linear in the accelerations i and i , which are generated from the time derivatives, we can always obtain the equation M i (q i ; _ q i ;t) q i =F i (q i ; _ q i ;t); (2.33) where M i is a symmetric n i by n i matrix. In general, the symmetric matrix M i may be singular while the n i -vector F i consists of the 'given' generalized forces on the i th rigid body, where the given generalized forces are known functions of their arguments. The equations of motion of the system, which consist of theN rigid bodies, are then governed by M q := 2 6 6 6 6 6 6 6 6 6 4 M 1 0 0 0 M 2 0 . . . . . . . . . . . . 0 0 M N 3 7 7 7 7 7 7 7 7 7 5 2 6 6 6 6 6 6 6 6 6 4 q 1 q 2 . . . q N 3 7 7 7 7 7 7 7 7 7 5 = 2 6 6 6 6 6 6 6 6 6 4 F 1 F 2 . . . F N 3 7 7 7 7 7 7 7 7 7 5 :=F; (2.34) where the symmetric block diagonalK byK matrixM could be singular and the vectorF is aK-vector. The generalized acceleration vector a, whenM is not singular, is obtained as q =M 1 F :=a(q; _ q;t): (2.35) 18 The dimensions of the matrix M and the vector F are determined by K = N X i=1 n i ; (2.36) where K is the total number of generalized coordinates used. The minimum value of K, which is denoted by K min , determines the total number of degrees of freedom of the system. In deriving these equations, we have assumed the components of the generalized co- ordinate vector q i are independent for each of the N rigid bodies. When K = K min , the equations of motion for theN rigid body system are simply given by equation (2.34). When the components of any of the generalized coordinate vectorsq i are dependent on one another, equation (2.34) is incomplete because it ignores the requiredKK min constraint equations. Thus, to properly model the N rigid body system we must impose coordinate constraint requirements to equation (2.34). The coordinate constraint equations have the form ' k;m (q;t) = 0; k = 1; 2;:::;KK min ; (2.37) which make-up a general set of consistent KK min holonomic constraint equations for the N rigid body system. These additional requirements on the N rigid body system encompass what we shall call the modeling constraints. 19 2.5 Modeling a ComplexN Rigid Body System Thus far, we have assumed the generalized coordinates for the N rigid body system may be generally chosen so that KK min , which may or may not introduce the coordinate constraint equations in (2.37). In fact, thei th rigid body in the system may not represent a single rigid body, but instead a subsystem of multiple rigid bodies with an arbitrary number of degrees of freedom that has the form of equation (2.33). Additionally, the symmetric K by K matrix M in equation (2.34) may be singular since, in general, K K min , which may causeM to become rank decient. However, we have not yet described the completely modeled N rigid body system because the i th rigid body will generally have additional modeling requirements. To construct the desired complex N rigid body system, we now dene the interactions, or additional modeling constraints, between each of the rigid bodies in the system. Here, it is assumed a complex N rigid body system is a collection of rigid bodies that interact with one another wherein all the interactions are governed by any combination of the modeling constraints k;m (q;t) = 0; k = 1; 2;:::;h m (2.38) k;m (q; _ q;t) = 0; k = 1; 2;:::;h m ; (2.39) which are holonomic and nonholonomoic constraint equations, respectively. Equations (2.39) are called nonholonomic if they cannot be integrated and put in the form of equation (2.38). Thus, in addition to the constraints that arise due to the choice of generalized coordinates in equation (2.37), we also impose those modeling constraints that construct the desired complex N rigid body system. 20 The total number of modeling constraints for the complex N rigid body system is found by summing the modeling constraint equations (2.37), (2.38), and (2.39), and is given by h =KK min +h m +h m : (2.40) When the set of h modeling constraint equations are suciently smooth, we can obtain the general constraint matrix equation A m (q; _ q;t) q =b m (q; _ q;t); (2.41) by appropriately dierentiating the h modeling constraints with respect to time. Note that the general set of constraint equations in (2.41) may include constraints that are (1) nonlinear functions of position and velocity, (2) explicitly dependent on time, and (3) functionally dependent. Equation (2.41) thus includes all the consistent modeling constraints necessary to appropriately model general complex N rigid body systems. 2.6 Explicit Equations of Motion for ComplexN Rigid Body Systems To accommodate the modeling constraint requirements on the N rigid body system an additional force K-vector, F m , is required. This is the additional force that is applied to the N rigid body system in equation (2.34) so that the required motion, which satises the modeling constraints in equation (2.41), is obtained. First, we shall assume the 21 symmetric K by K matrix M in equation (2.34) is nonsingular. The equation of motion for the properly modeled complex N rigid body system is then expressed as M q =F (q; _ q;t) +F m (q; _ q;t); (2.42) and the corresponding acceleration a m is given by q =M 1 (F +F m ) :=a m (q; _ q;t): (2.43) The constraint force F m associated with the modeling constraints in equation (2.41) is determined explicitly by Udwadia and Kalaba [27] and is given by F m =M 1 2 (A m M 1 2 ) + (b m A m a): (2.44) The '+' symbol denotes the Moore-Penrose generalized inverse of the matrix B m = A m M 1 2 . The explicit acceleration of the complex N rigid body system is then q =a +M 1 2 B + m (b m A m a): (2.45) The explicit determination of the acceleration in equation (2.45) conforms to Gauss's prin- ciple of least constraint [8], which requires the global minimization of the scalar quantity G = ( qa) T M( qa): (2.46) 22 Thus, the acceleration in equation (2.45) is the unique acceleration that minimizes the scalar in equation (2.46) and yields the acceleration of the appropriately modeled complex N rigid body system at each instant of time. Since the theory allows K > K min , the symmetric K by K matrix M could be, in general, singular. Equation (2.45) is undened when M is singular. The explicit equations of motion for constrained mechanical systems with singular matrices M are obtained in [29]. When, in particular, the matrix M A = 2 6 4 M A m 3 7 5 has full rank, the following alternative approach can be used. Consider the general modeling constraint equation (2.41) and pre-multiply to it the matrix A + m so that A + m A m q =A + m b: (2.47) By adding equation (2.34) to equation (2.47), we obtain the augmented system M q := (M +A + m A m ) q =F +A + m b m := F: (2.48) Since the K by K matrix M is symmetric and positive denite (see Appendix A), the acceleration is a = M 1 F: (2.49) Note that equation (2.48) is simply the unconstrained equation of motion of the complex N rigid body system. When the modeling constraints (2.41) are applied to this system, the additional K-vector F m needs to be determined so that M q = F + F m (2.50) 23 exactly satises the modeling constraints, where q = M 1 ( F + F m ) := a m : (2.51) The additional K-vector F m is explicitly given by F m = M 1 2 (A m M 1 2 ) + (b m A m a): (2.52) The explicit acceleration of the modeled multibody system when M is singular is then found by q = a + M 1 2 B + m (b m A m a); (2.53) where B m =A m M 1 2 . When M is symmetric and positive denite (nonsingular), equa- tion (2.53) is equivalent to the explicit acceleration found by equation (2.45). Thus, the explicit acceleration of the modeled multibody system is uniquely found when the matrix M A has full rank whetherM is singular or not. To motivate the development thus far let us consider the following example. 2.7 Example: Two Spacecraft Linked at an Arbitrary Point Consider a multibody spacecraft system in orbit around a central body with a uniform gravity eld. The multibody spacecraft system consists of N = 2 spacecraft elastically linked together by a nonlinear spring at arbitrary points P i , i = 1; 2; as shown in Figure (2.3). The two spacecraft are free to rotate independently about the line P 1 P 2 , which is xed in direction relative to V 1 and V 2 . The relative distancekk is free to vary so that relative motion between the spacecraft V 1 and V 2 occurs when the spring stretches 24 1 η 2 η 1 μ 1 γ 2 γ 2 μ 2 ρ 1 ρ 2 P 1 P 2 P r 1 P r i R 1 i R 2 i R 3 nl l k k , 1 V 2 V 1 1 , J m 2 2 , J m Figure 2.3: A two spacecraft system|denoted by rigid bodies V 1 and V 2 |elastically linked together at the pointsP 1 andP 2 by a spring whose constantsk l andk nl denote the linear and cubically nonlinear restoring forces exerted by the spring. The two spacecraft are free to rotate independently about the line P 1 P 2 . or compresses along the line P 1 P 2 that is xed relative to both V 1 andV 2 . The relative distance vector is given by = 1 +r 1 P 2 r 2 P ; (2.54) where the position vectors r i P ; i = 1; 2; are xed in their corresponding rotating body- xed coordinate frames. This multibody system has 8 degrees of freedom since relative motion only occurs in translation along the line P 1 P 2 and in rotation about the line P 1 P 2 . Without loss of generality, let us assume the desired generalized coordinate vector for each spacecraft is taken as q i = i 1 i 2 i 3 i 1 i 2 i 3 i 4 T ; i = 1; 2; (2.55) 25 so that n 1 = n 2 = 7. The orbital position of the i th spacecraft is represented by the coordinates i 1 ; i 2 ; and i 3 , and the attitude of the i th spacecraft is represented by the coordinates i 1 , i 2 , i 3 , and i 4 . Recall that these generalized coordinates came about by the transformations ^ f i and ^ g i . For the sake of generality, we'll assume ^ f i ; ^ g i ; and h i ! are known, which implies the transformations f i and g i given in equations (2.24) and (2.28) are also known. The uniform gravitational potential of the i th spacecraft is U i g ( i ) = g m i k i k ; i = 1; 2; (2.56) where kk denotes the 2-norm operation and g |usually denoted the gravitational parameter|is the product of the gravitational constant and the total mass of the central body. Let us assume the elastic potential of the nonlinear spring is given by U e (q) = 1 2 k l (kkL e ) 2 + 1 4 k nl (kkL e ) 4 ; (2.57) where the distance between V 1 and V 2 when the spring is unstretched is L e . The gener- alized forces on the i th spacecraft are thus i = @U i g @ i @U e @ i + i nc ; i = 1; 2; (2.58) where i nc are the non-conservative forces. Similarly, the generalized torques on the i th spacecraft are i = @U e @ i + i nc ; i = 1; 2; (2.59) 26 where i nc are the non-conservative torques. In this particular multibody spacecraft prob- lem, we shall assume i nc = i nc = 0. Assuming the generalized coordinates in equation (2.55) are independent, the equa- tions of motion of the two spacecraftV 1 andV 2 are found by equations (2.31) and (2.32), which subsequently yields M q := 2 6 4 M 1 0 0 M 2 3 7 5 2 6 4 q 1 q 2 3 7 5 = 2 6 4 F 1 F 2 3 7 5 :=F: (2.60) Here, we shall assume the 14 by 14 matrixM is singular and the vectorF is a column 14- vector sinceK = 14. As a result, we have 2 coordinate constraint equations (KK min = 2) given by ' 1;m ( 1 ) = 0 (2.61) ' 2;m ( 2 ) = 0: (2.62) The nature of the modeling constraints in equations (2.61) and (2.62) will be described in detail in chapter 4. The set of equations in (2.60) describes the un-modeled motion of the two interacting spacecraft. The interaction between the spacecraft is governed by the elastic force generated from the spring that connects the two spacecraft. The multibody spacecraft system (2.60) is not correctly modeled because, in addi- tion to the modeling constraints (2.61) and (2.62), there exist the additional interactions between the bodies V 1 and V 2 due to physical modeling constraints. To model these physical interactions, we have the requirement in which the orientation of the relative 27 distance vector is xed relative to the rotating body-xed coordinate frames of V 1 and V 2 . This requirement is given by the following two modeling constraint equations k;m (q 1 ;q 2 ) := ^ n 1 = 0; k = 1; 2; 3; (2.63) k;m (q 1 ;q 2 ) := ^ n 2 = 0; k = 4; 5; 6; (2.64) which demands that the body-xed unit 3-vectors ^ n 1 and ^ n 2 remain xed in orientation relative to the line P 1 P 2 . The unit vectors ^ n 1 and ^ n 2 are determined by the location of the spring connections toV 1 andV 2 . It then follows that equations (2.61), (2.62), (2.63), and (2.64)|a total ofh = 8 constraints|completely describe the modeling constraints for this multibody spacecraft problem. Notice that we have chosen more modeling constraints (h = 8) than necessary. The minimum number of modeling constraints necessary to model the system is 6 since we know the system has 8 degrees of freedom. This is an important feature in the formulation because it allows the freedom to choose any desirable number of consistent modeling constraint equations without any concern regarding dependency. Upon appropriately dierentiating theh = 8 modeling constraints with respect to time, we can write the modeling constraint matrix equation as in (2.41). The explicit acceleration of the appropriately modeled two spacecraft system is then found by equation (2.53). 2.8 Summary A general formulation for modeling complex N rigid body systems has been developed. The multibody system is modeled by rst considering the general motion of the i th rigid body in the multibody system. The individual rigid bodies are then put together by dening the interactions, or modeling constraints, that build the multibody system. The 28 conguration space for thei th rigid body is freely chosen by then i generalized coordinates. Furthermore, the design of the multibody system is simple in its implementation since there is no restriction on the total number of the general modeling constraints except that they be consistent. The explicit equations of motion of the modeled multibody system given by M q = F + F m = F + M 1 2 B + m (b m A m a) are simple to construct and implement numerically. The forces that arise due to the presence of the modeling constraints are found explicitly in closed form and can be used to determine the eect they have on the motion of the modeled multibody system. The illustrative example in section 2.7 shows the actual model development of a realistic multibody spacecraft problem, which is seen as a set of elastically connected rigid bodies subjected to the modeling constraints A m q =b m : The whole procedure of modeling a complex multibody system is thus reduced to simply applying Lagrange's equation and subsequently dening those modeling constraints that appropriately congure the desired N rigid body system. In the next chapter, we will extend this formulation to develop a general methodology for nding, in closed form, the control forces necessary for the complex N rigid body system to exactly satisfy a set of control objectives. 29 Chapter 3 Explicit Control for Multibody Systems In the previous chapter, a general approach to modeling a complex system of N rigid bodies was developed. The equations of motion of the modeled multibody system con- sist of the given generalized force K-vector F and the additional modeling constraint force K-vector F m , which is appropriately designed, to model the multibody system. In this chapter, the general nonlinear control problem of the modeled multibody system is investigated and a general approach for explicitly generating nonlinear controllers is developed. 3.1 Multibody Control Formulation The goal in this section is to develop a general framework for controlling the modeled multibody system in equation (2.50). The basic idea is to nd a control force F c 2R K so that the controlled system M(q; _ q;t) q = F (q; _ q;t) + F m (q; _ q;t) + F c (q; _ q;t) (3.1) 30 satises some control objective for the modeled multibody system. The rst thing to note is that equation (3.1) is a general nonlinear control problem with no restrictions on its structure such as being a control-ane system or utilizing any approximations. Typical control objectives for multibody systems might include stabilization, tracking, and disturbance attenuation. While the nonlinear control problem considered here is general, the additional force F c is restricted only in the sense that it cannot violate any of the modeling constraints ' m , m , and m , which are inherent to the modeling force F m . Thus, the modeling constraints given in equation (2.41) must also be imposed to the design of the control force F c . To investigate the structure of the control force F c , let us rst consider the modeled multibody system given by equation (2.50). Substituting equations (2.49) and (2.52) into equation (2.50), we can write the explicit equations of motion for the modeled multibody system as M q = F + M 1 2 B + m (b m A m M 1 F ) (3.2) = M 1 2 B + m b m + (I M 1 2 B + m B m M 1 2 ) F: (3.3) = M 1 2 B + m b m + M 1 2 (I B + m B m ) M 1 2 F; (3.4) where the matrixI is aK byK identity matrix. In equation (3.4), observe that the total force of the modeled multibody system is made-up of the following two K-vectors F 1 = M 1 2 B + m b m (3.5) F 2 = M 1 2 (I B + m B m ) M 1 2 F = P F: (3.6) 31 By inspection of equation (3.5), we see the K-vector F 1 belongs to the column space of the K by h matrix M 1 2 B + m . The force F 1 arises solely due to the modeling constraints imposed on theN rigid body system. TheK-vectorF 2 in equation (3.6) shows the eect the modeling constraints have on theK-vector F . Consider now the following denitions. Denition 1 Moore-Penrose Inverse [28]. The unique matrix B + m that satises the following four conditions is called the Moore- Penrose inverse of the h by K matrix B m . 1. ( B + m B m ) T = B + m B m . 2. ( B m B + m ) T = B m B + m . 3. B m B + m B m = B m . 4. B + m B m B + m = B + m . Denition 2 Idempotent Matrix [28]. A square matrix W satisfying W =W 2 is dened to be an idempotent matrix. Denition 3 Projection Matrix [1]. A square matrix W is a projection matrix if and only if it is idempotent. Let the idempotent matrix W project any vector X n to a vector subspaceV s , whereV s is the space spanned by the columns of W . The projection is orthogonal if W is a symmetric matrix, or it is oblique if W is a nonsymmetric matrix. We now state the following results. 32 Corollary 1 TheK-vectorF 2 in equation (3.6) is the projection of theK-vector F by the projection matrix P. When B + m B m M = M B + m B m the projection is orthogonal, otherwise the projection is oblique. Proof. The matrix P is a projection matrix since we can show it is idempotent by P 2 = h M 1 2 (I B + m B m ) M 1 2 ih M 1 2 (I B + m B m ) M 1 2 i = M 1 2 (I 2 B + m B m + B + m B m B + m B m ) M 1 2 = M 1 2 (I B + m B m ) M 1 2 = P; (3.7) where we have used conditions 3 and 4 of the Moore-Penrose inverse. Additionally, we can show the matrix P is generally nonsymmetric by P T = M 1 2 (I (B + m B m ) T )M 1 2 = M 1 2 (IB + m B m )M 1 2 ; (3.8) where we have used condition 1 of the Moore-Penrose inverse. Thus, the matrix P is, in general, the oblique projection operator on the K-vector F . However, if it is true that M 1 2 (I B + m B m ) M 1 2 = M 1 2 (I B + m B m ) M 1 2 (I B + m B m ) M = M(I B + m B m ) B + m B m M = M B + m B m ; (3.9) 33 then the matrix P becomes the orthogonal projection operator on the K-vector F . Corollary 2 The K-vector F 2 in equation (3.6) belongs to the null space of A m M 1 . Proof. Pre-multiplying F 2 by the h by K matrix B m M 1 2 we have B m M 1 2 F 2 = B m M 1 2 h M 1 2 I B + m B m M 1 2 F i = B m I B + m B m M 1 2 F = 0; (3.10) where condition 3 of the Moore-Penrose inverse is used. The result follows since the matrix B m M 1 2 =A m M 1 . Corollary 3 The K-vectors F 1 and F 2 are M 1 -orthogonal. Proof. The vectors F 1 and F 2 are M 1 -orthogonal if the quadric F T 1 M 1 F 2 is identi- cally zero. Thus, we have F T 1 M 1 F 2 = h b T m ( B + m ) T M 1 2 i M 1 h M 1 2 I B + m B m M 1 2 F i = b T m ( B + m ) T I ( B m ) T ( B + m ) T M 1 2 F = 0; (3.11) where conditions 1 and 4 of the Moore-Penrose inverse are used. 34 The above results show that the motion of the modeled multibody system in equation (3.4) is governed by a total force that is composed of the two M 1 -orthogonalK-vectors F 1 and F 2 . As already stated, the vector F 1 is the part of the total force that lies in the column space of the matrix M 1 2 B + m and arises as a result of the modeling constraints imposed on the N rigid body system. The vector F 2 is the projection of the K-vector F to the null space of A m M 1 . Suppose now we add any K-vector ^ F c to equation (3.4) so that M q = M 1 2 B + m b m + P F + ^ F c ; (3.12) or M q = ^ F + ^ F c ; (3.13) where ^ F = F 1 +F 2 . However, we cannot add any vector we wish to equation (3.4) since equation (3.13) must satisfy the modeling constraints given by equation (2.41). The explicit equations of motion when equation (3.13) is constrained by equation (2.41) are thus found by using equation (2.52) so that M q = ^ F + ^ F c + M 1 2 B + m b m B m M 1 2 ( ^ F + ^ F c ) (3.14) = ^ F + ^ F c + M 1 2 B + m b m M 1 2 B + m B m M 1 2 ( ^ F + ^ F c ) (3.15) = M 1 2 B + m b m + P ^ F + ^ F c : (3.16) Corollary 4 The K-vector F 2 is equal to the K-vector P ^ F . 35 Proof. Substituting ^ F =F 1 +F 2 into P ^ F , we have P ^ F = P M 1 2 B + m b m + P 2 F: (3.17) The rst term on the right hand side of equation (3.17) becomes P M 1 2 B + m b m = h M 1 2 (I B + m B m ) M 1 2 i M 1 2 B + m b m = M 1 2 (I B + m B m ) B + m b m = 0; (3.18) where condition 4 of the Moore-Penrose inverse is used. The second term on the right hand side of equation (3.17) is simplied as P 2 F = P F (3.19) since P is idempotent. Thus, we nd that F 2 = P ^ F . As a result of Corollary 4, we have M q = F 1 + P ^ F + ^ F c (3.20) = F 1 +F 2 + P ^ F c : (3.21) Equation (3.21) is equivalent to M q = F + F m + P ^ F c ; (3.22) 36 where ^ F c is an arbitrary K-vector. Therefore, we nd that the permissible control forces of the modeled multibody system must have the form F c = P ^ F c ; (3.23) where F c are thoseK-vectors in the same direction asF 2 that belong to the null space of A m M 1 . 3.2 Trajectory Control Given the additional K-vector F c (the controller), we are now interested in designing a controller to satisfy some control objective for the modeled multibody system. The control objectives, or trajectory requirements, considered herein are expressible by any combination of the general constraint equations k;c (q;t) = 0; k = 1; 2;:::;h c (3.24) k;c (q; _ q;t) = 0; k = 1; 2;:::;h c ; (3.25) which again include constraint equations of the holonomic and nonholonomic variety. Upon appropriately dierentiating these constraints with respect to time, we can obtain the constraint matrix equation A c (q; _ q;t) q =b c (q; _ q;t): (3.26) 37 Thus, the control objectives for the modeled multibody system are represented by what we shall call the control constraints in equation (3.26). To satisfy these imposed control constraints, we obtain the arbitrary K-vector ^ F c explicitly by ^ F c = M 1 2 B + c (b c A c a m ); (3.27) where B c =A c M 1 2 . The control force applied to the modeled multibody system is thus provided by the K-vector F c so that F c = P M 1 2 B + c (b c A c a m ): (3.28) Note that the control constraints are prescribed independently from the modeling con- straints. In general, the control constraints c and c may prescribe trajectories that are physically unrealizable for the modeled multibody system; i.e., trajectories that do not lie on the modeling constraint manifold may be prescribed. In this case, the control force F c obtained in equation (3.28) will provide the control action that adjusts the desired trajectory (equation (3.26)) to actually occur along the modeling constraint manifold, which is given by the intersection of the modeling constraints ' m , m , and m . This is possible since the K-vector F c lies in the null space of A m M 1 and is thus consistent with all the modeling constraints. For example, a trajectory, (q c ; _ q c ), may be prescribed to the multibody system by equation (3.26) such that a desired end point is achieved, and the exact path taken, (q; _ q), does not always lie on the control constraints c and c as illustrated in Figure (3.1). This simplies the controller design for complex multibody systems since it is not necessary to prescribe the control constraints along the modeling 38 ) ( ), ( f f t q t q & ) ( ), ( 0 0 t q t q & q q & , c c q q & , m m m Ψ ∩ Φ ∩ φ c F c F c F c F ˆ c F ˆ c F ˆ Figure 3.1: A prescribed trajectory (q c ; _ q c ) and the actual trajectory (q; _ q) over a time duration t2 (t 0 ;t f ) along the modeling constraint manifold dened by ' m \ m \ m . constraint manifold. However, it is necessary to ensure that the desired control constraints c and c can properly intersect the modeling constraint manifold. When the control force F c = ^ F c , the K-vector ^ F c is projected to itself. In this case, the control force F c yields the control action that exactly satises the desired control objectives at each instant of time. This indicates that the prescribed trajectories c and c lie exactly on the modeling constraint manifold. When this occurs we note that the control force F c not only causes the multibody system to exactly track the desired trajectory, but it does so by minimizing the cost function [26] C(t) = F T c M 1 F c (3.29) at each instant of time. 39 3.3 Stabilization The control methodology developed thus far assumes the initial conditions q(0) and _ q(0) will exactly satisfy the control constraints. In reality, it is usually necessary to rst control the system to the desired trajectory and subsequently along the desired trajectory. Furthermore, we cannot expect that the initial states of the multibody system will lie exactly on the control constraints. To deal with errors between the actual trajectory and the desired trajectory, the controller F c must provide stability along the control constraints. A simple approach to create stability along the control constraints c and c is ac- complished by the modications c + _ c + c = 0 (3.30) _ c + c = 0 (3.31) for holonomic and nonholonomic constraints, respectively. By appropriately selecting the control parameters ; ; and , we can ensure the xed points ( c ; _ c ) = (0; 0) and c = 0 are asymptotically stable. In fact, instead of equations (3.30) and (3.31), we could just as well choose other suitable rst and second-order linear, or nonlinear, dierential equations that contain the desired asymptotically stable xed points, and possess the salient features sought for converging to the control constraints. For example, instead of equation (3.30), we could choose the nonlinear damped oscillator c + _ c + sin c = 0: (3.32) 40 In either case, we can obtain the modied control constraints in the form of an appropri- ately constructed control constraint matrix equation (3.26). The idea behind enforcing these modied control constraints with the control force F c is to cause the multibody system to approach the desired trajectory as illustrated in Figure (3.2), where the actual ) ( ), ( 0 0 t q t q & q q & , c c q q & , ) ( ), ( 0 0 t q t q c c & m m m Ψ ∩ Φ ∩ φ c c Ψ ∩ Φ c F ˆ c F ˆ c F c F c c F F ˆ = Figure 3.2: An actual trajectory (q; _ q) converging along the modeling constraint manifold ' m \ m \ m to the desired trajectory (q c ; _ q c ), which travels along the control constraint manifold dened by c \ c . dynamical path taken is governed by the nature of the modied control constraints. While Figure (3.2) shows the desired trajectories (q c ; _ q c ) do not initially lie on the modeling con- straint manifold, we note that they would lie on ' m \ m \ m when c and c are chosen so that ^ F c = P ^ F c . In linear systems, the origin of a stabilized system is always globally asymptotically stable since there exists only one xed point. However, in nonlinear systems there may be additional xed points, or trajectories, which are also asymptotically stable. Moreover, periodic solutions may exist somewhere in the phase space of the system. Thus, depending 41 on the multibody system under investigation further stability analysis may be required to determine the level of stability achieved when F c stabilizes the multibody system under the appropriately selected modied control constraints. The modied control constraints (3.30) and (3.31) will at least achieve local stabilization for the nonlinear system since the xed points ( c ; _ c ) = (0; 0) and c = 0 are asymptotically stable with the appropriate choice of , , and . 3.4 Underactuated Control The general set of controllers obtained by F c presume the multibody system under inves- tigation is fully actuated. In many mechanical systems, it is not always feasible to design a system such that all of its degrees of freedom are actuated. In this section, we present a method to underactuate the modeled multibody system (3.4) thereby creating fewer actuators than degrees of freedom. Consider now a control force vector F a;c 2R Kc , where K c K min denotes the dimension of F a;c . A fully actuated system will have an equal number of control inputs and independent generalized coordinates, i.e., K c = K min . In the following, we shall assume the multibody system given by equation (3.1) is underac- tuated whenK c <K min , and under this assumption we will haveK min K c second-order dynamic equations of the form k;m (q; _ q; q;t) = 0 k = 1; 2;:::;K min K c ; (3.33) which are linear in the accelerations q. Equation (3.33) is called a second-order nonholo- nomic constraint equation if there exist no functions such that = m (q; _ q; q;t). The 42 modeling constraints in equation (3.33) are obtainable, in certain cases, by partitioning equation (3.1) so that 2 6 4 M a 0 0 M u 3 7 5 2 6 4 q a q u 3 7 5 2 6 4 F a F u 3 7 5 2 6 4 F a;m F u;m 3 7 5 = 2 6 4 F a;c 0 3 7 5; (3.34) where the subscripts 'a' and 'u' denote the actuated and unactuated subsystems, respec- tively. It follows that the modeling constraints (3.33) are given by the unactuated part of equation (3.34). However, it is not always possible to obtain the necessary modeling constraints (3.33) directly from equation (3.34). They may lie in a transformed set of coordinates; i.e., coordinates relative to a body-xed coordinate frame. Departing from multibody appli- cations for the time being, let us consider the surface vessel vehicle (N = 1) in Figure 3.3. The vehicle is shown relative to an inertial coordinate frame (x;y), where the vectors u 1 x y θ 1 u 2 u J m, Figure 3.3: A surface vessel vehicle with mass m and inertia J. 43 (surge) and u 2 (sway) show the directions of its linear velocities. The angle determines the orientation of the vehicle. Typically these types of vehicles are designed with no side thruster such that they can only actuate in the direction of u 1 and about the angle . The vehicle is thus underactuated since it cannot actuate in the u 2 direction. The general dynamics of a surface vessel vehicle can be found in [7], where the equa- tions of motion are given by M _ +C() +D =F (3.35) _ q =T (): (3.36) The 3-vector = [u 1 ;u 2 ;#] T denotes the linear velocitiesu 1 ,u 2 , and the angular velocity # in the vehicles body-xed coordinate frame. The 3-vectorq = [x;y;] T gives the vehicles position relative to the inertial coordinate frame. The constant diagonal matrices M = diag(J 1 ;J 2 ;J 3 ) andD = diag(d 1 ;d 2 ;d 3 ) yield inertia properties and damping coecients, respectively. The matrixC contains terms associated with Coriolis and centripetal forces given by C = 2 6 6 6 6 6 4 0 0 J 2 u 2 0 0 J 1 u 1 J 2 u 2 J 1 u 1 0 3 7 7 7 7 7 5 ; (3.37) and the matrixT is simply the transformation matrix T = 2 6 6 6 6 6 4 cos sin 0 sin cos 0 0 0 1 3 7 7 7 7 7 5 : (3.38) 44 To control the vehicle, the control inputs are given by the 3-vector F = [F u 1 ; F u 2 ; F # ] T . Without loss of generality, let us assume the matrix D = 0 and that the vehicle is symmetric with respect to the axes u 1 and u 2 so that J 1 = J 2 = m and J 3 = J. The equations of motion for a fully actuated surface vessel vehicle then become m _ u 1 = F u 1 +mu 2 # (3.39) m _ u 2 = F u 2 mu 1 # (3.40) J _ # = F # : (3.41) Since the vehicle is underactuated (it cannot actuate along the u 2 direction), the com- ponent F u 2 is zero. Thus, the modeling constraint associated with the underactuation is given by m := _ u 2 +u 1 # = 0; (3.42) which is a second-order nonholonomic constraint. In inertial coordinates, equations (3.39)- (3.41) become M q := 2 6 6 6 6 6 4 m 0 0 0 m 0 0 0 J 3 7 7 7 7 7 5 2 6 6 6 6 6 4 x y 3 7 7 7 7 7 5 = 2 6 6 6 6 6 4 c;x c;y c; 3 7 7 7 7 7 5 :=F c ; (3.43) and the modeling constraint (3.42) becomes m := sin x + cos y = 0: (3.44) 45 The bar notation above M and F c is excluded here since the matrix M is nonsingular. The modeling constraint matrix equation associated with equation (3.44) is then A m q := sin cos 0 2 6 6 6 6 6 4 x y 3 7 7 7 7 7 5 = 0 :=b m : (3.45) If we are to control the vehicle under this modeling constraint, then by equation (3.23) we obtain the controlled system of equations by M q =P ^ F c ; (3.46) where F c =P ^ F c represents all the permissible controllers for the underactuated vehicle. The matrixP is computed as P =M 1 2 (IB + m B m )M 1 2 = 2 6 6 6 6 6 4 1 sin 2 sin cos 0 sin cos 1 cos 2 0 0 0 1 3 7 7 7 7 7 5 ; (3.47) and the 3-vector ^ F c is arbitrary. Let us now assume the underactuated vehicle, which starts from rest, is required to converge to and move along an ellipse of the form 1;c := x 2 e 2 x + y 2 e 2 y 1 = 0; (3.48) 46 wheree x ande y give its semimajor and semiminor axes. Equation (3.48) then represents a control constraint. To approach this control constraint, we obtain the modied control constraint equation 1;c + 1 _ 1;c + 1 1;c = 0; (3.49) which yields the asymptotically stable xed point ( 1;c ; _ 1;c ) = (0; 0) when the parameters 1 > 0 and 1 > 0. The control constraint matrix equation is then found as A c q := x=e 2 x y=e 2 y 0 2 6 6 6 6 6 4 x y 3 7 7 7 7 7 5 = _ x e 2 x _ y e 2 y 1 _ 1;c 1 1;c :=b c : (3.50) The force ^ F c that exactly satises the desired trajectory (3.49) is given as in equation (3.27), where the acceleration a m = 0. Thus, we can explicitly compute the equations of motion of the controlled system (3.46) as x = x e 2 x (1 sin 2 ) + y e 2 y sin cos (3.51) y = x e 2 x sin cos + y e 2 y (1 cos 2 ) (3.52) = 0; (3.53) where = bc x 2 =e 4 x +y 2 =e 4 y . We immediately see this system of equations is undened at the position (x;y) = (0; 0), which arises because the origin is a degenerate form of the ellipse in equation (3.48). Evaluating equations (3.51) - (3.53) at ( q; _ q) = (0; 0), and after some 47 simplication, we nd that this controlled system contains non-isolated xed points given by the relations x 2 e 2 x + y 2 e 2 y 1 = 0 (3.54) x e 2 x cos + y e 2 y sin = 0: (3.55) Equation (3.54) is exactly the ellipse we are attempting to control the vehicle along (equation (3.48)) while equation (3.55) yields orientations that occur when the vehicle axis u 2 is normal to the ellipse at position (x;y). We also see by inspection of equation (3.53) that this controller does not provide any angle control in the coordinate , and it can only actuate the vehicle along the u 1 direction making it impossible to move in a closed elliptical orbit. To move along the ellipse, we need to control the orientation of the vehicle so that it can appropriately align the actuating axis u 1 to apply the necessary centripetal forces. Here, we presume the nominal orientation of the vehicle along the ellipse path is given by 2;c := @ 1;c @x sin @ 1;c @y cos = 0; (3.56) which yields orientations so that the axis u 1 is pointed normal to the ellipse path. The modied control constraint that approaches this control objective is then given by 2;c + 2 _ 2;c + 2 2;c = 0; (3.57) 48 and the control constraint matrix equation becomes A c q := 2 6 4 x=e 2 x y=e 2 y 0 sin=e 2 x cos=e 2 y x cos=e 2 x +y sin=e 2 y 3 7 5 2 6 6 6 6 6 4 x y 3 7 7 7 7 7 5 = 2 6 4 _ x 2 =e 2 x _ y 2 =e 2 y 1 _ 1;c 1 1;c b 2;c 2 _ 2;c 2 2;c 3 7 5 :=b c ; (3.58) where b 2;c = x _ 2 sin 2 _ x _ cos =e 2 x y _ 2 cos + 2 _ y _ sin =e 2 y . The explicit equa- tions of motion of the controlled vehicle are again obtained by equation (3.46). The numerical integration of the resulting controlled equations is carried out to study the vehicle trajectories using a variable time step Runge-Kutta scheme with a relative error tolerance 10 12 and an absolute error tolerance 10 15 . We assume the vehicle starts at rest with a mass m = 30 kg and an inertia J = 50 kg m 2 . The desired ellipse path is chosen to have a semimajor axis e x = 1:5 m and a semiminor axis e y = 1 m, which corresponds to an eccentricity e = 0:745. The control parameters for this example are 1 = 2 = 2= p 10 and 1 = 2 = 1=10. Figure 3.4 shows various vehicle trajectories with initial positions interior and exterior to the desired ellipse path. Depending on the initial orientation (0) for a given initial position x(0) and y(0), the vehicle may or may not converge to and move along the ellipse. A single trajectory, which starts at the position (x;y) = (0; 1:5) m, is shown to diverge from the ellipse path while the other trajectories appear to converge to and along the ellipse path. In any case, stability along the ellipse path remains to be ascertained. 49 −2 −1 0 1 2 −2 −1.5 −1 −0.5 0 0.5 1 1.5 2 x [m] y [m] Figure 3.4: Surface vessel vehicle trajectories (dotted lines) near the desired ellipse path with eccentricity e = 0:745. The control forces F u 1 and F u 2 in the body-xed coordinate frame, and the control torque F # required to converge to and move along the ellipse for the initial conguration (x(0);y(0);(0)) = (1; 0; 0:8) are shown in Figure 3.5 for approximately 3 orbits. Note that the control force F u 2 is exactly zero as required by the modeling constraint (3.42). Figure 3.5(a) gives the control action required to approach the ellipse path while Figure 3.5(b) gives the control action to maintain motion along the ellipse path. Finally, Figure 3.6 shows the satisfaction of the control constraints 1;c and 2;c while the underactuated vehicle moves along the ellipse. The error in the control constraint 1;c as shown in Figure 3.6(a) reveals that the vehicle is periodically perturbed from the ellipse path. By examining Figure 3.6(b), we see this occurs when the control constraint 2;c deviates from zero. The peaks in Figure 3.6(b) correspond to vehicle positions near the positions (x;y) = (1:5; 0)m while the zero crossings correspond to vehicle positions (x;y) = (0;1)m. 50 0 5 10 15 20 25 −5 −4 −3 −2 −1 0 1 time [s] Control Forces F u1 [N] F u2 [N] F ϑ [N m] (a) Control forcesFu 1 andFu 2 , and control torqueF # over the time intervalt2 (0; 25)s. 500 1000 1500 2000 2500 3000 −0.015 −0.01 −0.005 0 0.005 0.01 0.015 time [s] Control Forces F u1 [N] F u2 [N] F ϑ [N m] (b) Control forces Fu 1 and Fu 2 , and con- trol torque F # over the time interval t 2 (25; 3000)s. Figure 3.5: Control action required to approach and move along the desired ellipse path. 500 1000 1500 2000 2500 3000 −0.5 0 0.5 1 1.5 2 2.5 3 x 10 −6 time [s] Φ 1,c (a) Error in control constraint 1;c over the time interval t2 (50; 3000)s. 500 1000 1500 2000 2500 3000 −3 −2 −1 0 1 2 3 x 10 −3 time [s] Φ 2,c (b) Error in control constraint 2;c over the time interval t2 (50; 3000)s. Figure 3.6: Error in control constraints 1;c and 2;c while moving along the ellipse. 51 The errors in Figure 3.6(b) correspond to corrections in the vehicles steering control law imposed by equation (3.56). 3.5 Summary In this chapter, we have derived the general form of the controller F c for the modeled multibody system. Since the multibody system is constructed using the h modeling con- straints in equation (2.41), we must ensure that they are not violated when any additional control forces F c are applied to the multibody system. The permissible control forces given by F c = P ^ F c may be designed using the arbitraryK-vector ^ F c . This arbitrary vector is then projected to the null space of A m M 1 , which is in the same direction as F 2 . When the arbitrary vector ^ F c is projected to itself by P, the control constraints c and c were seen to lie exactly on the modeling constraint manifold. The idea of a control constraint that is imposed to the modeled multibody system was presented. The control constraints are interpreted as the control objective, and to satisfy the control objective we select the arbitrary K-vector ^ F c by again using the fundamental equation so that ^ F c = M 1 2 B + c (b c A c a m ): To create stability along the control constraints, we presented the modied control con- straints. In particular, the control constraints in equations (3.24) and (3.25) are modied to create asymptotic stability at the xed points ( c ; _ c ) = (0; 0) and c = 0. This is 52 accomplished by appropriately choosing a rst and second-order linear, or nonlinear, dif- ferential equation that contains these asymptotically stable xed points. Furthermore, we saw that, initially, the control constraints need not lie exactly on the modeling constraint manifold, and we may prescribe c and c independently of the modeling constraints' m , m , and m . Since the controlled system resulting from a stabilizing controller F c is, in general, nonlinear, notions of stability beyond local stabilization need to be ascertained. The application of the control methodology developed herein to the control of underac- tuated systems was presented. It was shown that a multibody system is underactuated by introducing additional modeling constraints of the form m (q; _ q; q;t) = 0. The underactu- ated control of a surface vessel vehicle required to move along an ellipse was developed to show the ease with which controllers can be designed and studied using the approach. It also points out that a proper conceptualization of the desired control constraint manifold is required when trying to achieve a particular control objective. In the following chapter, we will utilize the modeling methodology developed in chapter 2, and the control methodology developed in this chapter, to develop a general approach for the rotational dynamics and control of rigid bodies using quaternions. 53 Chapter 4 Rotational Dynamics and Control Using Quaternions This chapter is devoted to explicitly obtaining Lagrange's rotational equations of motion in terms of quaternions using the methodology developed in chapter 2. Additionally, the nonlinear control torque required to change the orientation of a rigid body is explicitly obtained using the control methodology developed in chapter 3. As discussed in chapter 2, the appropriate selection of an attitude parametrization in rigid body dynamics is critical since many attitude descriptions are singular for certain orientations. The unit quaternion (also called Euler parameters) provides a compact nonsingular (no geometric singularities) representation of rotation and has a wide range of application in rigid body dynamics due to its ability to deal with large angle motions. The development of Lagrange's equations of rotational motion in terms of quaternions has been carried out in [13] and [14]. However, it doesn't appear that Lagrange's formulation has received much use in conventional rotational dynamic or control analysis. 54 4.1 Rotational Kinetic Energy in Terms of Quaternions Consider a unit quaternion 2R 4 (the superscript i is omitted for brevity) written as a 4 by 1 column vector by = 1 ; T T = cos 2 ; e T sin 2 T : (4.1) By a unit quaternion, it is implied that the square of its norm is unity given by N () := T = 1: (4.2) Typically, 1 is termed its scalar part and = [ 2 ; 3 ; 4 ] T is termed its vector part. Geometrically e = [e 1 ;e 2 ;e 3 ] T is any unit vector dened relative to an inertial coordinate frame and is a rotation about this unit vector. The 3-vector e is known as the principal rotation axis and is the principal rotational angle. The unit quaternion in equation (4.1) thus describes the orientation of a rigid body. The orientation of the rigid body is determined by so that represents an identical orientation. This is shown by taking the principal rotation angle 0 = + 2 so that 0 = cos 0 2 ; e T sin 0 2 T = cos 2 +; e T sin 2 + T = cos 2 ; e T sin 2 T =: (4.3) 55 Thus, the unit quaternions and represent the exact same orientation since the principal rotation angles and +2 taken about the principal rotation axis e yield the same orientation. Consider now a rigid body that has an angular velocity, !2R 3 , with respect to an inertial coordinate frame. As described in chapter 2, the components of this angular velocity with respect to its body-xed coordinate frame whose origin is located at the body's center of mass are denoted ! , ! , and ! . In terms of a unit quaternion, ! is given by ! = 2E _ =2 _ E :=g(; _ ); (4.4) where the transformation matrix E is E = 2 6 6 6 6 6 4 2 1 4 3 3 4 1 2 4 3 2 1 3 7 7 7 7 7 5 : (4.5) To obtain _ (!;), we note that the general solution to the equationT E _ := 2E _ =! is _ =T + E ! + (IT + E T E )w; (4.6) where I is a 4 by 4 identity matrix and w is an arbitrary 4-vector. Since the matrix T + E = 1 2 E T , we have _ = 1 2 E T ! + (IE T E)w: (4.7) 56 When is a unit quaternion, we nd that E T E = (I T ), which yields _ = 1 2 E T ! + T w: (4.8) If is a unit quaternion, we also know that T _ = 0 by dierentiating equation (4.1). Pre- multiplying equation (4.8) by T , we nd the requirement T w = 0 since T E T = 0 13 . Thus, the inverse transformation of equation (4.4) is uniquely given by _ = 1 2 E T !: (4.9) Without loss of generality, we shall now assume the body-xed coordinate axes of the rigid body are aligned along the principal axes of inertia of the rigid body. The principle moments of inertia are given byJ ,J , andJ . The rotational kinetic energy of the rigid body is then given by the second term in equation (2.17) as T = 1 2 ! T J! = 2 _ T E T JE _ = 2 T _ E T J _ E: (4.10) The principal inertia matrix of the rigid body is represented byJ =diag(J ;J ;J ). The second and third equalities in equation (4.10) follow from equation (4.4). 4.2 Lagrange's Rotational Equations of Motion Given the kinetic energy of the rigid body in terms of the generalized coordinates (unit quaternions), we now proceed to obtain Lagrange's rotational equations of motion by 57 following the procedure outlined in chapter 2. Assuming the generalized coordinates are independent of one another, Lagrange's equation becomes d dt @T @ _ @T @ = ; (4.11) which yields four second-order nonlinear dierential equations where the 4-vector rep- resents the generalized torque. Noting that @T @ _ = 4E T JE _ , we obtain d dt @T @ _ = 4E T JE + 4 _ E T JE _ + 4E T J _ E _ : (4.12) Also, since @T @ = @ @ 2 T _ E T J _ E = 4 _ E T J _ E =4 _ E T JE _ ; (4.13) Lagrange's equation, assuming the components of the unit quaternion are all independent, then becomes 4E T JE + 8 _ E T JE _ + 4E T J _ E _ = : (4.14) By using equation (4.5) and noting that _ E _ = 0 31 , the last term on the left hand side of equation (4.14) becomes zero. Equation (4.14) is thus 4E T JE + 8 _ E T JE _ = ; (4.15) where the matrix M = 4E T JE multiplying the 4-vector is symmetric and singular. The equation of motion in (4.15) presumes the components of the 4-vector are independent of one another. However, the four components of are actually constrained by 1 modeling constraint because KK min = 1. Since represents a rotation, it must 58 satisfy equation (4.2); i.e., must be a unit quaternion. The modeling constraint is then given by ' m :=N () 1 = 0: (4.16) Dierentiating equation (4.16) twice with respect to time yields the modeling constraint matrix equation A m := T =N ( _ ) :=b m : (4.17) Since the 4 by 4 matrix M is singular, we nd the augmented system by pre-multiplying A + m to equation (4.17) and adding the result to equation (4.15) to obtain M := 4E T JE + T =8 _ E T JE _ N ( _ ) + := F; (4.18) as illustrated in chapter 2. The explicit equation of motion is then directly given by 4E T JE + T =8 _ E T JE _ N ( _ ) + + M 1 2 B + m (b m A m a); (4.19) where a is the acceleration obtained from equation (4.18) as a = 1 4 E T J 1 E + T 8 _ E T JE _ N ( _ ) + : (4.20) We will now proceed to assemble the various quantities required to determine the last term on the right hand side of equation (4.19). This is accomplished by computing A m a, then 59 b m A m a, and nally M 1 2 B + m . We begin by noting that E = 0 31 and T E T = 0 13 . Hence, we nd that A m a = 1 4 T E T J 1 E + T T 8 _ E T JE _ N ( _ ) + = T 8 _ E T JE _ N ( _ ) + = 2! T J!N ( _ ) + T ; (4.21) where in the second and third equalities above we have used equations (4.2) and (4.4). This yields b m A m a =N ( _ ) 2! T J! +N ( _ ) T =2! T J! T : (4.22) Next, the matrix B + m is found as B + m = T 1 2 E T J 1 2 E + T + = T + = ; (4.23) where T + = since the matrix T is idempotent [10]. Finally, we nd that M 1 2 B + m = 2E T J 1 2 E + T =; (4.24) where equation (4.2) is again used. 60 Using equations (4.22) and (4.24), equation (4.19) reduces to 4E T JE + T =8 _ E T JE _ N ( _ ) 2 ! T J! + I T : (4.25) Equation (4.25) is Lagrange's equation for rotational motion in terms of unit quaternions, which was also found in [13]. In equation (4.25),I is the 4 by 4 identity matrix. Note the appearance of the term 2 ! T J! on the right hand side of equation (4.25). This term arises because of the modeling constraint ' m in equation (4.16) and is simply twice the kinetic energy of rotation of the rigid body, which is a result obtained in [14]. Equation (4.25) also informs us|and this does not seem to be widely recognized|that the gen- eralized torque must appear in the form (I T ). Furthermore, post-multiplying the relation I = (I T ) + T by , we nd that = (I T ) + T : (4.26) Thus, T is the component of the generalized torque 4-vector in the direction of the unit 4-vector , and (I T ) is the projection of in the plane normal to . The generalized acceleration is now explicitly obtained as = 1 4 E T J 1 E + T 8 _ E T JE _ N ( _ ) 2 ! T J! + I T : (4.27) Expanding equation (4.27), we obtain =2E T J 1 E _ E T JE _ N ( _ ) + 1 4 E T J 1 E (4.28) 61 by virtue of equations (4.2) and (4.4). Also, we nd that E _ E T = 1 2 ~ !; (4.29) where ~ ! = 2 6 6 6 6 6 4 0 ! ! ! 0 ! ! ! 0 3 7 7 7 7 7 5 : (4.30) The rst term on the right hand side of equation (4.28) is then expressed as 2E T J 1 E _ E T JE _ = 1 2 E T J 1 ~ !J!: (4.31) Thus, equation (4.28) simplies to = 1 2 E T J 1 ~ !J!N ( _ ) + 1 4 E T J 1 E := a m : (4.32) This is the explicit equation for the generalized acceleration of a rotating rigid body in terms of unit quaternions. In equation (4.32),E is the 3 by 4 transformation matrix given by equation (4.5) and the 4 by 4 matrix E T J 1 E is singular. In the next section, we will utilize this formulation to investigate the explicit closed form control torque required to rotate a rigid body form any initial orientation to any desired orientation. 4.3 Control of Rotational Motion in Terms of Quaternions Consider a rigid body whose orientation at any timet is described by the unit quaternion (t). The body is free to pivot about its center of mass. It is at rest, and its initial 62 orientation is given by the quaternion(0) = initial . In what follows, we aim to determine the explicit control torque c that changes the initial orientation, initial , of the rigid body to a desired orientation described by the quaternion d , and the rigid body is brought to rest at d . The equation of motion for the generalized acceleration of the rigid body is given by ( = 0) = 1 2 E T J 1 ~ !J!N ( _ ) + (t); (0) = initial ; _ (0) = 0; (4.33) where (t) = M 1 P ^ c (see chapter 3) is the acceleration 4-vector that is provided by the control torque so that the body rotates to the desired orientation d . The rigid body achieves the desired orientation when it reaches the control constraint c (t) :=(t) d = 0; (4.34) which describes the control objective. At this point, we could make the problem more general by taking d (t), but because the goal is to achieve the stationary orientation d we'll progress with the assumption that d is any constant 4-vector in the set f d 2R 4 : N ( d ) = 1g: (4.35) Since the rigid body initially does not start on the control constraint in equation (4.34), we can approach this control constraint by using the modied constraint equation given by (3.30). The desired control constraint is then obtained as (, > 0) = _ ( d ); (4.36) 63 where = diag( 1 ; 2 ; 3 ; 4 ) and = diag( 1 ; 2 ; 3 ; 4 ). Note that equation (4.36) describes the desired independent paths for each of the quaternion components. How- ever, throughout the control maneuver the satisfaction of equation (4.2) must also be guaranteed because the quaternion (t) must have a unit norm so that it represents a physical rotation. Thus, only three of the four components of the quaternion (t) may be independently controlled at each instant of time since the fourth component of (t) is dependent on the other three components by virtue of the modeling constraint equation (4.2). In the following sections, we present two control strategies to handle this control problem. The rst strategy simply enforces the unit norm equation (4.2) along with the modied control constraints in equation (4.36) for i (t), i = 2; 3; 4; (i.e., the vector part (t) of the quaternion). The second strategy imposes all four of the modied control constraints (4.36) and thus attempts to satisfy a physically unrealizable trajectory as described in chapter 3. In either strategy, we obtain the appropriately constructed control constraint equation (3.26). The resulting generalized control torque is found as c = M P ^ c : (4.37) The components of this torque in the body-xed coordinate frame are found by rst dierentiating equation (4.4), which yields _ ! = 2 _ E _ u + 2E u = 2E u: (4.38) 64 Pre-multiplying equation (4.32) by 2E, we obtain 2E u =EE T J 1 ~ !J! + 1 2 EE T J 1 E c : (4.39) Here, we nd that the 3 by 3 matrix EE T =I. Thus, we obtain _ ! =J 1 ~ !J! + 1 2 J 1 E c : (4.40) We see that Euler's rotational equations of motion appear in equation (4.40), and we nd that the generalized torque c is related to a torque in the body-xed coordinate frame by B = 1 2 E c ; (4.41) where B = ; ; ; T are the control torques about the body-xed , , and principal directions. 4.4 Rotational Controller 1 As a consequence to the general reorientation control problem described in the previous section, let us now investigate the application of the explicit controller obtained in chapter 3 in nding the control acceleration (t). To design (t), we take the unit norm constraint 65 (equation (4.2)) along with the last 3 equations in (4.36) so that the matrix A c and the vector b c in equation (3.26) become A c = 2 6 6 6 6 6 6 6 6 6 4 1 2 3 4 0 1 0 0 0 0 1 0 0 0 0 1 3 7 7 7 7 7 7 7 7 7 5 and b c = 2 6 6 6 6 6 6 6 6 6 4 N ( _ ) 2 _ 2 2 ( 2 2;d ) 3 _ 3 3 ( 3 3;d ) 4 _ 4 4 ( 4 4;d ) 3 7 7 7 7 7 7 7 7 7 5 : (4.42) Let us assume initially that nowhere along the controlled trajectory does 1 = 0, so that the matrixA c in equation (4.42) is nonsingular. The control acceleration is then explicitly found by (t) = M 1 P ^ c = M 1 2 B + c (b c A c a m ) = 1 2 E T J 1 ~ !J! +N ( _ ) +: (4.43) Notice that P ^ c = ^ c since we are prescribing a control constraint that is exactly along the modeling constraint manifold. The 4-vector = 1 ; T T is simply 1 1 = N ( _ ) + T ^ _ + T ^ ( d ) := 1 1 (4.44) = ^ _ ^ ( d ) :=; (4.45) where the 3 by 3 matrices ^ = diag( 2 ; 3 ; 4 ) and ^ = diag( 2 ; 3 ; 4 ). Observe that equation (4.45) exactly represents the desired control constraints given by the last 3 equa- tions in (4.36). Thus, we have obtained the nonlinear controller in equation (4.43), which 66 causes the coordinates 2 , 3 , and 4 to independently converge to their respective desired values given by the corresponding components of d . The coordinate 1 is subsequently found so that the unit norm constraint (equation (4.2)) is exactly satised at each instant of time as required by the rst control constraint in the set (4.42). Furthermore, we can immediately arrive at stability results for the controlled system in equations (4.44) and (4.45). The 3 equations in (4.45) are asymptotically stable at the vector part of the quaternion d provided ^ > 0 and ^ > 0. We can further deduce that the entire system is asymptotically stable at the two xed points = 1;d ; T d T ; (4.46) where the coordinate 1;d = q 1 2 2;d 2 3;d 2 4;d . We see that two xed points exist| a consequence of the ambiguity that results when solving for the coordinate 1 using equation (4.2). Let us now consider the case when 1 = 0. The matrix A c in equation (4.42) now becomes singular, and as a result the control constraints (3.26) are, in general, not con- sistent for all vectors b c . It is only consistent for the vectors b c that lie in the column space of the matrix A c . In this case, consistency is necessary to ensure that the norm of the quaternion (t) is unity and hence represents a physical rotation. In fact, equation (4.44) points out that when 1 = 0, for nite accelerations 1 , we would require that N ( _ ) + T ^ _ + T ^ ( d ) = 0: (4.47) 67 Note that the equation A c =b c is consistent if and only if [10] A c A + c b c =b c : (4.48) When 1 = 0, the matrix A c A + c is A c A + c = 2 6 4 0 T 0 31 I 3 7 5 2 6 4 0 0 13 1 2 I 1 2 T 3 7 5 = 2 6 4 1 2 T 1 2 T 1 2 I 1 2 T 3 7 5: (4.49) The requirement in equation (4.48) then yields the condition N ( _ ) + T ^ _ + T ^ ( d ) = 0; (4.50) which is the same condition given in equation (4.47). In general, the satisfaction of this condition, which arises when 1 = 0 is not possible for arbitrary values of d , ^ , and ^ . Thus, crossing the hyperplane 1 = 0 at any timet would generally result in a physically unrealizable rotation since the control constraints may not be consistent at 1 = 0. Physically speaking, the hyperplane 1 = 0, which corresponds to principal rotation angles =;2;:::; appears to separate the two xed points found in equation (4.46). The two xed points are separated by the boundary f2R 4 : N () = 1 \ 1 = 0g; (4.51) 68 which is the space of rotations dened by the intersection of the hyperplane 1 = 0 and the unit 3-sphere S 3 . These results show the controller found by equation (4.43) causes rotational trajectories of the rigid body to occur on the unit 3-sphere, and the rigid body will asymptotically approach one of the xed points in equation (4.46). To control the rigid body on the boundary (4.51) the condition in equation (4.47) must be satised. For rest-to-rest maneuvers, the sign of the coordinate 1 at the initial time determines the specic xed point to which the controller will rotate the rigid body. Yet, crossing the hyperplane 1 = 0 to achieve the desired control objective is not necessary. If the quaternions initial and d have their rst components ( 1;intitial and 1;d ) of opposite sign, we merely need to command the controller to orient the rigid body to the diametrically opposite quaternion, d , which corresponds to precisely the same physical orientation as the quaternion d . 4.5 Rotational Controller 2 In many instances, it may be desirable to control the rigid body in the entire set S 3 =f2R 4 : N () = 1g (4.52) such that no restriction on its motion exists. Recall that the prior control strategy was brought about by assigning the control constraints as the unit norm requirement (equation (4.2)) in conjunction with the last 3 equation in (4.36). Consider now all the modied control constraints in equation (4.36) so that A c =I and b c = _ ( d ): (4.53) 69 The 4-vector ^ c that causes each component of to independently approach d is explicitly found as ^ c = M 1 2 B + c (b c A c a m ) = M(b c a m ): (4.54) The controller (t) that adjusts the desired control constraints to occur along the modeling constraint manifold is then given by (t) = M 1 P ^ c = M 1 E T E M(b c a m ) = 1 2 E T J 1 ~ !J! +E T Eb c : (4.55) The generalized acceleration of the controlled system found by equation (4.33) is thus =N ( _ ) +E T Eb c : (4.56) Although equation (4.56) appears simple, the behavior of this nonlinear dynamical system can be quite complex. We shall now concentrate on determining the xed points of equation (4.56) and also their stability. The xed points are given by _ = 0 and any 4-vector = that satises the equation E T E( d ) = 0; (4.57) 70 where we recall that both and d are unit quaternions. One obvious solution to equation (4.57) is the isolated xed point = d . In general, for some real scalar, the satisfaction of equation (4.57) requires that the vector ( d ) = (4.58) since is a unit quaternion and E T E is the orthogonal projection operator that projects any vector in R 4 onto a plane normal to the 4-vector . Re-arranging, equation (4.58) becomes (I) = d : (4.59) The general solution to equation (4.59) is =H + d + (IH + H)w; (4.60) whereH =I andw is an arbitrary 4-vector. Observe that the solutions to equation (4.59) depend upon the values of , , and d . Also, note that the matrix H is not invertible only when = i ; i = 1; 2; 3; 4. By inspection of equation (4.59), we see that this could only occur when d has one or more zero components. When H is invertible equation (4.60) becomes i = i i i;d ; i = 1; 2; 3; 4: (4.61) The xed points are found by rst determining those real values that satisfy 4 X i=1 i i 2 2 i;d = 1; (4.62) 71 and subsequently equation (4.61). A solution of equation (4.62) is always = 0, which corresponds to = d , a result we have already obtained by inspection of equation (4.57). When the diagonal matrix H is singular we can assume that = i for some i2f1; 2; 3; 4g. For consistency of equation (4.58), this requires that i;d = 0. Without loss of generality, let us assume that = 4 and that, initially, i 6= 4 , i = 1; 2; 3. The solution to equation (4.59) then becomes i = i i 4 i;d ; i = 1; 2; 3; 4 =w 4 ; (4.63) where a realw 4 is sought by using equation (4.2). Depending on the values of this could yield 2 additional isolated xed points. If = 4 = 2 and i 6= 4 ,i = 1; 2, then we have the solutions i = i i 4 i;d ; i = 1; 2; 3 =w 3 ; 4 =w 4 ; (4.64) where again the arbitrary values w 3 and w 4 are determined by seeking real solutions to equation (4.2). We would again require, for consistency, that 3;d = 4;d = 0. Depending on the values of this could yield, in addition to the isolated xed points, a circle of non- isolated xed points. Similarly, if = 4 = 3 = 2 , which requires 2;d = 3;d = 4;d , the solution to equation (4.59) is i = 1 1 4 1;d ; 2 =w 2 ; 3 =w 3 ; 4 =w 4 : (4.65) In this case, we could have, in addition to the isolated xed points and the circle of non- isolated xed points, a sphere of non-isolated xed points. Note that there are a minimum 72 number of 2 isolated xed points and a maximum number of 8 isolated xed points since equation (4.62) can have at most 8 real solutions. When the control parameters i are chosen so that they are all identical ( i = 1 ; i = 2; 3; 4), then =2 1 satises equation (4.62). This yields = d as another xed point of the controlled system, which is exactly the same physical orientation as d as we have shown earlier. In fact, when i are all identical the only solutions to equation (4.59) are = d , and so this particular controller has only 2 xed points. Furthermore, when both of the control parameters i and i are identical, i.e., i = 1 ;i = 2; 3; 4; and i = 1 ; i = 2; 3; 4; the controller in (4.55) is of the same form as the ad hoc controller proposed in [30], which we now see is a subset of the controllers developed herein. We now consider the stability of the controlled system in equation (4.56) at the xed point of interest ( _ = 0; = d ). We shall assume, for simplicity, that either all the values of i are identical, or that they are all distinct, so that we always have a set of isolated xed points. Let (t) = d +(t); (4.66) where (t) is a perturbed trajectory of the rotating rigid body. The equation governing the perturbation then becomes = N ( _ )( d +) +E T E( _ ) = N ( _ )( d +) + I ( d +)( d +) T ( _ ) (4.67) since E T E =I T . Consider now the Lyapunov function V L (; _ ) = 1 2 T + 1 2 _ T _ ; (4.68) 73 which is positive for ; _ 6= 0 since > 0. Dierentiating equation (4.68) with respect to time we get _ V L = _ T + _ T = _ T + _ T h N ( _ )( d +) i + _ T I ( d +)( d +) T ( _ ) = _ T _ N ( _ ) _ T ( d +) _ T ( d +)( d +) T ( _ ): (4.69) Since the perturbed trajectory (t) describes a rotational motion, we nd that ( d + ) T ( d +) = 1 and upon dierentiating with respect to time we have _ T ( d +) = 0. Equation (4.69) then reduces to _ V L = _ T _ ; (4.70) which is negative semi-denite. To show = d is asymptotically stable, we consider LaSalle's theorem as follows. In the vicinity of the equilibrium point = 0, let us choose a suitably small region D. This region is positively invariant since V L is positive denite and _ V 0. Moreover, in the region D, when _ V L = 0 we have _ = 0 by equation (4.70). By equation (4.59), the set W =f2D : _ V L = 0; N ( d +) = 1g (4.71) that contains all values of 2D and satises the equation I ( d +)( d +) T = (I T ) = 0 (4.72) 74 has only one invariant trajectory = 0. This is because the solution,, of equation (4.72) is of the general form = 1 (4.73) for some real scalar such thatN ( d +) = 1. The xed point = 0 is always an isolated solution of equation (4.73). The other values of are found from equations (4.63), (4.64), and (4.65), and the corresponding xed points are isolated since the parameters i are assumed distinct. We see the largest (only) invariant set in W is then given by = 0, making the xed point = d asymptotically stable. 4.6 Numerical Examples In this section, numerical examples are provided to illustrate the capability of the two control methodologies presented in the last two sections for a rigid body with principal inertias J = 100 kg m 2 ; J = 200 kg m 2 ; and J = 250 kg m 2 . The simulations are carried out for rest-to-rest maneuvers so that _ initial = 0 and the nal orientation is d . 4.6.1 Rotational Controller 1 First, we employ the rotational controller described in section (4.4). Using this approach, the equations of motion of the controlled system were found by assuming 1 6= 0 in equations (4.44) and (4.45). The control parameters ^ and ^ are chosen as ^ = diag 3 5 ; 9 20 ; 9 25 and ^ = diag 1 9 ; 1 16 ; 1 25 : (4.74) 75 The desired quaternion, d , is chosen as d = 1 p 3 ; 1 p 3 ; 0; 1 p 3 T : (4.75) In order to visualize the controlled rotational trajectories, we shall assume the rigid body rotates in the set 2R 3 : N () = 1; 3 = 0 : (4.76) This allows us to plot rotational trajectories on the unit 2-sphere, S 2 , given by 2 1 + 2 2 + 2 4 = 1. The numerical integration of this system is carried out over a duration of t = 30s using a variable time step Runge-Kutta scheme with a relative error tolerance of 10 12 and an absolute error tolerance of 10 15 . Figure 4.1 shows the the resulting (a) Rotational trajectories converging to the equilibrium points in equation (4.46). The small dark gray sphere denotes the de- sired orientation d and the small light gray cube denotes the undesirable orientation. (b) Rotational trajectories converging to the equilibrium points d and d . The small dark gray sphere denotes the desired orientation d and the small light gray cube denotes the diametrically opposite quater- nion d , which corresponds to the same physical orientation. Figure 4.1: Rotational Controller 1 { Rotational trajectories (dotted lines) on a unit 2-sphere. 76 rotational trajectories on the unit 2-sphere for various initial quaternions initial near the boundary given by equation (4.51). In Figure 4.1(a), rotational trajectories are seen to converge to the xed points in equation (4.46). The intersection of the unit 2-sphere and the hyperplane 1 = 0 is illustrated by the unit circle that is embedded in the plane. When the 4-vector initial is such that 1;initial > 0, the rotational trajectories converge to d , which is denoted by the small dark gray sphere. When 1;initial < 0, the rotational trajectories converge to the undesirable orientation denoted by the small light gray cube. In this situation, it is not necessary to cross the hyperplane 1 = 0 in order to reach our target orientation since we can modify the desired quaternion location to d as shown in Figure 4.1(b). Next, we consider an arbitrary rigid body reorientation from the initial quaternion initial = " 3 10 ; 1 5 ; 7 10 ; p 38 10 # T (4.77) to the desired quaternion d = 9 10 ; 1 5 ; 3 10 ; 3 10 T (4.78) using the control parameters in equation (4.74) for a duration of t = 30s. Figure 4.2 shows the components of the quaternion (t) as a function of time. Figure 4.3 gives the components of the control torques (equation (4.41)) about the body-xed axes required to accomplish the change in orientation. The extent to which the quaternion(t) remains a unit quaternion, as is required, during the control maneuver is shown in Figure 4.4. Noting the vertical scale, we see that ' m and _ ' m are satised to the same numerical order of accuracy as the local error tolerance used in integrating equations (4.44) and (4.45). 77 0 5 10 15 20 25 30 −0.4 −0.2 0 0.2 0.4 0.6 0.8 1 time [s] υ(t) υ 1 υ 2 υ 3 υ 4 Figure 4.2: Rotational Controller 1 { The rotation (components of the 4-vector, (t)) generated by the control acceleration that is explicitly given by equation (4.43) as a function of time. 0 5 10 15 20 25 30 −100 −80 −60 −40 −20 0 20 time [s] Control Torques [N m] ¯ Γ η ¯ Γ μ ¯ Γ γ Figure 4.3: Rotational Controller 1 { Components of the control torque about the body- xed axes given by equation (4.41). 78 0 5 10 15 20 25 30 −3 −2 −1 0 1 2 3 4 5 x 10 −15 time [s] ϕ m and ˙ ϕ m ϕ m ˙ ϕ m Figure 4.4: Rotational Controller 1 {The modeling constraints ' m and _ ' m as a function of time. Note the vertical scale. 4.6.2 Rotational Controller 2 We now apply the second rotational controller given by equation (4.55) to the system in equation (4.33). The application of this controller is highlighted by considering the following two examples. Identical i When the parameters of the diagonal matrix i are chosen so that i = 0 , i = 1; 2; 3; the controlled system has two xed points, namely, the asymptotically stable xed point = d and a second xed point = d as we have found by equations (4.61) and (4.62). In this example, we shall assume the desired quaternion is d = [1; 0; 0; 0] T : (4.79) 79 (a) The stable xed point = d is de- noted by the small dark gray sphere and the unstable xed point = d is denoted by the small light gray cube. (b) The stable xed point = d is de- noted by the small dark gray sphere and the unstable xed point = d is denoted by the small light gray cube. Figure 4.5: Rotational Controller 2 i = 1 8 ; i = 0; 1; 2; 3 { Rotational trajectories (dot- ted lines) on the unit 2-sphere. Here, the xed points of the controlled system are then = [1; 0; 0; 0] T . The control parameters are chosen as = diag p 2 2 ; 2 p 3 ; 2 p 2 p 7 ; 2 p 2 ! and = diag 1 8 ; 1 8 ; 1 8 ; 1 8 : (4.80) As before, we visualize the rotational trajectories of the controlled system on the unit 2-sphere. Figure 4.5 shows rotational trajectories for various initial orientations near the unstable xed points = d (see Appendix B). Figure 4.5(a) reveals that initial orientations starting near = [1; 0; 0; 0] T will travel to the stable xed point = [1; 0; 0; 0] T , which is the same physical orientation separated by a principal rotation angle of = 2. To travel along the shortest rotational path to the desired orientation ( d ), we can always modify the desired orientation to d as shown in Figure 4.5(b). 80 0 10 20 30 40 50 60 70 80 90 −1 −0.5 0 0.5 1 time [s] υ(t) υ 0 υ 1 υ 2 υ 3 Figure 4.6: Rotational Controller 2 i = 1 8 ; i = 0; 1; 2; 3 { The rotation (components of the 4-vector, (t)) generated by the control acceleration that is explicitly given by equation (4.55) as a function of time. The rotational trajectories (t) corresponding to the initial quaternion initial = [0:995; 0:050; 0; 0:0865] T are shown in Figure 4.6. Also, for the same initial quaternion, the control torques about the body-xed coordinate axes required to accomplish the ma- neuver are shown in Figure 4.7. Figure 4.8 shows the extent to which the quaternion remains a unit quaternion. Distinct i When the control parameters i are all distinct it is possible to have a maximum of 8 isolated xed points, but non-isolated xed points are not possible since there are no repeated values of i . In this example, the distinct parameters i are chosen as = diag 1 8 ; 1 3 ; 2 7 ; 1 2 : (4.81) 81 0 10 20 30 40 50 60 70 80 90 −8 −6 −4 −2 0 2 4 6 8 time [s] Control Torques [N m] ¯ Γ η ¯ Γ μ ¯ Γ γ Figure 4.7: Rotational Controller 2 i = 1 8 ; i = 0; 1; 2; 3 { Components of the control torque about the body-xed axes given by equation (4.41). 0 10 20 30 40 50 60 70 80 90 −4 −3 −2 −1 0 1 2 3 4 x 10 −15 time [s] ϕ m and ˙ ϕ m ϕ m ˙ ϕ m Figure 4.8: Rotational Controller 2 i = 1 8 ; i = 0; 1; 2; 3 { The modeling constraints' m and _ ' m as a function of time. Note the vertical scale. 82 The desired quaternion d and the diagonal matrix are again given in equation (4.79) and (4.80), respectively. Given d and we can nd all the xed points of the controlled system as outlined in section (4.5). In this problem, equations (4.61) and (4.62) yield the isolated xed points = [1; 0; 0; 0] T : (4.82) In addition, since d contains zero components and the matrix contains no repeated values, we have the possible scalars = i , i = 1; 2; 3. When = 1 , we nd the 2 isolated xed points = 3 5 ; 4 5 ; 0; 0] T T ; (4.83) and taking = 2 we nd the 2 isolated xed points = " 7 9 ; 0; 4 p 2 9 ; 0 # T : (4.84) Finally, for = 3 , we obtain the 2 isolated xed points = " 1 3 ; 0; 0; 2 p 2 3 # T : (4.85) Thus, this particular reorientation controller yields a total of 8 isolated xed points. All of these xed points are found to be unstable except the desired orientation = d (see Appendix B). Therefore, the rotational controller (4.55) provides asymptotic stability at the two xed points d , which represents the desired physical orientation. We thus see that the desired orientation d appears to be globally stable. 83 (a) Rotational trajectories on S 2 when 2 = 0. The 2 small dark gray spheres denote the xed points d . The 4 small light gray cubes denote the xed points in equations (4.84) and (4.85) (b) Rotational trajectories on S 2 when 3 = 0. The 2 small dark gray spheres denote the xed points d . The 4 small light gray cubes denote the xed points in equations (4.83) and (4.85) (c) Rotational trajectories on S 2 when 4 = 0. The 2 small dark gray spheres denote the xed points d . The 4 small light gray cubes denote the xed points in equations (4.84) and (4.85) Figure 4.9: Rotational Controller 2 0 = 1 8 ; 1 = 1 3 ; 2 = 2 7 ; 3 = 1 2 { Rotational trajec- tories (dotted lines) on the unit 2-sphere. The numerical integration of equation (4.56) is carried out over a duration of t = 30s again using a relative error tolerance 10 12 and an absolute error tolerance 10 15 . Here, we plot rotational trajectories on the unit 2-sphere by assuming the rigid body rotates in the sets given by 2R 3 : N () = 1; 2 = 0 (4.86) and 2R 3 : N () = 1; 4 = 0 ; (4.87) and also in the set (4.76). Figure 4.9 shows rotational trajectories with various initial orientations near the unstable xed points in equations (4.83), (4.84), and (4.85), which are denoted by the light gray cubes (see Appendix B). In Figures 4.9(a), 4.9(b), and 4.9(c) rotational trajectories are seen to converge to the two stable xed points = d denoted by the dark gray spheres. The unstable xed points are saddle nodes (see Appendix B), 84 0 5 10 15 20 25 30 −0.8 −0.6 −0.4 −0.2 0 0.2 0.4 0.6 0.8 1 time [s] υ(t) υ 0 υ 1 υ 2 υ 3 Figure 4.10: Rotational Controller 2 0 = 1 8 ; 1 = 1 3 ; 2 = 2 7 ; 3 = 1 2 { The rotation (components of the 4-vector, (t)) generated by the control acceleration that is explicitly given by equation (4.55) as a function of time. and in each of these gures heteroclinic trajectories that connect the saddle points are shown. The rotational trajectories (t) corresponding to the initial quaternion initial = [0:333;0:200; 0; 0:921] T are shown in Figure 4.10. The resulting control torques in the body-xed coordinate frame are given in Figure 4.11. The extent to which the quaternion (t) remains a unit quaternion throughout the maneuver is again of the order of the error tolerances used in the integration process as shown in Figure 4.12. 4.7 Summary In this chapter, a general approach for the rotational dynamics and control of rigid bod- ies using quaternions was developed. In a direct manner, we have explicitly obtained Lagrange's equations for the rotational motion of a rigid body using quaternions. An 85 0 5 10 15 20 25 30 −15 −10 −5 0 5 10 15 20 25 time [s] Control Torques [N m] ¯ Γ η ¯ Γ μ ¯ Γ γ Figure 4.11: Rotational Controller 2 0 = 1 8 ; 1 = 1 3 ; 2 = 2 7 ; 3 = 1 2 { Components of the control torque about the body-xed axes given by equation (4.41). 0 5 10 15 20 25 30 −8 −6 −4 −2 0 2 4 6 8 x 10 −15 time [s] ϕ m and ˙ ϕ m ϕ m ˙ ϕ m Figure 4.12: Rotational Controller 2 0 = 1 8 ; 1 = 1 3 ; 2 = 2 7 ; 3 = 1 2 { The modeling constraints ' m and _ ' m as a function of time. Note the vertical scale. 86 explicit equation giving the generalized accelerations of the rotating rigid body was also obtained. By applying the idea of control constraints to the modeled system, we have also developed two simple control strategies for changing the orientation of a rigid body to a desired orientation. The rst control strategy allows the stabilization to the desired orientation to occur exactly along a pre-selected trajectory that is desired for the vector part of the quaternion, , while ensuring that remains a unit quaternion. The second control strategy independently prescribes trajectories for each of the quaternion compo- nents using the arbitrary K-vector ^ c , which species an inconsistent and non-physical trajectory. The dynamics are then modied by the explicit control force c = P ^ c so that the unit norm constraint is exactly satised. The explicit nature of the equations obtained for the controlled system allow a detailed analysis of the resulting nonlinear dynamical system that describes the control action. For rest-to-rest maneuvers, the rst control strategy appears to provide global stability to the desired orientation. Though it does not pose any restriction to the control, it was shown there exists a hyperplane in the quaternion 4-space in which the rigid body cannot generally move. The nonlinear dynamical system resulting from the second control strategy shows that the controlled system could have a multiplicity of xed points whose location and stability depend on the control parameters i . Numerical examples were provided for both strategies to illustrate the accuracy with which the closed form control achieves the desired control objective of reorienting the rigid body to a desired orientation. While the second control strategy may generate a multiplicity of xed points, it appears 87 the desired orientation d is globally attracting since xed points corresponding to orientations other than the desired orientation are unstable. 88 Chapter 5 Spacecraft Precision Pointing As a direct application of the formulation developed thus far, let us consider the spacecraft precision pointing problem. Spacecraft pointing is generally a critical design goal in space vehicle design since accurate pointing is necessary for communication, and also to carry out mission objectives. The criteria for pointing may vary depending on the spacecraft system. Space systems that will require high precision pointing include satellite constellations with cross-link requirements and space borne laser systems. In this chapter, we develop in a unied manner the complete model and nonlinear controller necessary for a system of two spacecraft, which may exist in arbitrary orbits, to precisely point at one another. 5.1 Precision Pointing Description Consider two spacecraft (N = 2) in orbit around a rotating central body with a non- uniform gravity eld. It is assumed the two spacecraft move independently of one another in arbitrarily given orbits. Furthermore, each spacecraft is equipped with a payload at an arbitrary point P i , i = 1; 2; xed to the i th spacecraft as shown in Figure 5.1. The two payloads on each spacecraft are required to constantly and precisely track one another 89 1 η 2 η 1 μ 1 γ 2 γ 2 μ 2 ρ 1 ρ 2 P 1 P 2 P r 1 P r 1 V 2 V 1 1 , J m 2 2 , J m δ 1 ˆ n 2 ˆ n Ω x y b z z = b x b y Figure 5.1: Schematic of a two spacecraft system denoted by rigid bodies V 1 and V 2 in orbit around a rotating central body with a non-uniform gravity eld. A payload is xed to each spacecraft located by the position vectors r i P , i = 1; 2: The desired pointing axis ^ n i for the i th spacecraft is required to point along the relative distance vector . by continuously pointing the desired body axes ^ n i , i = 1; 2; at one another, where ^ n i is a unit vector dened relative to the body-xed coordinate frame of the i th spacecraft. The pointing occurs along the line P 1 P 2 , which is denoted by the relative distance vector . This problem is highly nonlinear in both its dynamics and controls aspects. It exhibits motion at three dierent time scales, namely, the time scales associated with the rotation of the central body, the orbital motion of the two spacecraft, and the attitude motion of the two spacecraft due to the payload pointing requirements. The complexity of the 90 problem grows because the two spacecraft require rotational and orbital coordinates to completely describe the ensuing coupled dynamics. In this problem, we shall assume the central body rotates about a xed axis with a constant angular rate . The problem is non-autonomous since the gravitational eld is non-uniform. Thus, we have a multibody spacecraft system whose dynamics are described by a highly nonlinear set of non-autonomous dierential equations, which are coupled in rotational and orbital motions. In the following, we will apply the modeling and control methodologies developed in this dissertation to obtain the explicit nonlinear controller that precisely accomplishes the pointing requirement. 5.2 Model Development First, let us assume the generalized coordinate vector q2R 7 for each spacecraft is given by q i = x i y i z i i 1 i 2 i 3 i 4 T ; i = 1; 2; (5.1) where the 3-vector i = [x i ;y i ;z i ] T represents the Cartesian coordinates relative to the inertial coordinate frame (x;y;z) and the 4-vector i = [ i 1 ; i 2 ; i 3 ; i 4 ] T represents the unit quaternion for thei th spacecraft. It is assumed the two spacecraft move in a non-uniform gravitational eld around the rotating planet Mars. The potential function of the i th spacecraft due to an oblate planet is given by (see Ref. [25]) U i g (r i ; i ; i ) = g m i r i 1U i 1 +U i 2 : (5.2) 91 The rst term in equation (5.2) is simply the gravitational potential corresponding to a point mass. The term U i 1 represents the additional gravitational eects due to zonal harmonics given by U i 1 = 1 X k=2 J k R e r i k P k (sin i ); (5.3) and U i 2 is the contribution due to tesseral and sectoral harmonics given by U i 2 = 1 X k=2 k X l=1 R e r i k P l k (sin i )(C k;l cosl i +S k;l sinl i ): (5.4) The various quantities in the potential function U i g are described in the following list. r i ; i ; i = radius, latitude, and longitude of thei th spacecraft relative to the body- xed coordinate frame of the central body. R e = mean equatorial radius of the central body. P l k (sin i ) = associated Legendre polynomial in sin i . J k ;C k;l ;S k;l = zonal, tesseral, and sectoral harmonic coecients. The rst degree harmonics (k = 1) in the series of equations (5.3) and (5.4) are zero by the assumption that the origin of the coordinate frame (^ x i ; ^ y i ; ^ z i ) coincides with the center of mass of the central body as shown in Figure 5.2. The ^ x i -axis is directed along the radius to thei th spacecraft while the axes ^ y i and ^ z i complete the rectangular coordinate system. The body-xed coordinate frame denoted by (x b ;y b ;z b ) attached to the planet Mars is also shown in Figure 5.2 relative to the inertial coordinate frame. Here, it is assumed that Mars rotates at a constant angular rate about the z b -axis, which is mutually aligned with the inertial z-axis. 92 i V Ω x y b z z = b x b y i x ˆ i y ˆ i z ˆ i λ i ϕ Figure 5.2: The body-xed coordinate axes (x b ;y b ;z b ) xed to the central body and the coordinate frame (^ x i ; ^ y i ; ^ z i ) relative to the inertial coordinate frame (x;y;z). For this problem, let us assume a non-uniform gravity eld up to the fourth harmonic (k = 4). The coecients accurate to this degree for the planet Mars are provided in Appendix C along with other physical properties of Mars including its angular rate . The gravitational force in the (^ x i ; ^ y i ; ^ z i ) coordinate frame on thei th spacecraft due to the non-uniform gravity potential is thus given by ^ i (r i ; i ; i ) = rU i g ; = @U i g @r i ^ r 1 r i @U i g @ i ^ 1 r i sin i @U i g @ i ^ ; (5.5) wherer() is the gradient operator in spherical polar coordinates. The transformation from (^ x i ; ^ y i ; ^ z i ) coordinates to body-xed coordinates (x i b ;y i b ;z i b ) for the i th spacecraft is given by i b =T i 1 ^ i ; (5.6) 93 where T i 1 = 2 6 6 6 6 6 4 cos i cos i sin i sin i cos i cos i sin i cos i sin i sin i sin i 0 cos i 3 7 7 7 7 7 5 : (5.7) Since we have assumed the generalized coordinates i are represented by Cartesian coor- dinates, we convert the spherical coordinates to Cartesian coordinates so that ^ i (r i ; i ; i ) 7! ^ i (x i b ;y i b ;z i b ) (5.8) T i 1 ( i ; i ) 7! T i 1 (x i b ;y i b ;z i b ); (5.9) and we obtain the gravitational force in the coordinate frame (x b ;y b ;z b ) as i b =T i 1 ^ i : (5.10) By assuming that the body-xed coordinate frame (x b ;y b ;z b ) at initial time t = 0 is aligned with the inertial coordinate frame, we simply obtain the transformation from body-xed coordinates i b to inertial coordinates i as a function of time so that i =T 2 i b ; (5.11) where T 2 = 2 6 6 6 6 6 4 cos t sin t 0 sin t cos t 0 0 0 1 3 7 7 7 7 7 5 : (5.12) 94 The gravitational force i g resolved in the inertial coordinate frame at each instant of time is thus i g =T 2 i b : (5.13) Given the above gravity model, we can now proceed to obtain the equations of mo- tion of the two spacecraft system. The transformations f i and g i from the generalized coordinates i and i are given by f i = I i (5.14) g i = 2E i _ i ; (5.15) where the matrix I is a 3 by 3 identity matrix and the 3 by 4 matrix E is given in equation (4.5). Applying equations (2.31) and (2.32), we obtain the equations of motion of the system as M q := 2 6 6 6 6 6 6 6 6 6 4 M 1 0 0 0 0 M 1 0 0 0 0 M 2 0 0 0 0 M 2 3 7 7 7 7 7 7 7 7 7 5 2 6 6 6 6 6 6 6 6 6 4 1 1 2 2 3 7 7 7 7 7 7 7 7 7 5 = 2 6 6 6 6 6 6 6 6 6 4 1 g 8( _ E 1 ) T J 1 E 1 _ 1 2 g 8( _ E 2 ) T J 2 E 2 _ 2 3 7 7 7 7 7 7 7 7 7 5 :=F; (5.16) where the 3 by 3 matrix M i = diag(m i ;m i ;m i ) and the 4 by 4 singular matrix M i = 4(E i ) T J i E i . The generalized forces i = i g (equation (5.13)) and the generalized torques i = 0. The terms8( _ E i ) T J i E i _ i in the 14-vector F are a result obtained by the application of Lagrange's equation as we have seen in chapter 4. 95 The total number of coordinates used in equation (5.16) isK = 14 while the minimum number of coordinates isK min = 12, which reveals we must have two modeling constraints associated with the chosen generalized coordinates. Since the rotational coordinates of the two spacecraft are described by unit quaternions, the required two modeling constraints are given by ' 1;m = N ( 1 ) 1 = 0 (5.17) ' 2;m = N ( 2 ) 1 = 0: (5.18) This yields the modeling constraint matrix equation A m q := 2 6 4 0 13 ( 1 ) T 0 13 0 14 0 13 0 14 0 13 ( 2 ) T 3 7 5 2 6 6 6 6 6 6 6 6 6 4 1 1 2 2 3 7 7 7 7 7 7 7 7 7 5 = 2 6 4 N ( _ 1 ) N ( _ 2 ) 3 7 5 :=b m ; (5.19) which must be imposed to equation (5.16) so that the quaternions 1 and 2 exactly maintain the unit norm requirement (see chapter 4). However, in equation (5.16), the 14 by 14 matrix M is singular because the matrices M 1 and M 2 are singular. We then require that the equations of motion in (5.16) take the form M q :=M +A + m A m =F +A + m b m := F; (5.20) where the acceleration a = M 1 F: (5.21) 96 The explicit equations of motion of the modeled two spacecraft system then become q = a + M 1 2 B + m (b m A m a) := a m ; (5.22) which describes the general motion of each spacecraft moving in a time-varying non- uniform gravity eld. 5.3 Multibody Spacecraft Pointing Control Having obtained the explicit equations of motion for the modeled two spacecraft system, we now dene the pointing requirement, or control objective. The control objective for this problem is dened such that the two spacecraftV 1 andV 2 continuously and precisely point their desired body-xed unit vectors ^ n i ,i = 1; 2; in the direction of one another along the relative distance vector. The vector, which gives the relative distance between the two payloads, is given by = 2 +r 2 P 1 r 1 P : (5.23) Thus, the control objective is achieved when the desired unit vectors ^ n i , i = 1; 2; are aligned along so that the two payloads exactly point at one another, which can be expressed as ^ n 1 kk = ^ n 2 + kk = 0: (5.24) 97 The control constraints that achieve the control objective are then given by k;c := T 1 ^ n 1 2 +T 2 r 2 P 1 T 1 r 1 P k 2 +T 2 r 2 P 1 T 1 r 1 P k = 0; k = 1; 2; 3; (5.25) k;c := T 2 ^ n 2 + 2 +T 2 r 2 P 1 T 1 r 1 P k 2 +T 2 r 2 P 1 T 1 r 1 P k = 0; k = 4; 5; 6; (5.26) which yields a total of 6 control constraints. In equations (5.25) and (5.26), the trans- formation matrixT i for the i th spacecraft from its body-xed coordinate frame to the inertial coordinate frame is given by T i = 2 6 6 6 6 6 4 i 2 1 + i 2 2 i 2 3 i 2 4 2( i 2 i 3 i 1 i 4 ) 2( i 2 i 4 + i 1 i 3 ) 2( i 2 i 3 + i 1 i 4 ) i 2 1 i 2 2 + i 2 3 i 2 4 2( i 3 i 4 i 1 i 2 ) 2( i 2 i 4 i 1 i 3 ) 2( i 3 i 4 + i 1 i 2 ) i 2 1 i 2 2 i 2 3 + i 2 4 3 7 7 7 7 7 5 : (5.27) We then have a set of control constraints by equations (5.25) and (5.26) that exhibit a coupling between orbital and rotational coordinates. However, these control constraints only describe those trajectories that exist when the pointing requirement is satised. When the pointing requirement is not satised, (the unit vectors ^ n i do not point along in the correct direction) we must obtain a trajectory that approaches the control constraint manifold given by the intersection of the 6 control constraints in equations (5.25) and (5.26). This is accomplished by modifying the control constraints to obtain (see chapter 3) k;c + k _ k;c + k k;c = 0; k = 1; 2;:::; 6; (5.28) 98 where k and k are chosen so that the xed point ( _ k;c ; k;c ) = (0; 0) becomes asymp- totically stable. After obtaining the control constraints in the form of equation (5.28), we then have the control constraint matrix equation A c q =b c : (5.29) The explicit control force F c that causes the two spacecraft to converge to and precisely satisfy the pointing requirement is then given by (see chapter 3) F c = P M 1 2 B + c (b c A c a m ); (5.30) and the explicit equations of motion of the controlled two spacecraft system are q = a m + M 1 F c : (5.31) 5.4 Numerical Examples We now demonstrate the ecacy of the formulation developed thus far by considering the following scenario. The two spacecraft, denoted by V 1 and V 2 in Figure 5.1, are in dierent orbits about the rotating planet Mars with the non-uniform gravity model 99 developed in section 5.2. The principal inertia matrix, J i , and mass, m i , for the i th spacecraft are given by J 1 = diag(2300; 4500; 3600) kg m 2 (5.32) m 1 = 2200 kg (5.33) J 2 = diag(1700; 2000; 600) kg m 2 (5.34) m 2 = 1200 kg: (5.35) The payload xed to spacecraft V 1 is located at r 1 P = [2; 1; 0] T m, and the payload xed to spacecraft V 2 is located at r 1 P = [0;0:5; 1] T m relative to their body-xed coordinate frames. The desired pointing axes ^ n i , i = 1; 2; as shown in Figure 5.1, are chosen as ^ n 1 = [1; 0; 0] T and ^ n 2 = [0; 0; 1] T . The spacecraftV 1 is in an elliptic equatorial orbit with an eccentricity e = 0:5, and the spacecraft V 2 is in a circular polar orbit. The initial position and velocity vectors of the two spacecraft are then chosen as 1 (0) = [R e + 150; 0; 0] T km (5.36) _ 1 (0) = [0; 4:256; 0] T km s 1 (5.37) 2 (0) = [0;R e + 300; 0] T km (5.38) _ 2 (0) = 0; 0; r g y 2 (0) T km s 1 ; (5.39) whereR e and g are given in Appendix C. The orbital trajectories for the two spacecraft are shown in Figure 5.3 relative to the rotating planet Mars. 100 Figure 5.3: Orbital trajectories of spacecraft V 1 and V 2 . Spacecraft V 1 is in an elliptic equatorial orbit with an eccentricitye = 0:5 and spacecraftV 2 is in a polar circular orbit. To highlight the capability of the method, we select the initial orientation of the two spacecraft so that their desired pointing axes ^ n i are not correctly aligned along the relative distance vector . We shall also assume the two spacecraft have an arbitrary rotational velocity so that the initial quaternion and quaternion rate vectors are 1 (0) = [1; 0; 0; 0] T (5.40) _ 1 (0) = 1 2 (E 1 ) T ! 1 (0) (5.41) 2 (0) = [1; 0; 0; 0] T (5.42) _ 2 (0) = 1 2 (E 2 ) T ! 2 (0); (5.43) 101 0 25 50 75 100 125 150 −1 −0.5 0 0.5 1 1.5 2 time [s] Φ c Φ1,c Φ2,c Φ3,c (a) Error in the control constraints associ- ated with spacecraft V 1 given by equation (5.25). 0 25 50 75 100 125 150 −1 −0.5 0 0.5 1 1.5 time [s] Φ c Φ4,c Φ5,c Φ6,c (b) Error in the control constraints associ- ated with spacecraft V 2 given by equation (5.26) Figure 5.4: Pointing error over the time period t2 (0; 150)s. where ! 1 (0) = [0:1;0:2; 0:15] T rad s 1 and ! 2 (0) = [0:25; 0:1;0:05] T rad s 1 . The control parameters in equation (5.28) are chosen as i = 0:005 i = 1; 2;:::; 6; (5.44) i = 2 p i i = 1; 2;:::; 6: (5.45) The numerical integration of the controlled two spacecraft system in equation (5.31) is carried out for a time periodt2 (0; 38600)s using a variable time step Runge-Kutta scheme with a relative error tolerance of 10 12 and an absolute error tolerance of 10 15 . This time duration allows spacecraft V 1 to complete approximately 2 orbits while spacecraft V 2 completes approximately 5.658 orbits. First, we examine the systems response over the time period t2 (0; 150)s. We see in Figures 5.4(a,b) that the control constraints given by equations (5.25) and (5.26) approach zero. This indicates that the payload vectors ^ n 1 and ^ n 2 , which are dened relative to 102 0 25 50 75 100 125 150 −50 0 50 100 150 200 time [s] Control Torques [N m] ¯ Γ 1 η ¯ Γ 1 μ ¯ Γ 1 γ (a) Control torques 1 c about the 1 , 1 , and 1 principal axes. 0 25 50 75 100 125 150 −60 −50 −40 −30 −20 −10 0 10 time [s] Control Torques [N m] ¯ Γ 2 η ¯ Γ 2 μ ¯ Γ 2 γ (b) Control torques 2 c about the 2 , 2 , and 2 principal axes. Figure 5.5: Control torques, i c , generated by equation (5.30) required to approach the payload pointing requirement over the time interval t2 (500; 38600). their respective body-xed coordinate frames, are aligning along the relative distance vector , and the two spacecraft system is successfully meeting the control objective. The control torques i c = [ i ; i ; i ] T ,i = 1; 2; (equation (4.41)) about the body-xed coordinate frame required to approach the payload pointing requirements are shown in Figures 5.5(a,b). In addition to the control torques, there are also control forces i c = [ i x ; i y ; i z ] T ,i = 1; 2; required to reach the pointing requirement as shown in Figures 5.6(a,b). Next, we examine the systems response over the time period t2 (500; 38600)s. The extent to which the control constraints are satised over this time interval is shown by Figures 5.7(a,b), and we see they are satised with accuracy commensurate with the tolerances used in the integration. In Figures 5.8(a,b), we nd that the unit norm requirement on the quaternions i , i = 1; 2; is also satised to the accuracy used in integrating the controlled system (5.31). Finally, the control torques required to maintain the pointing requirement over the remaining time interval are shown in Figures 5.9(a,b). 103 0 25 50 75 100 125 150 −4 −3 −2 −1 0 1 2 x 10 −3 time [s] Control Forces [N] ¯ Λ 1 x ¯ Λ 1 y ¯ Λ 1 z (a) Control forces 1 c in the x 1 , y 1 , and z 1 inertial directions. 0 25 50 75 100 125 150 −2 −1 0 1 2 3 4 x 10 −3 time [s] Control Forces [N] ¯ Λ 2 x ¯ Λ 2 y ¯ Λ 2 z (b) Control forces 2 c in the x 2 , y 2 , and z 2 inertial directions. Figure 5.6: Control forces, i c , generated by equation (5.30) required to approach the payload pointing requirement. 0.5 1 1.5 2 2.5 3 3.5 x 10 4 −5 −2.5 0 2.5 5 x 10 −14 time [s] Φ c Φ1,c Φ2,c Φ3,c (a) Error in the control constraints associ- ated with spacecraft V 1 given by equation (5.31) for i = 1; 2; 3. 0.5 1 1.5 2 2.5 3 3.5 x 10 4 −5 −2.5 0 2.5 5 x 10 −14 time [s] Φ c Φ4,c Φ5,c Φ6,c (b) Error in the control constraints associ- ated with spacecraft V 2 given by equation (5.31) for i = 4; 5; 6. Figure 5.7: Pointing error over the time period t2 (500; 38600)s. 104 0 1 2 3 4 x 10 4 −5 −2.5 0 2.5 5 x 10 −15 time [s] ϕ 1,m and ˙ ϕ 1,m ϕ 1,m ˙ ϕ 1,m (a) Error in the modeling constraint associ- ated with spacecraft V 1 given by equation (5.17). 0 1 2 3 4 x 10 4 −5 −2.5 0 2.5 5 x 10 −15 time [s] ϕ 2,m and ˙ ϕ 2,m ϕ 2,m ˙ ϕ 2,m (b) Error in the modeling constraint associ- ated with spacecraft V 2 given by equation (5.18). Figure 5.8: Error in T 1 and T _ over the time interval t2 (0; 38600)s. 1 2 3 x 10 4 −0.8 −0.6 −0.4 −0.2 0 0.2 0.4 0.6 time [s] Control Torques [N m] ¯ Γ 1 η ¯ Γ 1 μ ¯ Γ 1 γ (a) Control torques 1 c about the 1 , 1 , and 1 principal axes. 1 2 3 x 10 4 −0.4 −0.3 −0.2 −0.1 0 0.1 0.2 0.3 time [s] Control Torques [N m] ¯ Γ 2 η ¯ Γ 2 μ ¯ Γ 2 γ (b) Control torques 2 c about the 2 , 2 , and 2 principal axes. Figure 5.9: Control torques, i c , generated by equation (5.30) required to maintain the payload pointing requirement over the time interval t2 (500; 38600)s. 105 5.5 Summary In this chapter, we have developed in a unied manner the complete nonlinear model and controller to accomplish the precision pointing control for a system of two spacecraft orbiting in a non-uniform gravity eld about the rotating planet Mars. The two spacecraft may exist in arbitrarily given orbits and they may also have an arbitrarily given initial rotational state. The pointing is precise in that two payloads xed to each spacecraft at arbitrary locations will point exactly along their relative distance vector for the modeled dynamics. The problem is complex in that there are motions at three dierent times scales, which include the rotational motion of the central body, the orbital motion of the two spacecraft, and the rotational motion of the two spacecraft due to the pointing requirement. The development shows the capability to achieve precise control authority for the general motion of multibody spacecraft problems, where no approximations are considered. 106 Chapter 6 Multibody Spacecraft Tumbling Control In this chapter, we develop the complete nonlinear model and controller for a multibody spacecraft system in which orbital, rotational, and elastic motions are coupled. The tradi- tional method of analyzing spacecraft dynamics is to decouple the attitude motion of the spacecraft about its center of mass and that of its orbital motion. In general, this approach works well, but for complex space systems, such as proposed systems with massive ap- pendages and booms, or for multibody spacecraft systems that are tethered, understand- ing the coupled dynamics and control becomes important, especially for precision-driven space missions. Perturbative forces on these complex systems cause them to torque, ex, compress, expand, and behave in ways that aect the entire dynamics of the system. This requires one to study the motion of the system as a whole if high control authority is required. The multibody spacecraft system considered in this chapter is modeled as two rigid bodies that are elastically connected by a nonlinear spring as illustrated by the example provided in chapter 2. After obtaining the equations of motion that describe the general motion of the elastically connected multibody spacecraft system, the nonlinear controller that causes the system to precisely tumble is developed. 107 6.1 Multibody Spacecraft Tumbling Description Consider the multibody spacecraft system described in section 2.7. The system was mod- eled using two rigid bodies (with masses m 1 andm 2 and principal inertia tensors J 1 and J 2 ) in a dumbbell-like structure. The two spacecraft are connected at arbitrary locations r 1 P and r 2 P by a nonlinear elastic spring with stiness coecients k l and k nl as shown in Figure 6.1, where r 1 P and r 2 P are xed relative to their rotating body-xed coordinate δ Ω x y b z z = b x b y 2 ρ 1 ρ 1 η 2 η 1 μ 1 γ 2 γ 2 μ 2 P 1 P 2 P r 1 P r nl l k k , 1 V 2 V 1 1 , J m 2 2 , J m Figure 6.1: Schematic of the multibody spacecraft system in orbit relative to the rotating central body. frames. The distance between the two spacecraft when the spring is unstretched isL e , and denotes the relative distance between the points P 1 and P 2 given by equation (2.54). The system is free to translate, rotate, and vibrate in space. Yet, relative motion between the spacecraft V 1 and V 2 occurs when the spring stretches or compresses along the line P 1 P 2 , and also in rotation about the line P 1 P 2 . In this problem, we shall assume the 108 system is in a circular polar orbit around the rotating planet Mars whose gravitational eld is non-uniform. As we have seen in chapter 5, we then have a complex multibody spacecraft system whose modeled dynamics are described by equations of motion that are highly nonlinear, and also non-autonomous since the central body rotates as a function of time. Multiple time scales exist due to the rotational motion of the central body, the orbital motion of the two spacecraft, the rotational motion of the two spacecraft, and the elastic motion relative to the two spacecraft. The description of this particular multibody spacecraft system, which has 8 degrees of freedom, is provided in section 2.7. The tumbling trajectory requirement on the multibody spacecraft system is chosen so that the tumbling is a function of its orbital position. The tumbling requirement imposed here is (t) =& 1 (t)z 1 (t) +& 2 (t)z 2 (t); i = 1; 2; (6.1) where & 1 (t) = 2wn kcmk m 1 m 1 +m 2 , & 2 (t) = 2wn kcmk m 2 m 1 +m 2 , and cm = [x cm ;y cm ;z cm ] T de- notes the distance to the center of mass of the system from the origin of the inertial coordinate frame. The angle denotes the principal rotation angle of spacecraft V 1 and V 2 . This particular tumbling trajectory is illustrated in Figure 6.2, where the system is initially orbiting in the x-z orbital plane. The goal is to control the system so that its center of mass remains in a circular orbit in the x-z orbital plane while executing the tumbling motion described by equation (6.1). By assuming y cm = 0, equation (6.1) can be written as = 2w n z cm k cm k = 2w n z cm (x 2 cm +z 2 cm ) 1 2 : (6.2) 109 A B C D cm ρ x z 0 2 > = Θ cm cm n x z w & & π 0 2 < = Θ cm cm n x z w & & π 0 = Θ & 0 = Θ & 0 = Θ 0 = Θ n w π 2 = Θ n w π 2 - = Θ Figure 6.2: Schematic of the orbital and tumbling trajectory of the multibody spacecraft system when y cm = 0. The dark gray circle and the light gray circle denote spacecraft V 1 and V 2 , respectively. The inertial coordinate system is oriented so that the y-axis is directed into the orbital plane. The spacecraft system is orbiting in the counter-clockwise direction and the principal rotational states and _ are listed at intervals A, B, C, and D, which are located at 90 degree increments along the orbit. Here, we presume the angle denotes the rotation about the y-axis so that the principal rotation axis e (see equation (4.1)) is aligned with the inertial y-axis. Beginning at position B in Figure 6.2, we assume the system is orbiting in the counter- clockwise direction. By inspection of equation (6.1), we see that = 2w n when the system is located atx cm = 0 along the positivez-axis. The integerw n denotes the number of desired tumbles incurred from position A to position B. Dierentiating equation (6.2) with respect to time we obtain d dt = @ @x cm dx cm dt + @ @z cm dz cm dt = 2w n " _ z cm (x 2 cm +z 2 cm ) 1 2 x cm z cm _ x cm (x 2 cm +z 2 cm ) 3 2 z 2 cm _ z cm (x 2 cm +z 2 cm ) 3 2 # : (6.3) 110 At position B, equation (6.3) informs us that _ = 0 since x cm = 0. Moving to position C the system is at z cm = 0 and thus = 0 and _ = 2w n _ z cm =jx cm j, where _ < 0 at position C since _ z cm < 0. Therefore, as the system orbits from position B to position C it completes w n tumbles from = 2w n to = 0 wherein its principal angular velocity reaches a maximum _ = 2w n _ z cm =jx cm j at the equatorial plane. As the system continues through its orbit it reaches =2w n and _ = 0 at position D since it is located at x cm = 0 along the negative z-axis. At this point, the system has traveled in the left half plane of Figure 6.2 from position B to position D, where _ 0 since in this region _ z cm 0. To determine the tumbling motion as the system moves in the right half plane of Figure 6.2, we rst note that equation (6.1) possesses symmetry in that the principal rotation angle is symmetric about the z-axis. By examining the tumbling motion in the left half plane, it is then clear the system will rotate from =2w n to = 2w n as it travels from position D to position B. At position A, we have = 0 and again a maximum principal angular velocity of _ = 2w n _ z cm =jx cm j, but in the positive sense since _ z cm > 0. We see the tumbling motion of the system is then determined by equation (6.1) so that it completes 2w n tumbles in each half orbital plane as shown in Figure 6.2. Furthermore, the system tumbles in opposing directions in each plane and repeats this motion every orbital period under the assumption y cm = 0. 6.2 Model Development The general model development of the multibody spacecraft system as shown in Figure 6.1 was discussed in some detail in section 2.7, where the main variation here is the 111 conversion to a non-unform gravity eld for the rotating planet Mars. First, we denote the generalized coordinate vector as usual by q i = x i y i z i i 1 i 2 i 3 i 4 T ; (6.4) which describes the inertial position and the unit quaternion of the i th spacecraft. The transformations from the generalized coordinates i = [x i ;y i ;z i ] T and i = [ i 1 ; i 2 ; i 3 ; i 4 ] T are given by f i = I i (6.5) g i = 2E i _ i ; (6.6) where I is a 3 by 3 identity matrix. The equations of motion of this un-modeled and un-controlled system orbiting in a non-uniform gravity eld are thus found by equations (2.31) and (2.32) as M q := 2 6 6 6 6 6 6 6 6 6 4 M 1 0 0 0 0 M 1 0 0 0 0 M 2 0 0 0 0 M 2 3 7 7 7 7 7 7 7 7 7 5 2 6 6 6 6 6 6 6 6 6 4 1 1 2 2 3 7 7 7 7 7 7 7 7 7 5 = 2 6 6 6 6 6 6 6 6 6 4 1 g @Ue @ 1 8( _ E 1 ) T J 1 E 1 _ 1 @Ue @ 1 2 g @Ue @ 2 8( _ E 2 ) T J 2 E 2 _ 2 @Ue @ 2 3 7 7 7 7 7 7 7 7 7 5 :=F: (6.7) The 14 by 14 matrix M is the singular matrix given in equation (5.16). The generalized forces i consist of the gravitational forces i g (equation (5.13)) and the elastic forces due to the potentialU e given in equation (2.57). The generalized torques i are also generated by the potential U e , and the additional terms8( _ E i ) T J i E i _ i are a result of Lagrange's 112 equation. The system has a total of K = 14 coordinates and a minimum of K min = 12 coordinates. The requiredKK min modeling constraints associated with the generalized coordinates are simply the unit norm requirements on the quaternions i ,i = 1; 2; so that ' 1;m := N ( 1 ) 1 = 0; (6.8) ' 2;m := N ( 2 ) 1 = 0: (6.9) In addition to ' 1;m and ' 2;m , this particular multibody spacecraft system requires the modeling constraints m due to the physical connection of the two spacecraft. These additional modeling constraints were introduced in section 2.7 and are given here as k;m := ^ n 1 = 0; k = 1; 2; 3; (6.10) k;m := ^ n 2 = 0; k = 4; 5; 6: (6.11) Equations (6.10) and (6.11) ensure the orientation of the relative distance vector is xed relative to both of the rotating body-xed coordinate frames of V 1 and V 2 . The unit 3-vectors ^ n 1 and ^ n 2 are determined by the locations of the spring connections to the spacecraft. The relative distance vector was also described in section 2.7 and is given here as = 1 +r 1 P 2 r 2 P : (6.12) Dierentiating equations (6.10) and (6.11) twice with respect to time, we can write all six modeling constraints m as m := ( _ ! i ^ n i ) + 2(! i ^ n i ) _ + ^ n i = 0; i = 1; 2; (6.13) 113 where _ = _ 1 +! 1 r 1 P _ 2 ! 2 r 2 P and = 1 + _ ! 1 r 1 P 2 _ ! 2 r 2 P . Resolving equation (6.13) in the appropriate coordinate frame of the i th spacecraft, we can obtain equations (6.8 - 6.13)|a total ofh = 8 modeling constraints|as our modeling constraint matrix equation A m q =b m ; (6.14) where A m is an 8 by 14 matrix and b m is an 8-vector. Since the matrix M in equation (6.7) is singular, the augmented equations of motion become M q :=M +A + m A m =F +A + m b m := F; (6.15) so that M is symmetric and positive denite and a = M 1 F: (6.16) The explicit equations of motion are then q = a + M 1 2 B + m (b m A m a) := a m ; (6.17) which describes the general motion of the elastically connected multibody spacecraft sys- tem. 6.3 Tumbling Control To achieve the tumbling trajectory illustrated in Figure 6.2, we shall require that the system orbits in a polar orbit around Mars so that the orbital position of the systems 114 center of mass remains circular. The control constraints that achieve this control objective are 1;c := k cm kr cm = 0 (6.18) 2;c := y 1 = 0; (6.19) 3;c := y 2 = 0; (6.20) where r cm denotes the desired radius distance to the center of mass of the multibody system. Additionally, we must impose the necessary control constraints associated with the desired tumbling trajectory given by equation (6.1). These control constraints are obtained directly through the Euler parameters so that k;c := 1 cos =2 e sin =2 T = 0; k = 4; 5; 6; 7; (6.21) k;c := 2 cos =2 e sin =2 T = 0; k = 8; 9; 10; 11; (6.22) where we choose the principal rotation axis as e = [0; 1; 0] T so that the rotation occurs about the body axes i ; i = 1; 2; (see Figure 6.1), which are mutually aligned with the inertial y-axis. Finally, we have the tumbling requirement 12;c := (z 1 z 2 ) cos (x 2 x 1 ) sin = 0; (6.23) which denes the desired tumbling trajectory in terms of the systems position coordinates 1 and 2 . Thus, we have dened a total of 12 control constraints that completely describe the desired multibody spacecraft tumbling in Figure 6.2. 115 We can now dene the modied control constraints as ( k > 0; k > 0) k;c + k _ k;c + k k;c = 0; k = 1; 2;:::; 12; (6.24) so that trajectories, which occur along the modeling constraint manifold given by (6.14), will asymptotically converge to the control constraint manifold A c q =b c ; (6.25) where A c is a 12 by 14 matrix and b c is a 12-vector. The explicit control force F c = [( 1 c ) T ; ( 1 c ) T ; ( 2 c ) T ; ( 2 c ) T ] T that causes the elastically connected multibody spacecraft system to converge to and precisely satisfy the tumbling requirement is then F c = P M 1 2 B + c (b c A c a m ): (6.26) The explicit equations of motion giving the generalized accelerations of the controlled multibody system are q = a m + M 1 F c : (6.27) 6.4 Numerical Examples In this section, we provide a numerical example demonstrating the accuracy with which the nonlinear controller (6.26) stabilizes the multibody spacecraft system along the desired tumbling trajectory. For this example, we shall assume the principal inertia matrix, J i , and mass,m i , of thei th spacecraft are given by equations (5.32)-(5.35). The unstretched 116 length of the nonlinear spring connecting the two spacecraft is L e = 0:002 km while the linear and cubically nonlinear spring coecients are k l = 10 N m 1 and k nl = 1 N m 3 , respectively. The spring connections on spacecraft V 1 and V 2 are located by the body- xed vectors (see Figure 6.1) r 1 P = [0:001; 0; 0:001] T km (6.28) r 2 P = [0:001; 0; 0:001] T km; (6.29) which requires that the unit vectors ^ n 1 = [1; 0; 0] T and ^ n 2 = [1; 0; 0] T . Figure 6.3 illustrates the connection of the two spacecraft and the desired tumbling axis. In this 1 μ 1 γ 2 γ 2 μ 2 P 1 P 2 P r 1 P r nl l k k , 1 V 2 V 1 1 ˆ n = η 2 2 ˆ n = η Θ Θ Figure 6.3: Schematic of spacecraftV 1 andV 2 showing spring connections atP 1 andP 2 . The tumbling axis is depicted along the axes i , i = 1; 2: 117 example, we select w n = 1 so that the system completes one tumble at 1=4 orbit inter- vals as described in Figure (6.2). We assume the multibody spacecraft system is initially located by the inertial position vectors 1 (0) = [R e + 300; 0; 0] T km (6.30) 2 (0) = [R e + 300:004; 0; 0] T km; (6.31) whereR e is the mean equatorial radius of Mars given in Appendix C. The initial orbital velocities of the two spacecraft are given by _ 1 (0) = [0;5:941 10 4 ; 3:404] T km s 1 (6.32) _ 2 (0) = [0;5:941 10 4 ; 3:404] T km s 1 ; (6.33) so that its orbit slightly deviates from a polar circular orbit with an inclinationI = 89:99deg. The initial orientation of the system is given by 1 (0) = [1; 0; 0; 0] T (6.34) 2 (0) = [1; 0; 0; 0] T ; (6.35) and the initial rotational rate is set to zero so that _ 1 (0) = [0; 0; 0; 0] T (6.36) _ 2 (0) = [0; 0; 0; 0] T : (6.37) 118 This yields an initial state (q(0); _ q(0)); which does not exist along the desired tumbling trajectory dened by the 12 control constraints c . The control parameters for all 12 control constraints are chosen for simplicity as k = 0:2 and k = 0:01 fork = 1; 2;:::; 12. The numerical integration of the controlled system (6.27) is carried out over a time duration t2 (0; 13645)s, which correlates to approximately 2 orbits, around the rotating planet Mars using a variable time step Runge-Kutta scheme with a relative error toler- ance 10 12 and an absolute error tolerance 10 15 . The extent to which the 8 modeling constraints' m and m |the quaternion unit norm and the multibody spacecraft physical requirements|are satised is shown in Figure 6.4. All 8 of the modeling constraints are 0 2000 4000 6000 8000 10000 12000 14000 −2 −1.5 −1 −0.5 0 0.5 1 1.5 2 x 10 −12 time [s] ϕ m and Φ m Figure 6.4: Error in the 8 modeling constraints ' m and m over the time duration t2 (0; 13645)s. graphically shown in Figure 6.4 for ease of presentation. We note that, as required, all 8 of the modeling constraints are satised to the local error tolerances used for integrating (6.27). Figures 6.5(a,b) show the error in the 12 control constraints c over the time durations t2 (0; 150)s and t2 (250; 13645)s. Figure 6.5(a) shows the initial errors con- verging to zero, and Figure 6.5(b) shows all 12 control constraints are satised to the 119 0 25 50 75 100 125 150 −12 −10 −8 −6 −4 −2 0 2 x 10 −3 time [s] Φc (a) c over the time duration t2 (0; 150)s. 2000 4000 6000 8000 10000 12000 14000 −4 −3 −2 −1 0 1 2 3 4 x 10 −11 time [s] Φc (b) c over the time duration t 2 (250; 13645)s. Figure 6.5: Error in the 12 control constraints c showing their satisfaction throughout the integration. local error tolerances used throughout the integration. The initial errors exist because we have specied the initial conditions (q(0); _ q(0)) so that the multibody spacecraft system is not traveling along the desired tumbling trajectory. The control forces and torques that cause the multibody system to converge to the tumbling trajectory are provided in Figures 6.6(a,b) and Figures 6.7(a,b). The control forces and torques required to maintain the tumbling trajectory so that the multibody system remains in a polar orbit are shown in Figures 6.8(a,b) and Figures 6.9(a,b). The control inputs are applied to both spacecraft since we have assumed the system is fully actuated. The oscillations sustained relative to the two spacecraft due to the nonlinear spring are shown in Figure 6.10 in millimeters. Finally, the angular velocities ! i ; i = 1; 2; throughout the tumbling trajectory are shown in Figures 6.11(a,b). In Figure 6.11(b), we see the maximum and minimum angular velocities of the multibody spacecraft system correspond to the require- ment _ := 2w n _ z cm =kx cm k =5:786 10 3 , where w n = 1 and the values _ z cm and x cm are easily computed from the initial conditions. 120 0 25 50 75 100 125 150 −20 0 20 40 60 80 100 120 140 time [s] Control Forces [N] ¯ Λ 1 x ¯ Λ 1 y ¯ Λ 1 z (a) 1 c over the time duration t2 (0; 150)s. 0 50 100 150 −20 0 20 40 60 80 time [s] Control Forces [N] ¯ Λ 2 x ¯ Λ 2 y ¯ Λ 2 z (b) 2 c over the time duration t2 (0; 150)s. Figure 6.6: Control forces, i c , required to converge to the desired tumbling trajectory. 0 25 50 75 100 125 150 −1 0 1 2 3 4 5 6 time [s] Control Torques [N m] ¯ Γ 1 η ¯ Γ 1 μ ¯ Γ 1 γ (a) 1 c over the time duration t2 (0; 150)s 0 25 50 75 100 125 150 −0.5 0 0.5 1 1.5 2 2.5 time [s] Control Torques [N m] ¯ Γ 2 η ¯ Γ 2 μ ¯ Γ 2 γ (b) 2 c over the time duration t2 (0; 150)s. Figure 6.7: Control torques, i c , required to converge to the desired tumbling trajectory. 121 2000 4000 6000 8000 10000 12000 14000 −15 −10 −5 0 5 10 15 20 time [s] Control Forces [N] ¯ Λ 1 x ¯ Λ 1 y ¯ Λ 1 z (a) 1 c over the time duration t 2 (250; 13645)s. 2000 4000 6000 8000 10000 12000 14000 −10 −5 0 5 10 time [s] Control Forces [N] ¯ Λ 2 x ¯ Λ 2 y ¯ Λ 2 z (b) 2 c over the time duration t 2 (250; 13645)s. Figure 6.8: Control forces, i c , required to maintain the desired tumbling trajectory over 2 orbits. 2000 4000 6000 8000 10000 12000 14000 −0.04 −0.03 −0.02 −0.01 0 0.01 0.02 0.03 0.04 time [s] Control Torques [N m] ¯ Γ 1 η ¯ Γ 1 μ ¯ Γ 1 γ (a) 1 c over the time duration t 2 (250; 13645)s 2000 4000 6000 8000 10000 12000 14000 −0.015 −0.01 −0.005 0 0.005 0.01 0.015 time [s] Control Torques [N m] ¯ Γ 2 η ¯ Γ 2 μ ¯ Γ 2 γ (b) 2 c over the time duration t 2 (250; 13645)s. Figure 6.9: Control torques, i c , required to maintain the desired tumbling trajectory over 2 orbits. 122 0 2000 4000 6000 8000 10000 12000 14000 −20 −10 0 10 20 30 time [s] δ−L e [mm] Figure 6.10: Vibrational motion between spacecraft V 1 andV 2 throughout the tumbling trajectory. 6.5 Summary In this chapter, we have shown the capability of the modeling and control methodologies developed in this dissertation by applying it to a complex multibody spacecraft system that is required to execute a precise tumbling trajectory while it maintains a specied orbit around the rotating planet Mars. The tumbling trajectory was dened to be a function of the systems orbital position so that the principal rotation angle ( 1 ; 2 ) of spacecraftV 1 andV 2 easily determines the orientation of the multibody spacecraft system based on its orbital location. The coupled vibrational, rotational, and orbital motions are considered in their entirety without any linearizations or approximations, and the generalized control force F c is found in closed form also without any approximation. The complexity of this problem is highlighted by the fact that there are four time scales associated with the dynamics, namely, the spacecraft orbital motions, rotational motions, elastic motions, and the rotational motion of the central body. The modeling and control eort were seen to be see simple in implementation and the numerical results showed the precision with 123 0 25 50 75 100 125 150 −1 0 1 2 3 4 5 6 7 x 10 −3 time [s] ω i [rad s −1 ] ω i η ω i μ ω i γ (a) ! i over the time duration t2 (0; 150)s. 2000 4000 6000 8000 10000 12000 14000 −6 −4 −2 0 2 4 6 x 10 −3 time [s] ω i [rad s −1 ] ω i η ω i μ ω i γ (b) ! i over the time duration t 2 (250; 13645)s. Figure 6.11: Angular velocity,! i , of thei th spacecraft throughout the tumbling trajectory. which the control satised a set of 12 control objectives dened by the control constraints c . 124 Chapter 7 Conclusions This dissertation has provided a new and systematic way of modeling and controlling general multibody systems. The methodology is developed by a judicious application of a fundamental advance in analytical dynamics by Udwadia and Kalaba [27]. We have achieved a synthesis of dynamics and control for multibody systems in that the modeling and control design is carried out using a unied approach. In chapter 2, a general approach to the modeling of N rigid body systems is developed. The development allows the generalized coordinates q i , which describe linear and rotational motions of the i th rigid body, to be independently described by a total of n i coordinates. The multibody system utilizes any desirable number of K generalized coordinates so that K = N X i=1 n i : In addition to this, the generalized coordinatesq i need not represent the minimum number of (independent) coordinates. This particular attribute may cause the matrix M in the equations of motion (2.33) to become singular. Motivated by results obtained in [29], an alternative approach was developed to deal with equations of motion that may have 125 singular K by K matrices M. This is accomplished by obtaining the explicit augmented matrix M, which is symmetric and positive denite. The modeling of the entire multibody system is carried out by dening those consistent modeling constraints so that we obtain an appropriately constructed modeling constraint matrix equation A m q =b m : The general variety of modeling constraints' m (q;t), m (q;t), m (q; _ q;t), and m (q; _ q; q;t) are introduced to enable the modeling of a broad range of modeling requirements. We note that these modeling constraints include all the general holonomic and nonholonomic constraint types. The explicit equations of motion governing the general motion of the modeled multibody system is given by M q = F + F m : = F + M 1 2 B + m (b m +A m a): The forces that arise due to the modeling constraints are explicitly given by F m , and the eect that each modeling constraint has on the multibody system is easily evaluated. To control the multibody system, the permissible set of control forces F c = P ^ F c is derived so that we can eectively generate nonlinear controllers for the general nonlinear system M q = F + F m + F c : The control forces F c are permissible in the sense that they will exactly satisfy the mod- eling constraints of the appropriately modeled system. This points out the importance 126 of distinguishing between modeling and control constraints in multibody systems since any control action would become meaningless if the modeling constraints are violated. The design of the controller is accomplished using the arbitrary K-vector ^ F c , and to sat- isfy a general set of control objectives the control constraints c (q;t) and c (q; _ q;t) are introduced. To ensure stability along these control constraints, they are modied so that A c q =b c yields the asymptotically stable xed points ( c ; _ c ) = (0; 0) and c = 0. Thus, the controller design is reduced to choosing those rst and second-order dierential equations that possess the desired dynamic behavior for reaching the control objective. The explicit control action applied to the modeled multibody system is then F c = P M 1 2 B + c (b c A c a m ): Though the controller F c can be designed to achieve stability at c and c , a proper con- ceptualization of the control problem is necessary to ensure that the modeling constraint manifold and the control constraint manifold can asymptotically converge, or converge with some acceptable error criteria. In describing the general motion of multibody systems it is noted that rotational parameterizations utilizing the minimum number of independent coordinates will yield singularities for certain orientations. The four parameter unit quaternion provides a rotational parametrization that does not suer from singularities. In chapter 4, a direct method for obtaining the generalized accelerations in terms of quaternions is developed, 127 and two new methods are provided for explicitly determining, in closed form, the nonlinear control torque needed to change the orientation of a rigid body. The explicit nature of the equations obtained for the controlled rigid body allows a detailed analysis of the nonlinear dynamical system that describes the control action. Rotational Controller 1 contained a hyperplane with which the rigid body cannot generally move. However, this poses no consequence while attempting to reach the desired orientation, d , since a proper representation of the desired quaternion yields asymptotically stable control from any given initial quaternion to any desired quaternion. The controlled nonlinear dynamical system resulting from Rotational Controller 2 could have a multiplicity of xed points depending on the control parameters i . The stability of the desired quaternion using this controller is shown to be asymptotically stable while unstable xed points arise at locations other than the desired quaternion, and the physical desired orientation appears to be globally attracting. The application of these theoretical developments to realistic problems in multibody spacecraft are developed in chapters 5 and 6. The precision pointing of two spacecraft, which exist in arbitrarily given orbits around a rotating central body with a non-uniform gravity eld, is accomplished so that the two spacecraft will converge to and maintain a precise payload pointing requirement. The development of a model and controller is also demonstrated for a complex multibody spacecraft system required to precisely tumble so that its orientation is determined by its orbital location. Here, the multibody spacecraft system is seen as two rigid bodies elastically connected by a nonlinear spring at arbitrary locations relative to their body-xed coordinate frames. In both applications, the theoretical developments carried out in this dissertation are seen to handle (with considerable ease and accuracy) the modeling and control of complex multibody systems. 128 The entire process utilizes no approximations and provides a new and simple method for modeling and controlling highly nonlinear, multiscale, multibody systems. 129 Bibliography [1] K. M. Abadir and J. R. Magnus. Matrix Algebra. Cambridge University Press, 2005. [2] F. Amirouche. Fundamentals of Multibody Dynamics: Theory and Applications. Birkhauser, 2005. [3] P. Appell. Sur une forme generale des equations de la dynamique. C. R. Acad. Sci., Paris, 129:459{460, 1899. [4] B. Brogliato, R. Lozano, and B. Maschke. Dissipative Systems Analysis and Control. Springer, 2007. [5] W. H. Clohessy and R. S. Wiltshire. Terminal guidance for satellite rendezvous. Journal of Aerospace Science, 27, 1960. [6] P. A. M. Dirac. Lectures in Quantum Mechanics. Yeshiva University, New York, 1964. [7] T. I. Fossen. Guidance and Control of Ocean Vehicles. John Wiley & Sons Ltd., 1994. [8] C. F. Gauss. Uber ein neues allgemeines grundgesatz der mechanik. J. Reine Ange- wandte Mathematik, 4:232{235, 1829. [9] J. W. Gibbs. On the fundamental formulae of dynamics. Amer. J. Math., 2:49{64, 1879. [10] F. Graybill. Matrices with Application to Statistics. Duxbury Publishing Company, 2001. [11] T. R. Kane and D. A. Levinson. Formulation of equations of motion for complex spacecraft. Journal of Guidance and Control, 1980. [12] J. L. Lagrange. M ecanique Analytique. Paris, 1788. [13] H. S. Morton. Hamiltonian and lagrangian formulations of rigid-body rotational dynamics based on the euler parameters. Journal of the Astronautical Sciences, 41:569{591, 1993. 130 [14] P. E. Nikravesh, R. A. Wehage, and O. K. Kwon. Euler parameters in computational kinematics and dynamics. 107:358{369, 1985. Parts 1 and 2. [15] PDS Geosciences Node. Mars global surveyor: Radio science data products. mors1024: MGS85F2 and MGS95I gravity models. [16] J. C. Orden, C. G. Goicolea, and J. Cuadrado. Multibody Dynamics: Computational Methods and Applications. Springer, 2007. [17] L. A. Pars. A Treatise on Analytical Dynamics. Oxbow Press, Woodbridge, Con- necticut, 1979. [18] M. B. Quadrelli, E. Mettler, and J. K. Langmaier. Dynamics and controls of a conceptual jovian moon tour spacecraft. In AIAA/AAS Astrodynamics Specialist Conference and Exhibit, 2004. [19] A. Sanyal, J. Shen, and N. H. McClamroch. Dynamics and control of an elastic dumbbell spacecraft in a central gravitational eld. In Proceedings of 42 nd IEEE Conference on Decision and Control, 2003. [20] A. Sanyal, J. Shen, and N. H. McClamroch. Control of a dumbbell spacecraft using attitude and shape control inputs only. In Proceedings of 2004 American Control Conference, 2004. [21] H. Schaub and J. L. Junkins. Analytical Mechanics of Space Systems. AIAA Educa- tion Series, 2003. [22] H. Schuab, S. R. Vadali, J. L. Junkins, and K. T. Alfriend. Spacecraft formation ying control using mean orbit elements. Journal of the Astronautical Sciences, 79, 2001. [23] A. A. Shabana. Dynamics of Multibody Systems. Cambridge University Press, 1998. [24] M. D. Shuster. A survey of attitude representations. Journal of the Astronautical Sciences, 41(4):439{517, 1993. [25] G.W. Spier. Design and Implementation of Models for the Double Precision Trajec- tory Program (DPTRAJ). Jet Propulsion Laboratory, 1971. Technical Memorandum. [26] F. E. Udwadia. A new perspective on the tracking control of nonlinear structural and mechanical systems. Proc. Roy. Soc. London Ser. A, 459:1783{1800, 2003. [27] F. E. Udwadia and R. E. Kalaba. A new perspective on constrained motion. Proc. Roy. Soc. London Ser. A, 439(1906):407{410, 1992. [28] F. E. Udwadia and R. E. Kalaba. Analytical Dynamics: A New Approach. Cambridge University Press, 1996. 131 [29] F. E. Udwadia and P. Phohomsiri. Explicit equations of motion for constrained mechanical systems with singular mass matrices and application to multi-body dy- namics. Proc. Roy. Soc. London Ser. A, 462:2097{2117, 2006. [30] B. Wie, H. Weiss, and A. Arapostathis. Quaternion feedback regulator for spacecraft eigenaxis rotations. Journal of Guidance, Control, and Dynamics, 12(3):375{380, 1989. 132 Appendix A The Necessary and Sucient Condition for M Being Positive Denite When the K by K matrix M is singular as outlined in chapter 2, we obtain the augmented system given by equation (2.48). We now illustrate the necessary and sucient condition for a positive denite matrix M. Note that M is a symmetric matrix since the matrices M and A + m A m are symmetric. Corollary 5 The symmetric K by K matrix M is positive denite if the matrix M A = 2 6 4 M A m 3 7 5 has full rank. Proof. In the following two parts, we rst show the matrix M has full rank when the matrix M A has full rank. i When the matrix M A does not have full rank, there exists a K-vector w6= 0 such that 2 6 4 M A m 3 7 5w = 2 6 4 Mw A m w 3 7 5 = 0: (A-1) It follows that Mw +A + m A m w = M +A + m A m w = Mw = 0: (A-2) Thus, if the matrix M A does not have full rank, then the matrix M does not have full rank. 133 ii When the matrix M does not have full rank, there exists a K-vector w6= 0 such that M +A + m A m w =Mw +A + m A m w = 0: (A-3) If M =A + m A m , then we have M A w = 2 6 4 Mw A m w 3 7 5 = 2 6 4 A + m A m w A m w 3 7 5: (A-4) Since rank(A m ) = rank(A + m A m ), it follows that M A does not have full rank. If M6=A + m A m , then we must have Mw = 0 (A-5) A + m A m w = 0: (A-6) Pre-multiplying equation (A-6) by A m yields A m A + m A m w =A m w = 0; (A-7) where the third condition of the Moore-Penrose inverse is applied. Therefore, we have 2 6 4 Mw A m w 3 7 5 =M A w = 0: (A-8) Thus, if the matrix M does not have full rank, then the matrix M A does not have full rank. 134 By combining parts (i) and (ii) above, we have proven that the matrix M does not have full rank when the matrix M A does not have full rank. In turn, this yields the necessary and sucient condition that the matrix M has full rank when the matrix M A has full rank. Furthermore, the matrix M is positive denite whenM A has full rank since the matrices M and A + m A m are nonnegative. 135 Appendix B Stability of the Fixed Points for Rotational Controller 2 In chapter 4, we proved that the controlled system, given by the nonlinear equations of motion in (4.56), is asymptotically stable at the xed point = d . To investigate stability at the additional xed points, let us consider equation (4.56) as a system of 8 rst-order dierential equations in phase space so that _ x =f(x); (B-1) where the state vector x = [ 1 ; _ 1 ; 2 ; _ 2 ; 3 ; _ 3 ; 4 ; _ 4 ] T . To facilitate the linearization of this system, we rst reduce the number of rst order dierential equations from 8 to 6, so that the linearized dynamical system is restricted to the modeling constraint manifold N () = 1. This is accomplished by substitution of the modeling constraints T 1 = 0 and T _ = 0 into equation (B-1), so that we have 6 independent rst-order equations of motion. Linearization of the 6 th order dynamical system about the xed point x then yields _ x = @f @x x=x x :=J F x: (B-2) Using the control parameters chosen in equation (4.80), the eigenvalues 1 ; 2 ;:::; 6 of J F are shown in Table 1. Here, the xed point = d is a stable node, which is a result we have already obtained. The xed point = d is an unstable saddle node since its eigenvalues are real and negative. Next, we compute the eigenvalues for the distinct control parameters given by equation (4.81). Table 2 reveals that the two xed 136 points = d are of the stable node type and the remaining 6 isolated xed points are unstable saddle nodes. Fixed Point Eigenvalues = [1; 0; 0; 0] T = [0:095;0:121;0:134;0:935;1:034;1:319] T = [1; 0; 0; 0] T = [0:106; 0:010; 0:083;1:175;1:254;1:498] T Table 1: Eigenvalues of J F for Rotational Controller 2 with identical i . Fixed Point Eigenvalues = [1; 0; 0; 0] T = [0:534;0:534;0:577;0:577;0:707;0:707] T = [1; 0; 0; 0] T = [0:034;0:077;0:207;1:034;1:077;1:207] T = [ 3 5 ; 4 5 ; 0; 0] T = [0:133; 0:043;0:130;1:001;1:112;1:284] T = [ 7 9 ; 0; 4 p 2 9 ; 0] T = [0:064;0:043;0:173;0:990;1:112;1:242] T = [ 1 3 ; 0; 0; 2 p 2 3 ] T = [0:305; 0:173; 0:130;1:091;1:242;1:284] T Table 2: Eigenvalues of J F for Rotational Controller 2 with distinct i . . 137 Appendix C Non-Uniform Gravity Field Coecients and Physical Properties of Mars To accurately model the non-uniform gravitational eld of Mars, we utilize the spheri- cal harmonic coecients and physical properties obtained from the Mars Global Surveyor [15]. The gravitational parameter, g , the mean equatorial radius,R e , and the angular rate of Mars are given in Table 3. Parameter Value g 4:2828380415705753 10 4 km 3 s 2 R e 3:3962 10 3 km 7:08823595918567 10 5 rad s 1 Table 3: Mars physical properties. Table 4 gives the zonal harmonic coecients and Table 5 gives the tesseral and sectoral harmonic coecients up to the fourth degree. The rst harmonics (k = 1) are not shown since they are zero. Zonal Coecients J 2 = 195:556398928615 10 5 J 3 = 3:14165422991657 10 5 J 4 =1:53315531432299 10 5 Table 4: Zonal harmonic coecients up to degree four (k = 4) for the planet Mars. . 138 Tesseral Coecients Sectoral Coecients C 2;1 =3:94480634264613 10 9 S 2;1 = 2:23593963786895 10 9 C 2;2 =5:43352869066853 10 5 S 2;2 = 3:20170504988189 10 5 C 2;3 = 0 S 2;3 = 0 C 2;4 = 0 S 2;4 = 0 C 3;1 = 4:21857224736813 10 6 S 3;1 = 2:71533881368686 10 5 C 3;2 =5:41812865445241 10 6 S 3;2 = 2:89865730476064 10 6 C 3;3 = 4:92733035642495 10 6 S 3;3 = 3:50193614477309 10 6 C 3;4 = 0 S 3;4 = 0 C 4;1 = 4:00779032010079 10 6 S 4;1 = 3:55154320241793 10 6 C 4;2 =2:29674292177281 10 7 S 4;2 =2:00384134071156 10 6 C 4;3 = 3:85511535680243 10 7 S 4;3 =1:61149835552153 10 8 C 4;4 = 2:05367185921771 10 9 S 4;4 =2:71737558012244 10 7 Table 5: Tesseral and sectoral harmonic coecients up to degree four (k = 4) for the planet Mars. . 139

Abstract (if available)

Abstract This dissertation develops in a unified manner a new and simple approach for the modeling and control of complex multibody systems. Complex multibody systems are those nonlinear mechanical systems consisting of individual subsystems that are describable by nonlinear differential equations. The interaction of the various subsystems is, in general, governed by nonlinear 'elements.' The nonlinear analysis of general multibody systems presents theoretically challenging problems in both its dynamics and controls aspects. The theoretical developments obtained in this dissertation permit the modeling of multibody systems so that no restrictions are imposed in its formulation, except that its physical model is continuous. The idea of permissible multibody control is developed, and an effective method is provided for generating nonlinear controllers that satisfy a general set of control objectives.

Linked assets

University of Southern California Dissertations and Theses

Conceptually similar

PDF

On the synthesis of controls for general nonlinear constrained mechanical systems

PDF

An analytical dynamics approach to the control of mechanical systems

PDF

Dynamic modeling and simulations of rigid-flexible coupled systems using quaternion dynamics

PDF

Nonlinear dynamics and nonlinear dynamical systems

PDF

Transient modeling, dynamic analysis, and feedback control of the Inductrack Maglev system

PDF

New approaches in modeling and control of dynamical systems

PDF

Nonlinear control of flexible rotating system with varying velocity

PDF

Dynamic analysis and control of one-dimensional distributed parameter systems

PDF

Flexible formation configuration for terrain following flight: formation keeping constraints

PDF

Robust control of periodically time-varying systems

PDF

Control of spacecraft with flexible structures using pulse-modulated thrusters

PDF

Data-driven H∞ loop-shaping controller design and stability of switched nonlinear feedback systems with average time-variation rate

PDF

Terrain-following motion of an autonomous agent as means of motion planning in the unknown environment

PDF

New approaches to satellite formation-keeping and the inverse problem of the calculus of variations

PDF

Analysis, design, and optimization of large-scale networks of dynamical systems

PDF

Studies into computational intelligence approaches for the identification of complex nonlinear systems

PDF

Organizing complex projects around critical skills, and the mitigation of risks arising from system dynamic behavior

PDF

Distributed adaptive control with application to heating, ventilation and air-conditioning systems

PDF

Modeling and dynamic analysis of coupled structure-moving subsystem problem

PDF

Periodic solutions of flexible systems under discontinuous control

Asset Metadata

Creator Schutte, Aaron Dane (author)

Core Title On the synthesis of dynamics and control for complex multibody systems

School Viterbi School of Engineering

Degree Doctor of Philosophy

Degree Program Aerospace Engineering

Degree Conferral Date 2009-05

Publication Date 03/13/2009

Defense Date 11/19/2008

Publisher University of Southern California (original), University of Southern California. Libraries (digital)

Tag general constrained systems,multibody dynamics and control,nonlinear systems,OAI-PMH Harvest,quaternions,rigid-body dynamics and control,spacecraft dynamics and control

Language English

Contributor Electronically uploaded by the author (provenance)

Advisor Udwadia, Firdaus E. (committee chair), Flashner, Henryk (committee member), Goodson, Troy (committee member), Gruntman, Michael A. (committee member), Yang, Bingen (committee member)

Creator Email aaron.d.schutte@aero.org,adschutte@gmail.com

Permanent Link (DOI) https://doi.org/10.25549/usctheses-m2020

Unique identifier UC1483276

Identifier etd-Schutte-2566 (filename),usctheses-m40 (legacy collection record id),usctheses-c127-216110 (legacy record id),usctheses-m2020 (legacy record id)

Legacy Identifier etd-Schutte-2566.pdf

Dmrecord 216110

Document Type Dissertation

Rights Schutte, Aaron Dane

Type texts

Source University of Southern California (contributing entity), University of Southern California Dissertations and Theses (collection)

Repository Name Libraries, University of Southern California

Repository Location Los Angeles, California

Repository Email cisadmin@lib.usc.edu